source: codes/icosagcm/devel/src/output/xios_mod.F90 @ 905

Last change on this file since 905 was 905, checked in by dubos, 5 years ago

devel : XIOS output for unstructured mesh
File size: 12.0 KB
Line 
1MODULE xios_mod
2 
3#ifdef CPP_USING_XIOS
4  USE xios
5#endif
6
7  USE prec, ONLY : rstd
8  USE field_mod, ONLY : t_field, field_T, field_U, field_Z
9  USE domain_mod, ONLY : t_domain, t_cellset, domain, ndomain, mesh_loc
10  USE grid_param, ONLY : grid_type, grid_unst, grid_ico
11  IMPLICIT NONE   
12  PRIVATE
13  SAVE
14
15  LOGICAL :: using_xios
16
17  INTEGER :: ncell_i, ncell_v, ncell_e
18!$OMP THREADPRIVATE(ncell_i, ncell_v, ncell_e)
19
20  PUBLIC :: using_xios, xios_init, &
21       xios_init_write_field,  xios_write_field_finalize, &
22       xios_write_field, xios_read_field
23
24#ifdef CPP_USING_XIOS
25
26  PUBLIC :: xios_timestep,    & 
27       xios_set_file_attr, xios_set_fieldgroup_attr, &
28       xios_set_filegroup_attr, xios_get_axis_attr, &
29       xios_send_field, xios_read_var, &
30       xios_update_calendar, xios_set_context
31 
32CONTAINS
33 
34 SUBROUTINE xios_init
35   USE mpipara, ONLY : comm_icosa
36   TYPE(xios_context) :: ctx_hdl
37
38   using_xios=.TRUE.
39   CALL xios_context_initialize("icosagcm",comm_icosa)
40   CALL xios_get_handle("icosagcm",ctx_hdl)
41   CALL xios_set_current_context(ctx_hdl)   
42
43 END SUBROUTINE xios_init
44
45 SUBROUTINE xios_init_write_field
46   USE disvert_mod, ONLY : presnivs
47   USE time_mod, ONLY : dt, itau_out
48   USE grid_param, ONLY : llm, nqtot
49   USE mpi_mod, ONLY : MPI_INTEGER
50
51   TYPE(xios_context) :: ctx_hdl
52   TYPE(xios_duration)      :: dtime
53   REAL(rstd) :: lev_value(llm)
54   REAL(rstd) :: lev_valuep1(llm+1)
55   REAL(rstd) :: nq_value(nqtot)
56   INTEGER :: l, ncell, ncell_tot, displ
57   REAL(rstd),ALLOCATABLE    :: lon(:), lat(:), bounds_lon(:,:), bounds_lat(:,:)
58   INTEGER, ALLOCATABLE      :: ind_glo(:)
59   TYPE(t_domain),POINTER :: d
60   
61   !$OMP BARRIER
62   !$OMP MASTER
63   !   CALL xios_context_initialize("icosagcm",comm_icosa)
64   CALL xios_get_handle("icosagcm",ctx_hdl)
65   CALL xios_set_current_context(ctx_hdl)
66   lev_value(:) = (/ (l,l=1,llm) /)     
67   lev_valuep1(:) = (/ (l,l=1,llm+1) /)     
68   nq_value(:) = (/ (l,l=1,nqtot) /)     
69   CALL xios_set_axis_attr("lev",n_glo=llm ,value=lev_value) ;
70   CALL xios_set_axis_attr("levp1",n_glo=llm+1 ,value=lev_valuep1) ;
71   CALL xios_set_axis_attr("nq",n_glo=nqtot, value=nq_value) ;
72   CALL xios_set_axis_attr("presnivs_mb",n_glo=llm, value=presnivs/100., unit="mb") ;
73
74   !------------------ primal cells ------------------
75   CALL collect_bounds(6, mesh_loc%primal_own)
76   ncell_i=ncell   
77   CALL xios_set_domaingroup_attr("i",ni_glo=ncell_tot, ibegin=displ, ni=ncell)
78   CALL xios_set_domaingroup_attr("i", data_dim=1, type='unstructured' , nvertex=6, i_index=ind_glo)
79   CALL xios_set_domaingroup_attr("i",lonvalue_1d=lon, latvalue_1d=lat, bounds_lon_1d=bounds_lon, bounds_lat_1d=bounds_lat)
80   DEALLOCATE(lon, lat, bounds_lon, bounds_lat,ind_glo) 
81   
82   !--------------------- dual cells ------------------
83   
84   CALL collect_bounds(3, mesh_loc%dual_own)
85   ncell_v=ncell
86   CALL xios_set_domain_attr("v",ni_glo=ncell_tot, ibegin=displ, ni=ncell)
87   CALL xios_set_domain_attr("v", data_dim=1, type='unstructured' , nvertex=3)
88   CALL xios_set_domain_attr("v",lonvalue_1d=lon, latvalue_1d=lat, bounds_lon_1d=bounds_lon, bounds_lat_1d=bounds_lat)
89   DEALLOCATE(lon, lat, bounds_lon, bounds_lat,ind_glo) 
90   
91   !---------------------- edges -----------------------
92   
93   CALL collect_bounds(2, mesh_loc%edge_own)
94   ncell_e=ncell
95   CALL xios_set_domain_attr("u",ni_glo=ncell_tot, ibegin=displ, ni=ncell)
96   CALL xios_set_domain_attr("u", data_dim=1, type='unstructured' , nvertex=2, i_index=ind_glo)
97   CALL xios_set_domain_attr("u",lonvalue_1d=lon, latvalue_1d=lat, bounds_lon_1d=bounds_lon, bounds_lat_1d=bounds_lat)
98   DEALLOCATE(lon, lat, bounds_lon, bounds_lat,ind_glo) 
99   
100   dtime%second=dt
101   CALL xios_set_timestep(dtime)
102   
103   CALL xios_set_fieldgroup_attr("standard_output", freq_op=itau_out*xios_timestep, freq_offset=(itau_out-1)*xios_timestep)
104   
105   CALL xios_close_context_definition()
106   !$OMP END MASTER
107   !$OMP BARRIER
108   
109 CONTAINS
110   
111   SUBROUTINE collect_bounds(nvert, cells)
112     USE mpipara, ONLY : comm_icosa, mpi_size, mpi_rank
113     INTEGER, INTENT(IN) :: nvert
114     TYPE(t_cellset)     :: cells(:)
115     INTEGER :: i, ind, n_beg, n_end, ierr, ncell_glo(0:mpi_size-1)
116     ncell = SUM(cells%ncell)
117     CALL MPI_ALLGATHER(ncell,1,MPI_INTEGER, &
118          ncell_glo,1, MPI_INTEGER, comm_icosa, ierr)
119     displ=0
120     DO i=1,mpi_rank
121        displ=displ+ncell_glo(i-1)
122     ENDDO
123     ncell_tot=sum(ncell_glo(:))
124     
125     ALLOCATE(lon(ncell), lat(ncell), ind_glo(ncell))
126     ALLOCATE(bounds_lon(nvert,ncell), bounds_lat(nvert,ncell))
127     
128     n_beg=0
129     DO ind=1,ndomain
130        n_end = n_beg + cells(ind)%ncell
131        ind_glo(n_beg+1:n_end) = cells(ind)%ind_glo(:)
132        lon(n_beg+1:n_end) = cells(ind)%lon(:)
133        lat(n_beg+1:n_end) = cells(ind)%lat(:)
134        bounds_lon(:,n_beg+1:n_end) = cells(ind)%bnds_lon(:,:)
135        bounds_lat(:,n_beg+1:n_end) = cells(ind)%bnds_lat(:,:)
136        n_beg = n_end
137     END DO
138   END SUBROUTINE collect_bounds
139   
140 END SUBROUTINE xios_init_write_field
141 
142 SUBROUTINE xios_write_field(name,field)
143   CHARACTER(LEN=*),INTENT(IN) :: name
144   TYPE(t_field), POINTER :: field(:)
145   TYPE(t_cellset), POINTER :: cells(:)
146   INTEGER :: ncells
147
148!$OMP BARRIER
149!$OMP MASTER
150
151   SELECT CASE(field(1)%field_type)
152   CASE(field_T)
153      cells => mesh_loc%primal_own
154      ncells = ncell_i
155   CASE(field_U)
156      cells => mesh_loc%edge_own
157      ncells = ncell_e
158   CASE(field_Z)
159      cells => mesh_loc%dual_own
160      ncells = ncell_v
161   END SELECT
162
163   IF (field(1)%ndim>4) THEN
164      PRINT *, "xios_write_field : dimension > 4 are not supported for now"
165   ELSE
166      CALL xios_write_field_gen(name, field, cells, & 
167           ncells, field(1)%dim3, field(1)%dim4)
168   END IF
169   
170!$OMP END MASTER
171!$OMP BARRIER
172     
173 END SUBROUTINE xios_write_field
174
175 SUBROUTINE xios_read_field(name,field)
176   CHARACTER(LEN=*),INTENT(IN) :: name
177   TYPE(t_field), POINTER :: field(:)
178   TYPE(t_cellset), POINTER :: cells(:)
179   INTEGER :: ncells
180
181!$OMP BARRIER
182!$OMP MASTER
183
184   SELECT CASE(field(1)%field_type)
185   CASE(field_T)
186      cells => mesh_loc%primal_own
187      ncells = ncell_i
188   CASE(field_U)
189      cells => mesh_loc%edge_own
190      ncells = ncell_e
191   CASE(field_Z)
192      cells => mesh_loc%dual_own
193      ncells = ncell_v
194   END SELECT
195
196   IF (field(1)%ndim>4) THEN
197      PRINT *, "xios_read_field : dimension > 4 are not supported for now"
198   ELSE
199      CALL xios_read_field_hex(name, field, cells, & 
200           ncells, field(1)%dim3, field(1)%dim4)
201   END IF
202   
203!$OMP END MASTER
204!$OMP BARRIER
205     
206 END SUBROUTINE xios_read_field
207
208 SUBROUTINE xios_write_field_gen(name, field, cells, ncell_tot, nlev, nq)
209   CHARACTER(LEN=*),INTENT(IN) :: name
210   TYPE(t_field) :: field(:)
211   TYPE(t_cellset), TARGET :: cells(:)
212   INTEGER,INTENT(IN) :: ncell_tot, nlev, nq
213   REAL(rstd) :: field_tmp(ncell_tot, nlev, nq)
214   TYPE(t_cellset), POINTER :: cellset
215   INTEGER :: ind, n_beg, n, ij, sgn
216   LOGICAL :: signed
217
218   IF(ALLOCATED(cells(1)%sgn)) THEN
219      signed=.TRUE.
220   ELSE
221      signed=.FALSE.
222      sgn=1
223   END IF
224
225   n_beg=0
226   DO ind=1,ndomain
227      cellset => cells(ind)
228      DO n=1, cellset%ncell
229         ij = cellset%ij(n)
230         IF(signed) sgn = cellset%sgn(n)
231         SELECT CASE(grid_type)
232         CASE(grid_ico)
233            SELECT CASE(field(1)%ndim)
234            CASE(2)
235               field_tmp(n_beg+n,1,1) = sgn*field(ind)%rval2d(ij)
236            CASE(3)
237               field_tmp(n_beg+n,:,1) = sgn*field(ind)%rval3d(ij,:)
238            CASE(4)
239               field_tmp(n_beg+n,:,:) = sgn*field(ind)%rval4d(ij,:,:)
240            END SELECT
241         CASE(grid_unst)
242            SELECT CASE(field(1)%ndim)
243            CASE(2)
244               field_tmp(n_beg+n,1,1) = sgn*field(ind)%rval2d(ij)
245            CASE(3)
246               field_tmp(n_beg+n,:,1) = sgn*field(ind)%rval3d(:,ij)
247            CASE(4)
248               field_tmp(n_beg+n,:,:) = sgn*field(ind)%rval4d(:,ij,:)
249            END SELECT
250         END SELECT
251      END DO
252      n_beg = n_beg + cellset%ncell
253   END DO
254   CALL xios_send_field(name,field_tmp) 
255   
256 END SUBROUTINE xios_write_field_gen
257
258SUBROUTINE xios_read_field_hex(name, field, cells, ncell_tot, nlev, nq)
259   CHARACTER(LEN=*),INTENT(IN) :: name
260   TYPE(t_field) :: field(:)
261   TYPE(t_cellset), TARGET :: cells(:)
262   INTEGER,INTENT(IN) :: ncell_tot, nlev, nq
263
264   REAL(rstd) :: field_tmp(ncell_tot,nlev,nq)
265   TYPE(t_cellset), POINTER :: cellset
266   INTEGER :: ind, n_beg, n_end, n, ij, sgn
267   LOGICAL :: signed
268
269   CALL xios_recv_field(name,field_tmp)
270
271   IF(ALLOCATED(cells(1)%sgn)) THEN
272      signed=.TRUE.
273   ELSE
274      signed=.FALSE.
275      sgn=1
276   END IF
277
278   n_beg=0
279   DO ind=1,ndomain
280      cellset => cells(ind)
281      n_end = n_beg + cellset%ncell
282      DO n=1, cellset%ncell
283         ij = cellset%ij(n)
284         IF(signed) sgn = cellset%sgn(n)
285         SELECT CASE(field(1)%ndim)
286         CASE(2)
287            field(ind)%rval2d(ij) = sgn*field_tmp(n_beg+n,1,1) 
288         CASE(3)
289            field(ind)%rval3d(ij,:) = sgn*field_tmp(n_beg+n,:,1)
290         CASE(4)
291            field(ind)%rval4d(ij,:,:) = sgn*field_tmp(n_beg+n,:,:)
292         END SELECT
293      END DO
294   END DO
295 END SUBROUTINE xios_read_field_hex
296
297 SUBROUTINE xios_read_var(name,field)
298   USE prec
299   USE transfert_mod
300   CHARACTER(LEN=*),INTENT(IN) :: name
301   REAL(rstd), INTENT(OUT) :: field
302   !$OMP MASTER
303   CALL xios_recv_field(name,field)
304   !$OMP END MASTER
305   CALL bcast_omp(field)
306 END SUBROUTINE
307 
308 SUBROUTINE xios_write_field_finalize
309!$OMP BARRIER
310!$OMP MASTER
311   CALL xios_context_finalize
312!$OMP END MASTER
313!$OMP BARRIER
314 END SUBROUTINE xios_write_field_finalize
315
316 SUBROUTINE xios_set_context
317  TYPE(xios_context) :: ctx_hdl
318
319!$OMP MASTER
320   CALL xios_get_handle("icosagcm",ctx_hdl)
321   CALL xios_set_current_context(ctx_hdl)
322!$OMP END MASTER
323
324  END SUBROUTINE xios_set_context
325
326
327#else
328 
329
330INTERFACE xios_send_field
331  MODULE PROCEDURE xios_send_field_scalar, xios_send_field_1d
332END INTERFACE  xios_send_field
333
334INTEGER,PARAMETER :: xios_timestep=1
335
336CONTAINS 
337 
338 
339  SUBROUTINE xios_init
340    using_xios=.FALSE.   
341  END SUBROUTINE xios_init
342 
343  SUBROUTINE xios_send_field_scalar(name,field)
344    CHARACTER(LEN=*),INTENT(IN) :: name
345    REAL,INTENT(IN) :: field
346  END SUBROUTINE xios_send_field_scalar 
347
348  SUBROUTINE xios_send_field_1d(name,field)
349    CHARACTER(LEN=*),INTENT(IN) :: name
350    REAL,INTENT(IN) :: field(:)
351  END SUBROUTINE xios_send_field_1d 
352 
353  SUBROUTINE xios_write_field(name,field)
354  USE field_mod
355   CHARACTER(LEN=*),INTENT(IN) :: name
356   TYPE(t_field), POINTER :: field(:)
357  END SUBROUTINE xios_write_field
358
359  SUBROUTINE xios_read_field(name,field)
360  USE field_mod
361   CHARACTER(LEN=*),INTENT(IN) :: name
362   TYPE(t_field), POINTER :: field(:)
363  END SUBROUTINE xios_read_field
364
365 SUBROUTINE xios_read_var(name,field)
366   USE prec
367   CHARACTER(LEN=*),INTENT(IN) :: name
368   REAL(rstd), INTENT(OUT) :: field
369 END SUBROUTINE
370 
371  SUBROUTINE xios_update_calendar(step)
372   INTEGER, INTENT(IN):: step 
373  END SUBROUTINE xios_update_calendar
374
375  SUBROUTINE xios_write_field_finalize
376  END SUBROUTINE xios_write_field_finalize
377 
378  SUBROUTINE xios_init_write_field
379  END SUBROUTINE xios_init_write_field 
380 
381  SUBROUTINE xios_set_context
382  END SUBROUTINE xios_set_context
383 
384  SUBROUTINE xios_set_fieldgroup_attr(name,enabled,freq_op)
385    CHARACTER(LEN=*) :: name
386    LOGICAL,OPTIONAL          :: enabled
387    INTEGER,OPTIONAL          :: freq_op
388  END SUBROUTINE xios_set_fieldgroup_attr
389
390  SUBROUTINE xios_set_filegroup_attr(name,enabled)
391    CHARACTER(LEN=*) :: name
392    LOGICAL,OPTIONAL          :: enabled
393  END SUBROUTINE xios_set_filegroup_attr
394
395  SUBROUTINE xios_set_file_attr(id,name,mode,enabled, output_freq)
396    CHARACTER(LEN=*) :: id
397    CHARACTER(LEN=*),OPTIONAL :: name, mode
398    LOGICAL,OPTIONAL          :: enabled
399    INTEGER,OPTIONAL          :: output_freq
400  END SUBROUTINE xios_set_file_attr
401
402  SUBROUTINE xios_get_axis_attr(name,n_glo,value)
403    CHARACTER(LEN=*) :: name
404    INTEGER,OPTIONAL          :: n_glo
405    REAL,OPTIONAL             :: value(:)
406  END SUBROUTINE xios_get_axis_attr
407
408  SUBROUTINE xios_set_axis_attr(id,n_glo,value)
409    CHARACTER(LEN=*) :: id
410    INTEGER,OPTIONAL          :: n_glo
411    REAL,OPTIONAL             :: value(:)
412  END SUBROUTINE xios_set_axis_attr
413
414#endif 
415
416END MODULE xios_mod
Note: See TracBrowser for help on using the repository browser.