source: XIOS/branchs/xios-2.0/src/test/toy_cmip6.f90 @ 1627

Last change on this file since 1627 was 1118, checked in by mhnguyen, 7 years ago

Adding toy_cmip6 into repository
File size: 12.0 KB
Line 
1PROGRAM toy_cmip6 
2
3  USE xios
4  USE mod_wait
5  IMPLICIT NONE
6  INCLUDE "mpif.h"
7
8  INTEGER,PARAMETER :: il_unit=10
9  INTEGER :: comm, rank, size_loc, ierr
10  INTEGER :: ni,ibegin,iend,nj,jbegin,jend
11  INTEGER :: i,j,l,ts,n, nb_pt, il_run
12  DOUBLE PRECISION :: sypd, timestep_in_seconds, simulated_seconds_per_seconds, elapsed_per_timestep
13  CHARACTER(len=*),PARAMETER :: id="client"
14  CHARACTER(1000):: duration, timestep
15  INTEGER :: start_year,start_month,start_day
16  TYPE(xios_date) :: cdate, edate
17  TYPE(xios_duration)  :: dtime
18  TYPE(xios_context) :: ctx_hdl
19  REAL :: ilon,jlat
20  DOUBLE PRECISION,ALLOCATABLE :: lon_glo(:,:),lat_glo(:,:),lval(:)
21  DOUBLE PRECISION,ALLOCATABLE :: field_A_glo (:,:,:), pressure_glo (:,:,:), height_glo (:,:,:)
22  DOUBLE PRECISION,ALLOCATABLE :: bounds_lon_glo(:,:,:),bounds_lat_glo(:,:,:)
23  DOUBLE PRECISION,ALLOCATABLE :: pressure (:,:,:), height (:,:,:)
24  DOUBLE PRECISION,ALLOCATABLE :: lon(:,:),lat(:,:),lonvalue(:,:)
25  DOUBLE PRECISION,ALLOCATABLE :: bounds_lon(:,:,:),bounds_lat(:,:,:) 
26  DOUBLE PRECISION,ALLOCATABLE :: field_atm_2D(:,:),field_atm_3D(:,:,:),field_srf_2D(:),field_srf_3D(:,:)
27  DOUBLE PRECISION,ALLOCATABLE :: field_atm_2D_miss(:,:)
28  DOUBLE PRECISION,ALLOCATABLE :: field_oce_2D(:,:),field_oce_3D(:,:,:)
29  INTEGER, ALLOCATABLE :: kindex(:)
30
31  INTEGER :: ni_glo, nj_glo,llm
32
33  NAMELIST /param_toy/ ni_glo, nj_glo,llm,timestep,duration,sypd,start_year,start_month,start_day
34
35!!! MPI Initialization
36
37  CALL MPI_INIT(ierr)
38
39  CALL init_wait
40
41!!! Lecture des parametres du run
42
43  OPEN(unit=il_unit, file='param.def',status='old',iostat=ierr)
44  READ (il_unit, nml=param_toy)
45  !PRINT *, ni_glo, nj_glo,llm,duration
46
47!!! XIOS Initialization (get the local communicator)
48
49  CALL xios_initialize(id,return_comm=comm)
50
51  CALL MPI_COMM_RANK(comm,rank,ierr)
52  CALL MPI_COMM_SIZE(comm,size_loc,ierr)
53
54
55!!! Initialisation et allocation des coordonnées globales et locales pour la grille réguliÚre
56
57  ALLOCATE (lon_glo(ni_glo,nj_glo),lat_glo(ni_glo,nj_glo))
58  ALLOCATE(bounds_lon_glo(4,ni_glo,nj_glo))
59  ALLOCATE(bounds_lat_glo(4,ni_glo,nj_glo))
60  ALLOCATE (field_A_glo(ni_glo,nj_glo,llm))
61  ALLOCATE (pressure_glo(ni_glo,nj_glo,llm))
62  ALLOCATE (height_glo(ni_glo,nj_glo,llm))
63  ALLOCATE (lval(llm))
64
65  DO j=1,nj_glo
66    DO i=1,ni_glo
67
68      ilon=i-0.5
69      jlat=j-0.5
70
71      lat_glo(i,j)= 90-(jlat*180./nj_glo)
72      lon_glo(i,j)= (ilon*360./ni_glo)
73      !print*, 'i/lon=',i,'lon=',lon_glo(i,j), 'j/lat=',j,'lat=',lat_glo(i,j)
74
75      bounds_lat_glo(1,i,j)= 90-((jlat-0.5)*180./nj_glo)
76      bounds_lon_glo(1,i,j)=((ilon-0.5)*360./ni_glo)
77
78      bounds_lat_glo(2,i,j)= 90-((jlat-0.5)*180./nj_glo)
79      bounds_lon_glo(2,i,j)=((ilon+0.5)*360./ni_glo)
80
81      bounds_lat_glo(3,i,j)= 90-((jlat+0.5)*180./nj_glo)
82      bounds_lon_glo(3,i,j)=((ilon+0.5)*360./ni_glo)
83
84      bounds_lat_glo(4,i,j)= 90-((jlat+0.5)*180./nj_glo)
85      bounds_lon_glo(4,i,j)=((ilon-0.5)*360./ni_glo)
86
87      WHERE (abs(bounds_lat_glo(:,i,j) - 90) < 0.000000001) bounds_lat_glo(:,i,j) = 90
88      WHERE (abs(bounds_lat_glo(:,i,j) + 90) < 0.000000001) bounds_lat_glo(:,i,j) = -90
89
90      DO l=1,llm
91         field_A_glo(i,j,l)=(i-1)+(j-1)*ni_glo+10000*l
92         ! pressure at half levels. First index value is high altitude, low pressure
93         pressure_glo(i,j,l)=((l-0.)/llm)*100000 + (jlat -nj_glo/2.)/nj_glo * 10000
94         height_glo(i,j,l)=(llm-l+0.5)/llm * 15000 + jlat * 100
95      ENDDO
96    ENDDO
97  ENDDO
98  ni=ni_glo ; ibegin=0
99
100  jbegin=0
101  DO n=0,size_loc-1
102    nj=nj_glo/size_loc
103    IF (n<MOD(nj_glo,size_loc)) nj=nj+1
104    IF (n==rank) exit
105    jbegin=jbegin+nj
106  ENDDO
107
108  iend=ibegin+ni-1 ; jend=jbegin+nj-1
109
110  ALLOCATE(lon(ni,nj),lat(ni,nj),lonvalue(ni,nj))
111  ALLOCATE(bounds_lon(4,ni,nj))
112  ALLOCATE(bounds_lat(4,ni,nj))
113  lon(:,:)=lon_glo(ibegin+1:iend+1,jbegin+1:jend+1)
114  lat(:,:)=lat_glo(ibegin+1:iend+1,jbegin+1:jend+1)
115  bounds_lon(:,:,:)=bounds_lon_glo(:,ibegin+1:iend+1,jbegin+1:jend+1)
116  bounds_lat(:,:,:)=bounds_lat_glo(:,ibegin+1:iend+1,jbegin+1:jend+1)
117
118
119  DO i=1,llm
120    lval(i)=i
121 ENDDO
122 
123
124 
125!###########################################################################
126! Contexte ATM
127!###########################################################################
128 ALLOCATE(field_atm_2D(0:ni+1,-1:nj+2),field_atm_3D(0:ni+1,-1:nj+2,llm))
129 ALLOCATE(field_atm_2D_miss(0:ni+1,-1:nj+2))
130 ALLOCATE(pressure(0:ni+1,-1:nj+2,llm))
131 ALLOCATE(height(0:ni+1,-1:nj+2,llm))
132  field_atm_2D(1:ni,1:nj)=field_A_glo(ibegin+1:iend+1,jbegin+1:jend+1,1)
133  field_atm_2D_miss(1:ni,1:nj)=field_A_glo(ibegin+1:iend+1,jbegin+1:jend+1,1)
134  field_atm_3D(1:ni,1:nj,:)=field_A_glo(ibegin+1:iend+1,jbegin+1:jend+1,:)
135  pressure(1:ni,1:nj,:)=pressure_glo(ibegin+1:iend+1,jbegin+1:jend+1,:)
136  height(1:ni,1:nj,:)=height_glo(ibegin+1:iend+1,jbegin+1:jend+1,:)
137
138  CALL xios_context_initialize("arpsfx",comm)
139  CALL xios_define_calendar("Gregorian", &
140     start_date=xios_date(start_year,start_month,start_day,0,0,0), &
141     time_origin=xios_date(1850,1,1,0,0,0))
142
143  !write(0,*) 'atm context initialized' ; call flush(0)
144  CALL xios_get_handle("arpsfx",ctx_hdl)
145  CALL xios_set_current_context(ctx_hdl)
146
147  CALL xios_set_axis_attr("axis_atm",n_glo=llm ,value=lval) ;
148
149  CALL xios_set_domain_attr("domain_atm",ni_glo=ni_glo, nj_glo=nj_glo, ibegin=ibegin, &
150       ni=ni,jbegin=jbegin,nj=nj, type='rectilinear')
151  CALL xios_set_domain_attr("domain_atm",data_dim=2, data_ibegin=-1, &
152       data_ni=ni+2, data_jbegin=-2, data_nj=nj+4)
153  CALL xios_set_domain_attr("domain_atm",lonvalue_2D=lon,latvalue_2D=lat)
154  CALL xios_set_domain_attr("domain_atm", nvertex=4, bounds_lon_2d=bounds_lon, bounds_lat_2d=bounds_lat)
155  print *,'latmax/min=',minval(lat),maxval(lat)
156  print *,'latmax/min bounds=',minval(bounds_lat),maxval(bounds_lat)
157  print *,'lonmax/min=',minval(lon),maxval(lon)
158  print *,'lonmax/min bounds=',minval(bounds_lon),maxval(bounds_lon)
159
160
161!!! Definition du timestep
162
163  CALL xios_get_start_date(cdate)
164  edate=cdate+xios_duration_convert_from_string(duration)
165  dtime=xios_duration_convert_from_string(timestep)
166  CALL xios_set_timestep(timestep=dtime)
167
168!!! Fin de la definition du contexte
169
170  !CALL xios_close_context_definition()
171
172!!! Calcul de temps elaps par seconde pour respecter le SYPD (hyp : pas de délai d'I/O)
173 
174  timestep_in_seconds=xios_date_convert_to_seconds(cdate+dtime) - xios_date_convert_to_seconds(cdate)
175  simulated_seconds_per_seconds=sypd * 365 
176  elapsed_per_timestep=timestep_in_seconds/simulated_seconds_per_seconds ! in seconds
177
178!###########################################################################
179! Contexte SRF
180!###########################################################################
181
182!!! Initialisation des coordonnées globales et locales pour la grille indexee (1 point sur 2)
183
184    nb_pt=ni*nj/2
185    ALLOCATE(kindex(nb_pt),field_srf_2D(nb_pt),field_srf_3D(nb_pt,llm))
186    DO i=1,nb_pt
187      kindex(i)=2*i-1
188    ENDDO
189    field_srf_2D(1:nb_pt)=RESHAPE(field_A_glo(ibegin+1:iend+1:2,jbegin+1:jend+1,1),(/ nb_pt /))
190    field_srf_3D(1:nb_pt,:)=RESHAPE(field_A_glo(ibegin+1:iend+1:2,jbegin+1:jend+1,:),(/ nb_pt,llm /))
191
192  !CALL xios_context_initialize("surface",comm)
193  !CALL xios_get_handle("surface",ctx_hdl)
194  !CALL xios_set_current_context(ctx_hdl)
195
196  CALL xios_set_axis_attr("axis_srf",n_glo=llm ,value=lval)
197  CALL xios_set_domain_attr("domain_srf",ni_glo=ni_glo, nj_glo=nj_glo, ibegin=ibegin, ni=ni,jbegin=jbegin,nj=nj, type='rectilinear')
198  CALL xios_set_domain_attr("domain_srf",data_dim=1, data_ibegin=0, data_ni=nb_pt)
199  CALL xios_set_domain_attr("domain_srf",data_i_index=kindex)
200  CALL xios_set_domain_attr("domain_srf",lonvalue_2D=lon,latvalue_2D=lat)
201  CALL xios_set_domain_attr("domain_srf", nvertex=4, bounds_lon_2d=bounds_lon, bounds_lat_2d=bounds_lat)
202
203!!! Definition du timestep
204
205  !dtime%second=timestep
206  !dtime=xios_duration_convert_from_string(timestep)
207  !CALL xios_set_timestep(timestep=dtime)
208
209!!! Fin de la definition du contexte SRF
210
211  CALL xios_close_context_definition()
212  !write(0,*) 'srf context def closed' ; call flush(0)
213
214!###########################################################################
215! Contexte OCE
216!###########################################################################
217  ALLOCATE(field_oce_2D(0:ni+1,-1:nj+2),field_oce_3D(0:ni+1,-1:nj+2,llm))
218  field_oce_2D(1:ni,1:nj)=field_A_glo(ibegin+1:iend+1,jbegin+1:jend+1,1)
219  field_oce_3D(1:ni,1:nj,:)=field_A_glo(ibegin+1:iend+1,jbegin+1:jend+1,:)
220
221  CALL xios_context_initialize("nemo",comm)
222  CALL xios_define_calendar("Gregorian", &
223     start_date=xios_date(start_year,start_month,start_day,0,0,0), &
224     time_origin=xios_date(1850,1,1,0,0,0))
225  CALL xios_get_handle("nemo",ctx_hdl)
226  CALL xios_set_current_context(ctx_hdl)
227
228  CALL xios_set_axis_attr("axis_oce",n_glo=llm ,value=lval) ;
229
230  CALL xios_set_domain_attr("domain_oce",ni_glo=ni_glo, nj_glo=nj_glo, ibegin=ibegin, ni=ni,jbegin=jbegin,nj=nj, type='curvilinear')
231  CALL xios_set_domain_attr("domain_oce",data_dim=2, data_ibegin=-1, data_ni=ni+2, data_jbegin=-2, data_nj=nj+4)
232  CALL xios_set_domain_attr("domain_oce",lonvalue_2D=lon,latvalue_2D=lat)
233
234
235!!! Definition du timestep
236
237  dtime=xios_duration_convert_from_string(timestep)
238  CALL xios_set_timestep(timestep=dtime)
239
240  CALL xios_close_context_definition()
241
242!####################################################################################
243!!! Boucle temporelle
244!####################################################################################
245  ts=1
246  cdate=cdate+dtime
247  DO while ( cdate <= edate )
248
249      CALL xios_get_handle("arpsfx",ctx_hdl)
250      CALL xios_set_current_context(ctx_hdl)
251
252!!! Mise a jour du pas de temps
253
254      CALL xios_update_calendar(ts)
255
256!!! On donne la valeur du champ atm
257
258      !print *,'sending field_atm_2d at timestep',ts
259      CALL xios_send_field("field_atm_scalar",field_atm_2D(1,1)+ts)
260      CALL xios_send_field("field_atm_1D",field_atm_3D(1,1,:)+ts)
261      CALL xios_send_field("field_atm_2D",field_atm_2D+ts)
262      CALL xios_send_field("field_atm_3D",field_atm_3D+ts)
263      CALL xios_send_field("pressure" ,pressure)
264      CALL xios_send_field("height"   ,height)
265      if (mod(ts,2)==0) then
266         CALL xios_send_field("field_sub",field_atm_2D+ts)
267      endif
268      !! On crée un champ avec des missings qui bougent
269      !! dans le temps : un bande verticale de 1 à ni-3
270      field_atm_2D_miss(:,:)= field_atm_2D(:,:)+ts
271      field_atm_2D_miss(mod(ts,ni-3)+1,:)=1.e+20
272      CALL xios_send_field("field_miss",field_atm_2D_miss)
273     
274!!! On change de contexte
275
276      !CALL xios_get_handle("surface",ctx_hdl)
277      !CALL xios_set_current_context(ctx_hdl)
278
279!!! Mise a jour du pas de temps
280
281      !CALL xios_update_calendar(ts)
282
283!!! On donne la valeur du champ srf
284
285      CALL xios_send_field("field_srf_2D",field_srf_2D)
286      CALL xios_send_field("field_srf_3D",field_srf_3D)
287
288!!! On change de contexte
289
290      CALL xios_get_handle("nemo",ctx_hdl)
291      CALL xios_set_current_context(ctx_hdl)
292
293!!! Mise a jour du pas de temps
294
295      CALL xios_update_calendar(ts)
296
297!!! On donne la valeur du champ oce
298
299      CALL xios_send_field("field_oce_scalar",field_oce_2D(1,1)+ts)
300      CALL xios_send_field("field_oce_grid_2D",field_oce_2D)
301      CALL xios_send_field("field_oce_grid_3D",field_oce_3D)
302
303      CALL wait_us(int(elapsed_per_timestep*1.e6))   ! micro-secondes
304      cdate=cdate+dtime
305      ts=ts+1
306
307   ENDDO
308
309!####################################################################################
310!!! Finalisation
311!####################################################################################
312
313!!! Fin des contextes
314
315    CALL xios_context_finalize()
316    CALL xios_get_handle("arpsfx",ctx_hdl)
317    CALL xios_set_current_context(ctx_hdl)
318    CALL xios_context_finalize()
319    !CALL xios_get_handle("surface",ctx_hdl)
320    !CALL xios_set_current_context(ctx_hdl)
321    !CALL xios_context_finalize()
322
323    DEALLOCATE(lon, lat, lonvalue, field_atm_2D, field_atm_3D)
324    DEALLOCATE(pressure,height,pressure_glo,height_glo)
325    DEALLOCATE(kindex, field_srf_2D, field_srf_3D)
326    DEALLOCATE(field_oce_2D, field_oce_3D)
327    DEALLOCATE(lon_glo,lat_glo,field_A_glo,lval)
328
329!!! Fin de XIOS
330
331    CALL MPI_COMM_FREE(comm, ierr)
332     
333    CALL xios_finalize()
334
335  CALL MPI_FINALIZE(ierr)
336
337END PROGRAM toy_cmip6 
338
339
340
341
342
343
Note: See TracBrowser for help on using the repository browser.