source: trunk/LMDZ.MARS/libf/phymars/write_output_mod.F90

Last change on this file was 3225, checked in by emillour, 10 months ago

Mars PCM:
Remove interactive checking with XIOS whether a field should be sent to it;
some yet unresolved issues arise when using this in mixed MPI-OpenMP mode...
EM

File size: 9.4 KB
RevLine 
[2932]1MODULE write_output_mod
[3092]2
[3160]3implicit none
[3092]4
[3160]5private
[2970]6
[3092]7INTERFACE write_output
8    MODULE PROCEDURE write_output_d0, write_output_d1, write_output_d2, &
9                     write_output_i0, write_output_i1, write_output_i2, &
10                     write_output_l0, write_output_l1, write_output_l2
11END INTERFACE write_output
[2932]12
[3160]13public write_output
[3092]14
[3160]15!-----------------------------------------------------------------------
16contains
17!-----------------------------------------------------------------------
[2932]18
[3160]19SUBROUTINE write_output_d0(field_name,title,units,field)
[3161]20! For a scalar
[3160]21
[2934]22#ifdef CPP_XIOS
[3160]23    use xios_output_mod, only: xios_is_active_field
24    use xios_output_mod, only: send_xios_field
[2934]25#endif
[3092]26
[3160]27implicit none
28
29include "dimensions.h"
30
31integer, parameter :: ngrid = 2 + (jjm - 1)*iim - 1/jjm
32character(*), intent(in) :: field_name
33character(*), intent(in) :: title
34character(*), intent(in) :: units
35real,         intent(in) :: field
36
[3225]37logical :: is_active ! For XIOS, should this field be sent or not
38
[3160]39call writediagfi(ngrid,field_name,title,units,0,(/field/))
[3092]40#ifdef CPP_XIOS
[3225]41!is_active=xios_is_active_field(field_name)
42is_active=.true.
[3055]43    ! only send the field to xios if the user asked for it
[3225]44    if (is_active) call send_xios_field(field_name,field)
[2932]45#endif
[3092]46
[3160]47END SUBROUTINE write_output_d0
[2932]48
[3160]49!-----------------------------------------------------------------------
[3092]50
[3160]51SUBROUTINE write_output_d1(field_name,title,units,field)
52! For a surface field
53
[2934]54#ifdef CPP_XIOS
[3161]55    use xios_output_mod, only: xios_is_active_field
56    use xios_output_mod, only: send_xios_field
[2934]57#endif
[3092]58
[3160]59implicit none
60
61include "dimensions.h"
62
63integer, parameter :: ngrid = 2 + (jjm - 1)*iim - 1/jjm
64character(*),       intent(in) :: field_name
65character(*),       intent(in) :: title
66character(*),       intent(in) :: units
67real, dimension(:), intent(in) :: field
68
[3225]69logical :: is_active ! For XIOS, should this field be sent or not
70
[3160]71call writediagfi(ngrid,field_name,title,units,2,field)
[3092]72#ifdef CPP_XIOS
[3225]73!is_active=xios_is_active_field(field_name)
74is_active=.true.
[3055]75    ! only send the field to xios if the user asked for it
[3225]76    if (is_active) call send_xios_field(field_name,field)
[2932]77#endif
[3092]78
[3160]79END SUBROUTINE write_output_d1
[2932]80
[3160]81!-----------------------------------------------------------------------
[3092]82
[3160]83SUBROUTINE write_output_d2(field_name,title,units,field)
84! For a "3D" horizontal-vertical field
85
[2934]86#ifdef CPP_XIOS
[3160]87    use xios_output_mod, only: xios_is_active_field
88    use xios_output_mod, only: send_xios_field
[2934]89#endif
[3161]90
[3160]91use comsoil_h,         only: nsoilmx
92use writediagsoil_mod, only: writediagsoil
[2932]93
[3160]94implicit none
95
96include "dimensions.h"
97
98integer, parameter :: ngrid = 2 + (jjm - 1)*iim - 1/jjm
99character(*),         intent(in) :: field_name
100character(*),         intent(in) :: title
101character(*),         intent(in) :: units
102real, dimension(:,:), intent(in) :: field
103
[3225]104logical :: is_active ! For XIOS, should this field be sent or not
105
[3160]106if (size(field,2) == nsoilmx) then
[2932]107    call writediagsoil(ngrid,field_name,title,units,3,field)
[3160]108else
109    call writediagfi(ngrid,field_name,title,units,3,field)
110endif
111
[3092]112#ifdef CPP_XIOS
[3225]113!is_active=xios_is_active_field(field_name)
114is_active=.true.
[3055]115    ! only send the field to xios if the user asked for it
[3225]116    if (is_active) call send_xios_field(field_name,field)
[2932]117#endif
[3092]118
[3160]119END SUBROUTINE write_output_d2
[2932]120
[3160]121!-----------------------------------------------------------------------
[3092]122
[3160]123SUBROUTINE write_output_i0(field_name,title,units,field)
[3161]124! For a scalar
[3160]125
[2970]126#ifdef CPP_XIOS
[3160]127    use xios_output_mod, only: xios_is_active_field
128    use xios_output_mod, only: send_xios_field
[2970]129#endif
[3092]130
[3160]131implicit none
132
133include "dimensions.h"
134
135integer, parameter :: ngrid = 2 + (jjm - 1)*iim - 1/jjm
136character(*), intent(in) :: field_name
137character(*), intent(in) :: title
138character(*), intent(in) :: units
139integer,      intent(in) :: field
140
[3225]141logical :: is_active ! For XIOS, should this field be sent or not
142
[3160]143call writediagfi(ngrid,field_name,title,units,0,(/real(field)/))
[3092]144#ifdef CPP_XIOS
[3225]145!is_active=xios_is_active_field(field_name)
146is_active=.true.
[3055]147    ! only send the field to xios if the user asked for it
[3225]148    if (is_active) call send_xios_field(field_name,(/real(field)/))
[2970]149#endif
[3092]150
[3160]151END SUBROUTINE write_output_i0
[2970]152
[3160]153!-----------------------------------------------------------------------
[3092]154
[3160]155SUBROUTINE write_output_i1(field_name,title,units,field)
156! For a surface field
157
[2970]158#ifdef CPP_XIOS
[3160]159    use xios_output_mod, only: xios_is_active_field
160    use xios_output_mod, only: send_xios_field
[2970]161#endif
[3092]162
[3160]163implicit none
164
165include "dimensions.h"
166
167integer, parameter :: ngrid = 2 + (jjm - 1)*iim - 1/jjm
168character(*),          intent(in) :: field_name
169character(*),          intent(in) :: title
170character(*),          intent(in) :: units
171integer, dimension(:), intent(in) :: field
172
[3225]173logical :: is_active ! For XIOS, should this field be sent or not
174
[3160]175call writediagfi(ngrid,field_name,title,units,2,real(field))
[3092]176#ifdef CPP_XIOS
[3225]177!is_active=xios_is_active_field(field_name)
178is_active=.true.
[3055]179    ! only send the field to xios if the user asked for it
[3225]180    if (is_active) call send_xios_field(field_name,real(field))
[2970]181#endif
[3092]182
[3160]183END SUBROUTINE write_output_i1
[2970]184
[3160]185!-----------------------------------------------------------------------
[3092]186
[3160]187SUBROUTINE write_output_i2(field_name,title,units,field)
188! For a "3D" horizontal-vertical field
189
[2970]190#ifdef CPP_XIOS
[3160]191    use xios_output_mod, only: xios_is_active_field
192    use xios_output_mod, only: send_xios_field
[2970]193#endif
194
[3160]195use comsoil_h,         only: nsoilmx
196use writediagsoil_mod, only: writediagsoil
197
198implicit none
199
200include "dimensions.h"
201
202integer, parameter :: ngrid = 2 + (jjm - 1)*iim - 1/jjm
203character(*),            intent(in) :: field_name
204character(*),            intent(in) :: title
205character(*),            intent(in) :: units
206integer, dimension(:,:), intent(in) :: field
207
[3225]208logical :: is_active ! For XIOS, should this field be sent or not
209
[3160]210if (size(field,2) == nsoilmx) then
[2970]211    call writediagsoil(ngrid,field_name,title,units,3,real(field))
[3160]212else
213    call writediagfi(ngrid,field_name,title,units,3,real(field))
214endif
[3092]215#ifdef CPP_XIOS
[3225]216!is_active=xios_is_active_field(field_name)
217is_active=.true.
[3055]218    ! only send the field to xios if the user asked for it
[3225]219    if (is_active) call send_xios_field(field_name,real(field))
[2970]220#endif
[3092]221
[3160]222END SUBROUTINE write_output_i2
[2970]223
[3160]224!-----------------------------------------------------------------------
[3092]225
[3160]226SUBROUTINE write_output_l0(field_name,title,units,field)
[3161]227! For a scalar
[3160]228 
[2970]229#ifdef CPP_XIOS
[3160]230    use xios_output_mod, only: xios_is_active_field
231    use xios_output_mod, only: send_xios_field
[2970]232#endif
233
[3160]234implicit none
[3092]235
[3160]236include "dimensions.h"
237
238integer, parameter :: ngrid = 2 + (jjm - 1)*iim - 1/jjm
239character(*), intent(in) :: field_name
240character(*), intent(in) :: title
241character(*), intent(in) :: units
242logical,      intent(in) :: field
243! Local argument used to convert logical to real array
244real, dimension(1) :: field_real
[3225]245logical :: is_active ! For XIOS, should this field be sent or not
[3160]246
247field_real = 0.
248if (field) field_real = 1.
249
250call writediagfi(ngrid,field_name,title,units,0,field_real)
[3092]251#ifdef CPP_XIOS
[3225]252!is_active=xios_is_active_field(field_name)
253is_active=.true.
[3055]254    ! only send the field to xios if the user asked for it
[3225]255    if (is_active) call send_xios_field(field_name,field_real)
[2970]256#endif
[3092]257
[3160]258END SUBROUTINE write_output_l0
[2970]259
[3160]260!-----------------------------------------------------------------------
[3092]261
[3160]262SUBROUTINE write_output_l1(field_name,title,units,field)
263! For a surface field
264
[2970]265#ifdef CPP_XIOS
[3160]266    use xios_output_mod, only: xios_is_active_field
267    use xios_output_mod, only: send_xios_field
[2970]268#endif
[3160]269
270implicit none
271
272include "dimensions.h"
273
274integer, parameter :: ngrid = 2 + (jjm - 1)*iim - 1/jjm
275character(*),          intent(in) :: field_name
276character(*),          intent(in) :: title
277character(*),          intent(in) :: units
278logical, dimension(:), intent(in) :: field
[2970]279! Local argument used to convert logical to real
[3160]280real, dimension(ngrid) :: field_real
[3225]281logical :: is_active ! For XIOS, should this field be sent or not
[2970]282
[3160]283field_real = 0.
284where (field) field_real = 1.
[3092]285
[3160]286call writediagfi(ngrid,field_name,title,units,2,field_real)
[3092]287#ifdef CPP_XIOS
[3225]288!is_active=xios_is_active_field(field_name)
289is_active=.true.
[3055]290    ! only send the field to xios if the user asked for it
[3225]291    if (is_active) call send_xios_field(field_name,field_real)
[2970]292#endif
[3092]293
[3160]294END SUBROUTINE write_output_l1
[2970]295
[3160]296!-----------------------------------------------------------------------
[3092]297
[3160]298SUBROUTINE write_output_l2(field_name,title,units,field)
299! For a "3D" horizontal-vertical field
300
[2970]301#ifdef CPP_XIOS
[3160]302    use xios_output_mod, only: xios_is_active_field
303    use xios_output_mod, only: send_xios_field
[2970]304#endif
[3160]305
306use comsoil_h,         only: nsoilmx
307use writediagsoil_mod, only: writediagsoil
308
309implicit none
310
311include "dimensions.h"
312
313integer, parameter :: ngrid = 2 + (jjm - 1)*iim - 1/jjm
314character(*),            intent(in) :: field_name
315character(*),            intent(in) :: title
316character(*),            intent(in) :: units
317logical, dimension(:,:), intent(in) :: field
[2970]318! Local argument used to convert logical to real
[3160]319real, allocatable, dimension(:,:) :: field_real
[3225]320logical :: is_active ! For XIOS, should this field be sent or not
[2970]321
[3160]322allocate(field_real(size(field,1),size(field,2)))
323field_real = 0.
324where (field) field_real = 1.
[2970]325
[3160]326if (size(field,2) == nsoilmx) then
[2970]327    call writediagsoil(ngrid,field_name,title,units,3,field_real)
[3160]328else
329    call writediagfi(ngrid,field_name,title,units,3,field_real)
330endif
[2970]331
[3092]332#ifdef CPP_XIOS
[3225]333!is_active=xios_is_active_field(field_name)
334is_active=.true.
[3055]335    ! only send the field to xios if the user asked for it
[3225]336    if (is_active) call send_xios_field(field_name,field_real)
[2970]337#endif
[3092]338
[3160]339deallocate(field_real)
[2976]340
[3160]341END SUBROUTINE write_output_l2
[2970]342
[2932]343END MODULE write_output_mod
Note: See TracBrowser for help on using the repository browser.