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

Last change on this file since 993 was 992, checked in by rpennel, 5 years ago

devel : add reading metrics from file at init (folllowing revs @882, @887, @888)

[ !! needs to rewrite new subroutine in xios_mod using cellset strcuture ]

File size: 15.5 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_init_write_field_input, &
22       xios_write_field_finalize, &
23       xios_write_field, xios_read_field
24
25#ifdef CPP_USING_XIOS
26
27  PUBLIC :: xios_timestep,    & 
28       xios_set_file_attr, xios_set_fieldgroup_attr, &
29       xios_set_filegroup_attr, xios_get_axis_attr, &
30       xios_send_field, xios_read_var, &
31       xios_update_calendar, xios_set_context, xios_set_context_input
32 
33CONTAINS
34 
35 SUBROUTINE xios_init
36   USE mpipara, ONLY : comm_icosa
37   TYPE(xios_context) :: ctx_hdl
38
39   using_xios=.TRUE.
40   CALL xios_context_initialize("icosagcm",comm_icosa)
41   CALL xios_context_initialize("icosagcm_input",comm_icosa)
42   CALL xios_set_context
43
44 END SUBROUTINE xios_init
45
46 SUBROUTINE xios_init_write_field
47   USE disvert_mod, ONLY : presnivs
48   USE time_mod, ONLY : dt, itau_out
49   USE grid_param, ONLY : llm, nqtot
50   USE mpi_mod, ONLY : MPI_INTEGER
51
52   TYPE(xios_context) :: ctx_hdl
53   TYPE(xios_duration)      :: dtime
54   REAL(rstd) :: lev_value(llm)
55   REAL(rstd) :: lev_valuep1(llm+1)
56   REAL(rstd) :: nq_value(nqtot)
57   INTEGER :: l, ncell, ncell_tot, displ
58   REAL(rstd),ALLOCATABLE    :: lon(:), lat(:), bounds_lon(:,:), bounds_lat(:,:)
59   INTEGER, ALLOCATABLE      :: ind_glo(:)
60   TYPE(t_domain),POINTER :: d
61   
62   CALL xios_set_context
63   !$OMP BARRIER
64   !$OMP MASTER
65   !   CALL xios_context_initialize("icosagcm",comm_icosa)
66   !   CALL xios_get_handle("icosagcm",ctx_hdl)
67   !   CALL xios_set_current_context(ctx_hdl)
68   lev_value(:) = (/ (l,l=1,llm) /)     
69   lev_valuep1(:) = (/ (l,l=1,llm+1) /)     
70   nq_value(:) = (/ (l,l=1,nqtot) /)     
71   CALL xios_set_axis_attr("lev",n_glo=llm ,value=lev_value) ;
72   CALL xios_set_axis_attr("levp1",n_glo=llm+1 ,value=lev_valuep1) ;
73   CALL xios_set_axis_attr("nq",n_glo=nqtot, value=nq_value) ;
74   CALL xios_set_axis_attr("presnivs_mb",n_glo=llm, value=presnivs/100., unit="mb") ;
75
76   !------------------ primal cells ------------------
77   CALL collect_bounds(6, mesh_loc%primal_own)
78   ncell_i=ncell   
79   CALL xios_set_domaingroup_attr("i",ni_glo=ncell_tot, ibegin=displ, ni=ncell)
80   CALL xios_set_domaingroup_attr("i", data_dim=1, type='unstructured' , nvertex=6, i_index=ind_glo)
81   CALL xios_set_domaingroup_attr("i",lonvalue_1d=lon, latvalue_1d=lat, bounds_lon_1d=bounds_lon, bounds_lat_1d=bounds_lat)
82   DEALLOCATE(lon, lat, bounds_lon, bounds_lat,ind_glo) 
83   
84   !--------------------- dual cells ------------------
85   
86   CALL collect_bounds(3, mesh_loc%dual_own)
87   ncell_v=ncell
88   CALL xios_set_domain_attr("v",ni_glo=ncell_tot, ibegin=displ, ni=ncell)
89   CALL xios_set_domain_attr("v", data_dim=1, type='unstructured' , nvertex=3)
90   CALL xios_set_domain_attr("v",lonvalue_1d=lon, latvalue_1d=lat, bounds_lon_1d=bounds_lon, bounds_lat_1d=bounds_lat)
91   DEALLOCATE(lon, lat, bounds_lon, bounds_lat,ind_glo) 
92   
93   !---------------------- edges -----------------------
94   
95   CALL collect_bounds(2, mesh_loc%edge_own)
96   ncell_e=ncell
97   CALL xios_set_domain_attr("u",ni_glo=ncell_tot, ibegin=displ, ni=ncell)
98   CALL xios_set_domain_attr("u", data_dim=1, type='unstructured' , nvertex=2, i_index=ind_glo)
99   CALL xios_set_domain_attr("u",lonvalue_1d=lon, latvalue_1d=lat, bounds_lon_1d=bounds_lon, bounds_lat_1d=bounds_lat)
100   DEALLOCATE(lon, lat, bounds_lon, bounds_lat,ind_glo) 
101   
102   dtime%second=dt
103   CALL xios_set_timestep(dtime)
104   
105   CALL xios_set_fieldgroup_attr("standard_output", freq_op=itau_out*xios_timestep, freq_offset=(itau_out-1)*xios_timestep)
106   
107   CALL xios_close_context_definition()
108   !$OMP END MASTER
109   !$OMP BARRIER
110   
111 CONTAINS
112   
113   SUBROUTINE collect_bounds(nvert, cells)
114     USE mpipara, ONLY : comm_icosa, mpi_size, mpi_rank
115     INTEGER, INTENT(IN) :: nvert
116     TYPE(t_cellset)     :: cells(:)
117     INTEGER :: i, ind, n_beg, n_end, ierr, ncell_glo(0:mpi_size-1)
118     ncell = SUM(cells%ncell)
119     CALL MPI_ALLGATHER(ncell,1,MPI_INTEGER, &
120          ncell_glo,1, MPI_INTEGER, comm_icosa, ierr)
121     displ=0
122     DO i=1,mpi_rank
123        displ=displ+ncell_glo(i-1)
124     ENDDO
125     ncell_tot=sum(ncell_glo(:))
126     
127     ALLOCATE(lon(ncell), lat(ncell), ind_glo(ncell))
128     ALLOCATE(bounds_lon(nvert,ncell), bounds_lat(nvert,ncell))
129     
130     n_beg=0
131     DO ind=1,ndomain
132        n_end = n_beg + cells(ind)%ncell
133        ind_glo(n_beg+1:n_end) = cells(ind)%ind_glo(:)
134        lon(n_beg+1:n_end) = cells(ind)%lon(:)
135        lat(n_beg+1:n_end) = cells(ind)%lat(:)
136        bounds_lon(:,n_beg+1:n_end) = cells(ind)%bnds_lon(:,:)
137        bounds_lat(:,n_beg+1:n_end) = cells(ind)%bnds_lat(:,:)
138        n_beg = n_end
139     END DO
140   END SUBROUTINE collect_bounds
141   
142 END SUBROUTINE xios_init_write_field
143 
144 SUBROUTINE xios_init_write_field_input
145   USE disvert_mod,        ONLY : presnivs
146   USE time_mod,           ONLY : dt, itau_out
147   USE grid_param,         ONLY : llm, nqtot
148   USE mpi_mod,            ONLY : MPI_INTEGER
149   USE icosa,              ONLY : getin
150   USE mpipara,            ONLY : comm_icosa, mpi_rank, mpi_size
151   USE spherical_geom_mod, ONLY : xyz2lonlat
152   USE genmod
153 !USE genmod
154 !USE mpipara
155 !USE xios
156 !USE grid_param
157 !USE domain_mod
158 !USE dimensions
159 !USE spherical_geom_mod
160 !USE geometry
161 !USE mpi_mod
162 !USE time_mod
163 !USE metric, ONLY : vup,vdown, cell_glo
164 !USE icosa,ONLY  : getin
165 !IMPLICIT NONE
166  TYPE(xios_context) :: ctx_hdl
167  TYPE(xios_duration)      :: dtime
168  INTEGER :: ncell, ncell_tot, ncell_glo(0:mpi_size-1), displ
169  INTEGER :: ind, i,j,k,l,ij, ierr
170  REAL(rstd),ALLOCATABLE    :: lon(:), lat(:), bounds_lon(:,:), bounds_lat(:,:)
171  INTEGER, ALLOCATABLE      :: ind_glo(:)
172  TYPE(t_domain),POINTER :: d
173  CHARACTER(len=255) :: etat0_type
174  LOGICAL :: read_metric_
175
176   CALL xios_set_context_input
177!$OMP BARRIER
178!$OMP MASTER
179   
180   ncell=0
181   DO ind=1,ndomain
182     d=>domain(ind)
183       
184     DO j=d%jj_begin,d%jj_end
185       DO i=d%ii_begin,d%ii_end
186         IF (domain(ind)%own(i,j)) ncell=ncell+1
187       ENDDO
188     ENDDO
189   ENDDO     
190   ncell_i=ncell
191   
192   CALL MPI_ALLGATHER(ncell,1,MPI_INTEGER,ncell_glo,1,MPI_INTEGER,comm_icosa,ierr)
193
194   displ=0
195   DO i=1,mpi_rank
196     displ=displ+ncell_glo(i-1)
197   ENDDO
198
199   ncell_tot=sum(ncell_glo(:))
200   
201   ALLOCATE(lon(ncell), lat(ncell), bounds_lon(0:5,ncell), bounds_lat(0:5,ncell), ind_glo(ncell)) 
202   
203   ncell=0
204   DO ind=1,ndomain
205     d=>domain(ind)
206       
207     DO j=d%jj_begin,d%jj_end
208       DO i=d%ii_begin,d%ii_end
209         IF (domain(ind)%own(i,j)) THEN
210           ncell=ncell+1
211           CALL xyz2lonlat(d%xyz(:,i,j),lon(ncell),lat(ncell))
212           lon(ncell)=lon(ncell)*180/Pi
213           lat(ncell)=lat(ncell)*180/Pi
214           DO k=0,5
215             CALL xyz2lonlat(d%vertex(:,k,i,j),bounds_lon(k,ncell), bounds_lat(k,ncell))
216             bounds_lat(k,ncell)=bounds_lat(k,ncell)*180/Pi
217             bounds_lon(k,ncell)=bounds_lon(k,ncell)*180/Pi
218           ENDDO
219           ind_glo(ncell)=domain(ind)%assign_cell_glo(i,j)-1 
220         ENDIF
221       ENDDO
222     ENDDO
223   ENDDO         
224
225   CALL xios_set_domain_attr("i",ni_glo=ncell_tot, ibegin=displ, ni=ncell)
226   CALL xios_set_domain_attr("i", data_dim=1, type='unstructured' , nvertex=6, i_index=ind_glo)
227   CALL xios_set_domain_attr("i",lonvalue_1d=lon, latvalue_1d=lat, bounds_lon_1d=bounds_lon, bounds_lat_1d=bounds_lat)
228   
229   DEALLOCATE(lon, lat, bounds_lon, bounds_lat,ind_glo) 
230
231   dtime%second=1
232   CALL xios_set_timestep(dtime)
233!$OMP END MASTER
234
235   CALL getin('etat0',etat0_type)
236   CALL getin('read_metric', read_metric_) 
237
238!$OMP MASTER
239 
240   CALL xios_set_file_attr('start', enabled=.FALSE.)
241   IF (TRIM(etat0_type)=='start_file' .AND. read_metric_) THEN
242     CALL xios_set_file_attr('start', enabled=.TRUE.)
243   ENDIF
244     
245
246   CALL xios_close_context_definition()
247!$OMP END MASTER
248!$OMP BARRIER
249   
250 END SUBROUTINE xios_init_write_field_input 
251
252
253
254 SUBROUTINE xios_write_field(name,field)
255   CHARACTER(LEN=*),INTENT(IN) :: name
256   TYPE(t_field), POINTER :: field(:)
257   TYPE(t_cellset), POINTER :: cells(:)
258   INTEGER :: ncells
259
260!$OMP BARRIER
261!$OMP MASTER
262
263   SELECT CASE(field(1)%field_type)
264   CASE(field_T)
265      cells => mesh_loc%primal_own
266      ncells = ncell_i
267   CASE(field_U)
268      cells => mesh_loc%edge_own
269      ncells = ncell_e
270   CASE(field_Z)
271      cells => mesh_loc%dual_own
272      ncells = ncell_v
273   END SELECT
274
275   IF (field(1)%ndim>4) THEN
276      PRINT *, "xios_write_field : dimension > 4 are not supported for now"
277   ELSE
278      CALL xios_write_field_gen(name, field, cells, & 
279           ncells, field(1)%dim3, field(1)%dim4)
280   END IF
281   
282!$OMP END MASTER
283!$OMP BARRIER
284     
285 END SUBROUTINE xios_write_field
286
287
288
289 SUBROUTINE xios_read_field(name,field)
290   CHARACTER(LEN=*),INTENT(IN) :: name
291   TYPE(t_field), POINTER :: field(:)
292   TYPE(t_cellset), POINTER :: cells(:)
293   INTEGER :: ncells
294
295!$OMP BARRIER
296!$OMP MASTER
297
298   SELECT CASE(field(1)%field_type)
299   CASE(field_T)
300      cells => mesh_loc%primal_own
301      ncells = ncell_i
302   CASE(field_U)
303      cells => mesh_loc%edge_own
304      ncells = ncell_e
305   CASE(field_Z)
306      cells => mesh_loc%dual_own
307      ncells = ncell_v
308   END SELECT
309
310   IF (field(1)%ndim>4) THEN
311      PRINT *, "xios_read_field : dimension > 4 are not supported for now"
312   ELSE
313      CALL xios_read_field_hex(name, field, cells, & 
314           ncells, field(1)%dim3, field(1)%dim4)
315   END IF
316   
317!$OMP END MASTER
318!$OMP BARRIER
319     
320 END SUBROUTINE xios_read_field
321
322 SUBROUTINE xios_write_field_gen(name, field, cells, ncell_tot, nlev, nq)
323   CHARACTER(LEN=*),INTENT(IN) :: name
324   TYPE(t_field) :: field(:)
325   TYPE(t_cellset), TARGET :: cells(:)
326   INTEGER,INTENT(IN) :: ncell_tot, nlev, nq
327   REAL(rstd) :: field_tmp(ncell_tot, nlev, nq)
328   TYPE(t_cellset), POINTER :: cellset
329   INTEGER :: ind, n_beg, n, ij, sgn
330   LOGICAL :: signed
331
332   IF(ALLOCATED(cells(1)%sgn)) THEN
333      signed=.TRUE.
334   ELSE
335      signed=.FALSE.
336      sgn=1
337   END IF
338
339   n_beg=0
340   DO ind=1,ndomain
341      cellset => cells(ind)
342      DO n=1, cellset%ncell
343         ij = cellset%ij(n)
344         IF(signed) sgn = cellset%sgn(n)
345         SELECT CASE(grid_type)
346         CASE(grid_ico)
347            SELECT CASE(field(1)%ndim)
348            CASE(2)
349               field_tmp(n_beg+n,1,1) = sgn*field(ind)%rval2d(ij)
350            CASE(3)
351               field_tmp(n_beg+n,:,1) = sgn*field(ind)%rval3d(ij,:)
352            CASE(4)
353               field_tmp(n_beg+n,:,:) = sgn*field(ind)%rval4d(ij,:,:)
354            END SELECT
355         CASE(grid_unst)
356            SELECT CASE(field(1)%ndim)
357            CASE(2)
358               field_tmp(n_beg+n,1,1) = sgn*field(ind)%rval2d(ij)
359            CASE(3)
360               field_tmp(n_beg+n,:,1) = sgn*field(ind)%rval3d(:,ij)
361            CASE(4)
362               field_tmp(n_beg+n,:,:) = sgn*field(ind)%rval4d(:,ij,:)
363            END SELECT
364         END SELECT
365      END DO
366      n_beg = n_beg + cellset%ncell
367   END DO
368   CALL xios_send_field(name,field_tmp) 
369   
370 END SUBROUTINE xios_write_field_gen
371
372SUBROUTINE xios_read_field_hex(name, field, cells, ncell_tot, nlev, nq)
373   CHARACTER(LEN=*),INTENT(IN) :: name
374   TYPE(t_field), POINTER :: field(:)
375   TYPE(t_cellset), TARGET :: cells(:)
376   INTEGER,INTENT(IN) :: ncell_tot, nlev, nq
377
378   REAL(rstd) :: field_tmp(ncell_tot,nlev,nq)
379   TYPE(t_cellset), POINTER :: cellset
380   INTEGER :: ind, n_beg, n, ij, sgn
381   LOGICAL :: signed
382
383   CALL xios_recv_field(name,field_tmp)
384
385   IF(ALLOCATED(cells(1)%sgn)) THEN
386      signed=.TRUE.
387   ELSE
388      signed=.FALSE.
389      sgn=1
390   END IF
391
392
393   n_beg=0
394   DO ind=1,ndomain
395      cellset => cells(ind)
396      DO n=1, cellset%ncell
397         ij = cellset%ij(n)
398         IF(signed) sgn = cellset%sgn(n)
399         SELECT CASE(field(1)%ndim)
400         CASE(2)
401            field(ind)%rval2d(ij) = sgn*field_tmp(n_beg+n,1,1) 
402         CASE(3)
403            field(ind)%rval3d(ij,:) = sgn*field_tmp(n_beg+n,:,1)
404         CASE(4)
405            field(ind)%rval4d(ij,:,:) = sgn*field_tmp(n_beg+n,:,:)
406         END SELECT
407      END DO
408      n_beg = n_beg + cellset%ncell
409   END DO
410 END SUBROUTINE xios_read_field_hex
411
412 SUBROUTINE xios_read_var(name,field)
413   USE prec
414   USE transfert_mod
415   CHARACTER(LEN=*),INTENT(IN) :: name
416   REAL(rstd), INTENT(OUT) :: field
417   !$OMP MASTER
418   CALL xios_recv_field(name,field)
419   !$OMP END MASTER
420   CALL bcast_omp(field)
421 END SUBROUTINE
422 
423 SUBROUTINE xios_write_field_finalize
424!$OMP BARRIER
425!$OMP MASTER
426   CALL xios_context_finalize
427!$OMP END MASTER
428!$OMP BARRIER
429 END SUBROUTINE xios_write_field_finalize
430
431 SUBROUTINE xios_set_context
432  TYPE(xios_context) :: ctx_hdl
433
434!$OMP MASTER
435   CALL xios_get_handle("icosagcm",ctx_hdl)
436   CALL xios_set_current_context(ctx_hdl)
437!$OMP END MASTER
438
439  END SUBROUTINE xios_set_context
440
441 SUBROUTINE xios_set_context_input
442 IMPLICIT NONE   
443  TYPE(xios_context) :: ctx_hdl
444
445!$OMP MASTER
446   CALL xios_get_handle("icosagcm_input",ctx_hdl)
447   CALL xios_set_current_context(ctx_hdl)
448!$OMP END MASTER
449
450  END SUBROUTINE xios_set_context_input
451
452#else
453 
454
455INTERFACE xios_send_field
456  MODULE PROCEDURE xios_send_field_scalar, xios_send_field_1d
457END INTERFACE  xios_send_field
458
459INTEGER,PARAMETER :: xios_timestep=1
460
461CONTAINS 
462 
463 
464  SUBROUTINE xios_init
465    using_xios=.FALSE.   
466  END SUBROUTINE xios_init
467 
468  SUBROUTINE xios_send_field_scalar(name,field)
469    CHARACTER(LEN=*),INTENT(IN) :: name
470    REAL,INTENT(IN) :: field
471  END SUBROUTINE xios_send_field_scalar 
472
473  SUBROUTINE xios_send_field_1d(name,field)
474    CHARACTER(LEN=*),INTENT(IN) :: name
475    REAL,INTENT(IN) :: field(:)
476  END SUBROUTINE xios_send_field_1d 
477 
478  SUBROUTINE xios_write_field(name,field)
479  USE field_mod
480   CHARACTER(LEN=*),INTENT(IN) :: name
481   TYPE(t_field), POINTER :: field(:)
482  END SUBROUTINE xios_write_field
483
484  SUBROUTINE xios_read_field(name,field)
485  USE field_mod
486   CHARACTER(LEN=*),INTENT(IN) :: name
487   TYPE(t_field), POINTER :: field(:)
488  END SUBROUTINE xios_read_field
489
490 SUBROUTINE xios_read_var(name,field)
491   USE prec
492   CHARACTER(LEN=*),INTENT(IN) :: name
493   REAL(rstd), INTENT(OUT) :: field
494 END SUBROUTINE
495 
496  SUBROUTINE xios_update_calendar(step)
497   INTEGER, INTENT(IN):: step 
498  END SUBROUTINE xios_update_calendar
499
500  SUBROUTINE xios_write_field_finalize
501  END SUBROUTINE xios_write_field_finalize
502 
503  SUBROUTINE xios_init_write_field_input
504  END SUBROUTINE
505
506  SUBROUTINE xios_init_write_field
507  END SUBROUTINE xios_init_write_field 
508 
509  SUBROUTINE xios_set_context
510  END SUBROUTINE xios_set_context
511 
512  SUBROUTINE xios_set_fieldgroup_attr(name,enabled,freq_op)
513    CHARACTER(LEN=*) :: name
514    LOGICAL,OPTIONAL          :: enabled
515    INTEGER,OPTIONAL          :: freq_op
516  END SUBROUTINE xios_set_fieldgroup_attr
517
518  SUBROUTINE xios_set_filegroup_attr(name,enabled)
519    CHARACTER(LEN=*) :: name
520    LOGICAL,OPTIONAL          :: enabled
521  END SUBROUTINE xios_set_filegroup_attr
522
523  SUBROUTINE xios_set_file_attr(id,name,mode,enabled, output_freq)
524    CHARACTER(LEN=*) :: id
525    CHARACTER(LEN=*),OPTIONAL :: name, mode
526    LOGICAL,OPTIONAL          :: enabled
527    INTEGER,OPTIONAL          :: output_freq
528  END SUBROUTINE xios_set_file_attr
529
530  SUBROUTINE xios_get_axis_attr(name,n_glo,value)
531    CHARACTER(LEN=*) :: name
532    INTEGER,OPTIONAL          :: n_glo
533    REAL,OPTIONAL             :: value(:)
534  END SUBROUTINE xios_get_axis_attr
535
536  SUBROUTINE xios_set_axis_attr(id,n_glo,value)
537    CHARACTER(LEN=*) :: id
538    INTEGER,OPTIONAL          :: n_glo
539    REAL,OPTIONAL             :: value(:)
540  END SUBROUTINE xios_set_axis_attr
541
542#endif 
543
544END MODULE xios_mod
Note: See TracBrowser for help on using the repository browser.