Context Navigation

limit_read_mod.F90 @ 3496

Last change on this file since 3496 was 3435, checked in by Laurent Fairhead, 6 years ago

"Historic" :-) commit merging the physics branch used for DYNAMICO with the LMDZ trunk.
The same physics branch can now be used seamlessly with the traditional lon-lat LMDZ
dynamical core and DYNAMICO.
Testing consisted in running a lon-lat LMDZ bucket simulation with the NPv6.1 physics package
with the original trunk sources and the merged sources. Tests were succesful in the sense that
numeric continuity was preserved in the restart files from both simulation. Further tests
included running both versions of the physics codes for one year in a LMDZOR setting in which
the restart files also came out identical.

Caution:

as the physics package now manages unstructured grids, grid information needs to be transmitted

to the surface scheme ORCHIDEE. This means that the interface defined in surf_land_orchidee_mod.F90
is only compatible with ORCHIDEE version orchidee2.1 and later versions. If previous versions of
ORCHIDEE need to be used, the CPP key ORCHIDEE_NOUNSTRUCT needs to be set at compilation time.
This is done automatically if makelmdz/makelmdz_fcm are called with the veget orchidee2.0 switch

due to a limitation in XIOS, the time at which limit conditions will be read in by DYNAMICO will be

delayed by one physic timestep with respect to the time it is read in by the lon-lat model. This is caused
by the line

IF (MOD(itime-1, lmt_pas) == 0 .OR. (jour_lu /= jour .AND. grid_type /= unstructured)) THEN ! time to read

in limit_read_mod.F90

Work still needed on COSP integration and XML files for DYNAMICO

EM, YM, LF

Property copyright set to
Name of program: LMDZ
Creation date: 1984
Version: LMDZ5
License: CeCILL version 2
Holder: Laboratoire de m\'et\'eorologie dynamique, CNRS, UMR 8539
See the license file in the root directory
Property svn:eol-style set to native
Property svn:keywords set to Author Date Id Revision

File size: 17.9 KB

Line
1	!
2	! $Id: limit_read_mod.F90 3435 2019-01-22 15:21:59Z jyg $
3	!
4	MODULE limit_read_mod
5	!
6	! This module reads the fichier "limit.nc" containing fields for surface forcing.
7	!
8	! Module subroutines :
9	! limit_read_frac : call limit_read_tot and return the fractions
10	! limit_read_rug_alb : return rugosity and albedo, if coupled ocean call limit_read_tot first
11	! limit_read_sst : return sea ice temperature
12	! limit_read_tot : read limit.nc and store the fields in local modules variables
13	!
14	IMPLICIT NONE
15
16	REAL, ALLOCATABLE, DIMENSION(:,:), SAVE, PRIVATE :: pctsrf
17	!$OMP THREADPRIVATE(pctsrf)
18	REAL, ALLOCATABLE, DIMENSION(:), SAVE, PRIVATE :: rugos
19	!$OMP THREADPRIVATE(rugos)
20	REAL, ALLOCATABLE, DIMENSION(:), SAVE, PRIVATE :: albedo
21	!$OMP THREADPRIVATE(albedo)
22	REAL, ALLOCATABLE, DIMENSION(:), SAVE, PRIVATE :: sst
23	!$OMP THREADPRIVATE(sst)
24	LOGICAL,SAVE :: read_continents=.FALSE.
25	!$OMP THREADPRIVATE(read_continents)
26
27	CONTAINS
28	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
29	!!
30	!! Public subroutines :
31	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
32
33
34	SUBROUTINE init_limit_read(first_day)
35	USE mod_grid_phy_lmdz
36	USE surface_data
37	USE mod_phys_lmdz_para
38	#ifdef CPP_XIOS
39	USE XIOS
40	#endif
41	IMPLICIT NONE
42	INTEGER, INTENT(IN) :: first_day
43
44
45	IF ( type_ocean /= 'couple') THEN
46	IF (grid_type==unstructured) THEN
47	#ifdef CPP_XIOS
48	IF (is_omp_master) CALL xios_set_file_attr("limit_read",enabled=.TRUE.,record_offset=first_day)
49	#endif
50	ENDIF
51	ENDIF
52
53	END SUBROUTINE init_limit_read
54
55	SUBROUTINE limit_read_frac(itime, dtime, jour, pctsrf_new, is_modified)
56	!
57	! This subroutine is called from "change_srf_frac" for case of
58	! ocean=force or from ocean_slab_frac for ocean=slab.
59	! The fraction for all sub-surfaces at actual time step is returned.
60
61	USE dimphy
62	USE indice_sol_mod
63
64	! Input arguments
65	!****************************************************************************************
66	INTEGER, INTENT(IN) :: itime ! time step
67	INTEGER, INTENT(IN) :: jour ! current day
68	REAL , INTENT(IN) :: dtime ! length of time step
69
70	! Output arguments
71	!****************************************************************************************
72	REAL, DIMENSION(klon,nbsrf), INTENT(OUT) :: pctsrf_new ! sub surface fractions
73	LOGICAL, INTENT(OUT) :: is_modified ! true if pctsrf is modified at this time step
74
75	! End declaration
76	!****************************************************************************************
77
78	! 1) Read file limit.nc
79	CALL limit_read_tot(itime, dtime, jour, is_modified)
80
81	! 2) Return the fraction read in limit_read_tot
82	pctsrf_new(:,:) = pctsrf(:,:)
83
84	END SUBROUTINE limit_read_frac
85
86	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
87
88	SUBROUTINE limit_read_rug_alb(itime, dtime, jour, &
89	knon, knindex, &
90	rugos_out, alb_out)
91	!
92	! This subroutine is called from surf_land_bucket.
93	! The flag "ok_veget" must can not be true. If coupled run, "ocean=couple"
94	! then this routine will call limit_read_tot.
95	!
96	USE dimphy
97	USE surface_data
98
99	! Input arguments
100	!****************************************************************************************
101	INTEGER, INTENT(IN) :: itime ! numero du pas de temps courant
102	INTEGER, INTENT(IN) :: jour ! jour a lire dans l'annee
103	REAL , INTENT(IN) :: dtime ! pas de temps de la physique (en s)
104	INTEGER, INTENT(IN) :: knon ! nomber of points on compressed grid
105	INTEGER, DIMENSION(klon), INTENT(IN) :: knindex ! grid point number for compressed grid
106	! Output arguments
107	!****************************************************************************************
108	REAL, DIMENSION(klon), INTENT(OUT) :: rugos_out
109	REAL, DIMENSION(klon), INTENT(OUT) :: alb_out
110
111	! Local variables
112	!****************************************************************************************
113	INTEGER :: i
114	LOGICAL :: is_modified
115	!****************************************************************************************
116
117	IF (type_ocean == 'couple'.OR. &
118	(type_ocean == 'slab' .AND. version_ocean == 'sicINT')) THEN
119	! limit.nc has not yet been read. Do it now!
120	CALL limit_read_tot(itime, dtime, jour, is_modified)
121	END IF
122
123	DO i=1,knon
124	rugos_out(i) = rugos(knindex(i))
125	alb_out(i) = albedo(knindex(i))
126	END DO
127
128	END SUBROUTINE limit_read_rug_alb
129
130	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
131
132	SUBROUTINE limit_read_sst(knon, knindex, sst_out)
133	!
134	! This subroutine returns the sea surface temperature already read from limit.nc.
135	!
136	USE dimphy, ONLY : klon
137
138	INTEGER, INTENT(IN) :: knon ! nomber of points on compressed grid
139	INTEGER, DIMENSION(klon), INTENT(IN) :: knindex ! grid point number for compressed grid
140	REAL, DIMENSION(klon), INTENT(OUT) :: sst_out
141
142	INTEGER :: i
143
144	DO i = 1, knon
145	sst_out(i) = sst(knindex(i))
146	END DO
147
148	END SUBROUTINE limit_read_sst
149
150	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
151	!!
152	!! Private subroutine :
153	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
154
155	SUBROUTINE limit_read_tot(itime, dtime, jour, is_modified)
156	!
157	! Read everything needed from limit.nc
158	!
159	! 0) Initialize
160	! 1) Open the file limit.nc, if it is time
161	! 2) Read fraction, if not type_ocean=couple
162	! 3) Read sea surface temperature, if not type_ocean=couple
163	! 4) Read albedo and rugosity for land surface, only in case of no vegetation model
164	! 5) Close file and distribuate variables to all processus
165
166	USE dimphy
167	USE mod_grid_phy_lmdz
168	USE mod_phys_lmdz_para
169	USE surface_data, ONLY : type_ocean, ok_veget
170	USE netcdf
171	USE indice_sol_mod
172	USE phys_cal_mod, ONLY : calend, year_len
173	USE print_control_mod, ONLY: lunout, prt_level
174	#ifdef CPP_XIOS
175	USE XIOS, ONLY: xios_recv_field
176	#endif
177
178	IMPLICIT NONE
179
180	! In- and ouput arguments
181	!****************************************************************************************
182	INTEGER, INTENT(IN) :: itime ! numero du pas de temps courant
183	INTEGER, INTENT(IN) :: jour ! jour a lire dans l'annee
184	REAL , INTENT(IN) :: dtime ! pas de temps de la physique (en s)
185
186	LOGICAL, INTENT(OUT) :: is_modified ! true if pctsrf is modified at this time step
187
188	! Locals variables with attribute SAVE
189	!****************************************************************************************
190	! frequence de lecture des conditions limites (en pas de physique)
191	INTEGER,SAVE :: lmt_pas
192	!$OMP THREADPRIVATE(lmt_pas)
193	LOGICAL, SAVE :: first_call=.TRUE.
194	!$OMP THREADPRIVATE(first_call)
195	INTEGER, SAVE :: jour_lu = -1
196	!$OMP THREADPRIVATE(jour_lu)
197	! Locals variables
198	!****************************************************************************************
199	INTEGER :: nid, nvarid, ndimid, nn
200	INTEGER :: ii, ierr
201	INTEGER, DIMENSION(2) :: start, epais
202	REAL, DIMENSION(klon_glo,nbsrf) :: pct_glo ! fraction at global grid
203	REAL, DIMENSION(klon_glo) :: sst_glo ! sea-surface temperature at global grid
204	REAL, DIMENSION(klon_glo) :: rug_glo ! rugosity at global grid
205	REAL, DIMENSION(klon_glo) :: alb_glo ! albedo at global grid
206
207	REAL, DIMENSION(klon_mpi,nbsrf) :: pct_mpi ! fraction at global grid
208	REAL, DIMENSION(klon_mpi) :: sst_mpi ! sea-surface temperature at global grid
209	REAL, DIMENSION(klon_mpi) :: rug_mpi ! rugosity at global grid
210	REAL, DIMENSION(klon_mpi) :: alb_mpi ! albedo at global grid
211
212	CHARACTER(len=20) :: modname='limit_read_mod'
213	CHARACTER(LEN=99) :: abort_message, calendar, str
214
215	! End declaration
216	!****************************************************************************************
217
218	!****************************************************************************************
219	! 0) Initialization
220	!
221	!****************************************************************************************
222	IF (first_call) THEN
223	first_call=.FALSE.
224	! calculate number of time steps for one day
225	lmt_pas = NINT(86400./dtime * 1.0)
226
227	! Allocate module save variables
228	IF ( type_ocean /= 'couple' ) THEN
229	ALLOCATE(pctsrf(klon,nbsrf), sst(klon), stat=ierr)
230	IF (ierr /= 0) CALL abort_physic(modname, 'PB in allocating pctsrf and sst',1)
231	END IF
232
233	IF ( .NOT. ok_veget ) THEN
234	ALLOCATE(rugos(klon), albedo(klon), stat=ierr)
235	IF (ierr /= 0) CALL abort_physic(modname, 'PB in allocating rugos and albedo',1)
236	END IF
237
238	!$OMP MASTER ! Only master thread
239	IF (is_mpi_root) THEN ! Only master processus
240	ierr = NF90_OPEN ('limit.nc', NF90_NOWRITE, nid)
241	IF (ierr /= NF90_NOERR) CALL abort_physic(modname,&
242	'Pb d''ouverture du fichier de conditions aux limites',1)
243
244	!--- WARNING IF CALENDAR IS KNOWN AND DOES NOT MATCH THE ONE OF LMDZ
245	ierr=NF90_INQ_VARID(nid, 'TEMPS', nvarid)
246	ierr=NF90_GET_ATT(nid, nvarid, 'calendar', calendar)
247	IF(ierr==NF90_NOERR.AND.calendar/=calend.AND.prt_level>=1) THEN
248	WRITE(lunout,*)'BEWARE: gcm and limit.nc calendars differ: '
249	WRITE(lunout,*)' '//TRIM(calend)//' for gcm'
250	WRITE(lunout,*)' '//TRIM(calendar)//' for limit.nc file'
251	END IF
252
253	!--- ERROR IF FILE RECORDS NUMBER IS NOT EQUAL TO EXPECTED NUMBER OF DAYS
254	IF (grid_type==unstructured) THEN
255	ierr=NF90_INQ_DIMID(nid,"time_year",ndimid)
256	ELSE
257	ierr=NF90_INQUIRE(nid, UnlimitedDimID=ndimid)
258	ENDIF
259	ierr=NF90_INQUIRE_DIMENSION(nid, ndimid, len=nn)
260	WRITE(abort_message,'(a,2(i0,a))')'limit.nc records number (',nn,') does no'//&
261	't match year length (',year_len,')'
262	IF(nn/=year_len) CALL abort_physic(modname,abort_message,1)
263
264	!--- ERROR IF FILES AND LMDZ HORIZONTAL RESOLUTIONS DO NOT MATCH
265	IF (grid_type==unstructured) THEN
266	ierr=NF90_INQ_DIMID(nid, 'cell', ndimid)
267	ELSE
268	ierr=NF90_INQ_DIMID(nid, 'points_physiques', ndimid)
269	ENDIF
270	ierr=NF90_INQUIRE_DIMENSION(nid, ndimid, len=nn)
271	WRITE(abort_message,'(a,2(i0,a))')'limit.nc horizontal number of cells (',nn, &
272	') does not match LMDZ klon_glo (',klon_glo,')'
273	IF(nn/=klon_glo) CALL abort_physic(modname,abort_message,1)
274
275	ierr = NF90_CLOSE(nid)
276	IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Pb when closing file', 1)
277	END IF ! is_mpi_root
278	!$OMP END MASTER
279	!$OMP BARRIER
280	END IF
281
282	!****************************************************************************************
283	! 1) Open the file limit.nc if it is the right moment to read, once a day.
284	! The file is read only by the master thread of the master mpi process(is_mpi_root)
285	! Check by the way if the number of records is correct.
286	!
287	!****************************************************************************************
288
289	is_modified = .FALSE.
290	!ym IF (MOD(itime-1, lmt_pas) == 0 .OR. jour_lu /= jour ) THEN ! time to read
291	! not REALLY PERIODIC
292	IF (MOD(itime-1, lmt_pas) == 0 .OR. (jour_lu /= jour .AND. grid_type /= unstructured)) THEN ! time to read
293	! IF (MOD(itime-1, lmt_pas) == 0) THEN ! time to read
294	jour_lu = jour
295	is_modified = .TRUE.
296
297	IF (grid_type==unstructured) THEN
298
299	#ifdef CPP_XIOS
300	IF ( type_ocean /= 'couple') THEN
301
302	IF (is_omp_master) CALL xios_recv_field("foce_limin",pct_mpi(:,is_oce))
303	IF (is_omp_master) CALL xios_recv_field("fsic_limin",pct_mpi(:,is_sic))
304	! IF (read_continents .OR. itime == 1) THEN
305	IF (is_omp_master) CALL xios_recv_field("fter_limin",pct_mpi(:,is_ter))
306	IF (is_omp_master) CALL xios_recv_field("flic_limin",pct_mpi(:,is_lic))
307	! ENDIF
308	ENDIF! type_ocean /= couple
309
310	IF ( type_ocean /= 'couple') THEN
311	IF (is_omp_master) CALL xios_recv_field("sst_limin",sst_mpi)
312	ENDIF
313
314	IF (.NOT. ok_veget) THEN
315	IF (is_omp_master) CALL xios_recv_field("alb_limin",alb_mpi)
316	IF (is_omp_master) CALL xios_recv_field("rug_limin",rug_mpi)
317	ENDIF
318
319	IF ( type_ocean /= 'couple') THEN
320	CALL Scatter_omp(sst_mpi,sst)
321	CALL Scatter_omp(pct_mpi(:,is_oce),pctsrf(:,is_oce))
322	CALL Scatter_omp(pct_mpi(:,is_sic),pctsrf(:,is_sic))
323	! IF (read_continents .OR. itime == 1) THEN
324	CALL Scatter_omp(pct_mpi(:,is_ter),pctsrf(:,is_ter))
325	CALL Scatter_omp(pct_mpi(:,is_lic),pctsrf(:,is_lic))
326	! END IF
327	END IF
328
329	IF (.NOT. ok_veget) THEN
330	CALL Scatter_omp(alb_mpi, albedo)
331	CALL Scatter_omp(rug_mpi, rugos)
332	END IF
333	#endif
334
335
336	ELSE ! grid_type==regular
337
338	!$OMP MASTER ! Only master thread
339	IF (is_mpi_root) THEN ! Only master processus!
340
341	ierr = NF90_OPEN ('limit.nc', NF90_NOWRITE, nid)
342	IF (ierr /= NF90_NOERR) CALL abort_physic(modname,&
343	'Pb d''ouverture du fichier de conditions aux limites',1)
344
345	! La tranche de donnees a lire:
346	start(1) = 1
347	start(2) = jour
348	epais(1) = klon_glo
349	epais(2) = 1
350
351
352	!****************************************************************************************
353	! 2) Read fraction if not type_ocean=couple
354	!
355	!****************************************************************************************
356
357	IF ( type_ocean /= 'couple') THEN
358	!
359	! Ocean fraction
360	ierr = NF90_INQ_VARID(nid, 'FOCE', nvarid)
361	IF (ierr /= NF90_NOERR) CALL abort_physic(modname, 'Le champ <FOCE> est absent',1)
362
363	ierr = NF90_GET_VAR(nid,nvarid,pct_glo(:,is_oce),start,epais)
364	IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Lecture echouee pour <FOCE>' ,1)
365	!
366	! Sea-ice fraction
367	ierr = NF90_INQ_VARID(nid, 'FSIC', nvarid)
368	IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Le champ <FSIC> est absent',1)
369
370	ierr = NF90_GET_VAR(nid,nvarid,pct_glo(:,is_sic),start,epais)
371	IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Lecture echouee pour <FSIC>' ,1)
372
373
374	! Read land and continentals fraction only if asked for
375	IF (read_continents .OR. itime == 1) THEN
376	!
377	! Land fraction
378	ierr = NF90_INQ_VARID(nid, 'FTER', nvarid)
379	IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Le champ <FTER> est absent',1)
380
381	ierr = NF90_GET_VAR(nid,nvarid,pct_glo(:,is_ter),start,epais)
382	IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Lecture echouee pour <FTER>',1)
383	!
384	! Continentale ice fraction
385	ierr = NF90_INQ_VARID(nid, 'FLIC', nvarid)
386	IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Le champ <FLIC> est absent',1)
387
388	ierr = NF90_GET_VAR(nid,nvarid,pct_glo(:,is_lic),start,epais)
389	IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Lecture echouee pour <FLIC>',1)
390	END IF
391
392	END IF ! type_ocean /= couple
393
394	!****************************************************************************************
395	! 3) Read sea-surface temperature, if not coupled ocean
396	!
397	!****************************************************************************************
398	IF ( type_ocean /= 'couple') THEN
399
400	ierr = NF90_INQ_VARID(nid, 'SST', nvarid)
401	IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Le champ <SST> est absent',1)
402
403	ierr = NF90_GET_VAR(nid,nvarid,sst_glo,start,epais)
404	IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Lecture echouee pour <SST>',1)
405
406	END IF
407
408	!****************************************************************************************
409	! 4) Read albedo and rugosity for land surface, only in case of no vegetation model
410	!
411	!****************************************************************************************
412
413	IF (.NOT. ok_veget) THEN
414	!
415	! Read albedo
416	ierr = NF90_INQ_VARID(nid, 'ALB', nvarid)
417	IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Le champ <ALB> est absent',1)
418
419	ierr = NF90_GET_VAR(nid,nvarid,alb_glo,start,epais)
420	IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Lecture echouee pour <ALB>',1)
421	!
422	! Read rugosity
423	ierr = NF90_INQ_VARID(nid, 'RUG', nvarid)
424	IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Le champ <RUG> est absent',1)
425
426	ierr = NF90_GET_VAR(nid,nvarid,rug_glo,start,epais)
427	IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Lecture echouee pour <RUG>',1)
428
429	END IF
430
431	!****************************************************************************************
432	! 5) Close file and distribuate variables to all processus
433	!
434	!****************************************************************************************
435	ierr = NF90_CLOSE(nid)
436	IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Pb when closing file', 1)
437	ENDIF ! is_mpi_root
438
439	!$OMP END MASTER
440	!$OMP BARRIER
441
442	IF ( type_ocean /= 'couple') THEN
443	CALL Scatter(sst_glo,sst)
444	CALL Scatter(pct_glo(:,is_oce),pctsrf(:,is_oce))
445	CALL Scatter(pct_glo(:,is_sic),pctsrf(:,is_sic))
446	IF (read_continents .OR. itime == 1) THEN
447	CALL Scatter(pct_glo(:,is_ter),pctsrf(:,is_ter))
448	CALL Scatter(pct_glo(:,is_lic),pctsrf(:,is_lic))
449	END IF
450	END IF
451
452	IF (.NOT. ok_veget) THEN
453	CALL Scatter(alb_glo, albedo)
454	CALL Scatter(rug_glo, rugos)
455	END IF
456
457	ENDIF ! Grid type
458
459	ENDIF ! time to read
460
461	END SUBROUTINE limit_read_tot
462
463	END MODULE limit_read_mod

Note: See TracBrowser for help on using the repository browser.

Download in other formats:

Original Format