source: dynamico_lmdz/aquaplanet/LMDZ5/libf/phylmd/limit_read_mod.F90 @ 3867

Last change on this file since 3867 was 3867, checked in by ymipsl, 9 years ago
  • Switch to XIOS 2
  • Append main part of forcing configuration (LMDZ stand alone)
  • Create etat0 and limit from usual LMDZ input files, using XIOS 2 interpolation functionnalities
    • missing parametrization of gravity waves
    • missing aerosol
File size: 15.6 KB
Line 
1!
2! $Header$
3!
4MODULE 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
27CONTAINS
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 XIOS
38  USE etat0_limit_unstruct_mod, ONLY: create_etat0_limit
39  IMPLICIT NONE
40    INTEGER, INTENT(IN) :: first_day
41   
42   
43    IF ( type_ocean /= 'couple') THEN
44      IF (grid_type==unstructured) THEN
45        IF (.NOT. create_etat0_limit) CALL xios_set_file_attr("limit_read",enabled=.TRUE.,record_offset=first_day)
46      ENDIF 
47    ENDIF
48 
49  END SUBROUTINE init_limit_read
50 
51  SUBROUTINE limit_read_frac(itime, dtime, jour, pctsrf_new, is_modified)
52!
53! This subroutine is called from "change_srf_frac" for case of
54! ocean=force or from ocean_slab_frac for ocean=slab.
55! The fraction for all sub-surfaces at actual time step is returned.
56
57    USE dimphy
58    USE indice_sol_mod
59
60! Input arguments
61!****************************************************************************************
62    INTEGER, INTENT(IN) :: itime   ! time step
63    INTEGER, INTENT(IN) :: jour    ! current day
64    REAL   , INTENT(IN) :: dtime   ! length of time step
65 
66! Output arguments
67!****************************************************************************************
68    REAL, DIMENSION(klon,nbsrf), INTENT(OUT) :: pctsrf_new  ! sub surface fractions
69    LOGICAL, INTENT(OUT)                     :: is_modified ! true if pctsrf is modified at this time step
70
71! End declaration
72!****************************************************************************************
73
74! 1) Read file limit.nc
75    CALL limit_read_tot(itime, dtime, jour, is_modified)
76
77! 2) Return the fraction read in limit_read_tot
78    pctsrf_new(:,:) = pctsrf(:,:)
79   
80  END SUBROUTINE limit_read_frac
81
82!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
83
84  SUBROUTINE limit_read_rug_alb(itime, dtime, jour, &
85       knon, knindex, &
86       rugos_out, alb_out)
87!
88! This subroutine is called from surf_land_bucket.
89! The flag "ok_veget" must can not be true. If coupled run, "ocean=couple"
90! then this routine will call limit_read_tot.
91!
92    USE dimphy
93    USE surface_data
94
95! Input arguments
96!****************************************************************************************
97    INTEGER, INTENT(IN) :: itime                     ! numero du pas de temps courant
98    INTEGER, INTENT(IN) :: jour                      ! jour a lire dans l'annee
99    REAL   , INTENT(IN) :: dtime                     ! pas de temps de la physique (en s)
100    INTEGER, INTENT(IN) :: knon                      ! nomber of points on compressed grid
101    INTEGER, DIMENSION(klon), INTENT(IN) :: knindex  ! grid point number for compressed grid
102! Output arguments
103!****************************************************************************************
104    REAL, DIMENSION(klon), INTENT(OUT) :: rugos_out
105    REAL, DIMENSION(klon), INTENT(OUT) :: alb_out
106   
107! Local variables
108!****************************************************************************************
109    INTEGER :: i
110    LOGICAL :: is_modified
111!****************************************************************************************
112
113IF (type_ocean == 'couple'.OR. &
114         (type_ocean == 'slab' .AND. version_ocean == 'sicINT')) THEN
115       ! limit.nc has not yet been read. Do it now!
116       CALL limit_read_tot(itime, dtime, jour, is_modified)
117    END IF
118
119    DO i=1,knon
120       rugos_out(i) = rugos(knindex(i))
121       alb_out(i)  = albedo(knindex(i))
122    END DO
123
124  END SUBROUTINE limit_read_rug_alb
125
126!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
127
128  SUBROUTINE limit_read_sst(knon, knindex, sst_out)
129!
130! This subroutine returns the sea surface temperature already read from limit.nc.
131!
132    USE dimphy, ONLY : klon
133
134    INTEGER, INTENT(IN)                  :: knon     ! nomber of points on compressed grid
135    INTEGER, DIMENSION(klon), INTENT(IN) :: knindex  ! grid point number for compressed grid
136    REAL, DIMENSION(klon), INTENT(OUT)   :: sst_out
137
138    INTEGER :: i
139
140    DO i = 1, knon
141       sst_out(i) = sst(knindex(i))
142    END DO
143
144  END SUBROUTINE limit_read_sst
145
146!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
147!!
148!! Private subroutine :
149!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
150
151  SUBROUTINE limit_read_tot(itime, dtime, jour, is_modified)
152!
153! Read everything needed from limit.nc
154!
155! 0) Initialize
156! 1) Open the file limit.nc, if it is time
157! 2) Read fraction, if not type_ocean=couple
158! 3) Read sea surface temperature, if not type_ocean=couple
159! 4) Read albedo and rugosity for land surface, only in case of no vegetation model
160! 5) Close file and distribuate variables to all processus
161
162    USE dimphy
163    USE mod_grid_phy_lmdz
164    USE mod_phys_lmdz_para
165    USE surface_data, ONLY : type_ocean, ok_veget
166    USE netcdf
167    USE indice_sol_mod
168    USE XIOS
169   
170    IMPLICIT NONE
171   
172! In- and ouput arguments
173!****************************************************************************************
174    INTEGER, INTENT(IN) :: itime   ! numero du pas de temps courant
175    INTEGER, INTENT(IN) :: jour    ! jour a lire dans l'annee
176    REAL   , INTENT(IN) :: dtime   ! pas de temps de la physique (en s)
177
178    LOGICAL, INTENT(OUT) :: is_modified  ! true if pctsrf is modified at this time step
179
180! Locals variables with attribute SAVE
181!****************************************************************************************
182! frequence de lecture des conditions limites (en pas de physique)
183    INTEGER,SAVE                              :: lmt_pas
184!$OMP THREADPRIVATE(lmt_pas)
185    LOGICAL, SAVE                             :: first_call=.TRUE.
186!$OMP THREADPRIVATE(first_call) 
187    INTEGER, SAVE                             :: jour_lu = -1
188!$OMP THREADPRIVATE(jour_lu) 
189! Locals variables
190!****************************************************************************************
191    INTEGER                                   :: nid, nvarid
192    INTEGER                                   :: ii, ierr
193    INTEGER, DIMENSION(2)                     :: start, epais
194    REAL, DIMENSION(klon_glo,nbsrf)           :: pct_glo  ! fraction at global grid
195    REAL, DIMENSION(klon_glo)                 :: sst_glo  ! sea-surface temperature at global grid
196    REAL, DIMENSION(klon_glo)                 :: rug_glo  ! rugosity at global grid
197    REAL, DIMENSION(klon_glo)                 :: alb_glo  ! albedo at global grid
198
199    REAL, DIMENSION(klon_mpi,nbsrf)           :: pct_mpi  ! fraction at global grid
200    REAL, DIMENSION(klon_mpi)                 :: sst_mpi  ! sea-surface temperature at global grid
201    REAL, DIMENSION(klon_mpi)                 :: rug_mpi  ! rugosity at global grid
202    REAL, DIMENSION(klon_mpi)                 :: alb_mpi  ! albedo at global grid
203
204    CHARACTER(len=20)                         :: modname='limit_read_mod'     
205
206! End declaration
207!****************************************************************************************
208
209!****************************************************************************************
210! 0) Initialization
211!
212!****************************************************************************************
213    IF (first_call) THEN
214       ! calculate number of time steps for one day
215       lmt_pas = NINT(86400./dtime * 1.0)
216       
217       ! Allocate module save variables
218       IF ( type_ocean /= 'couple' ) THEN
219          ALLOCATE(pctsrf(klon,nbsrf), sst(klon), stat=ierr)
220          IF (ierr /= 0) CALL abort_physic(modname, 'PB in allocating pctsrf and sst',1)
221       END IF
222
223       IF ( .NOT. ok_veget ) THEN
224          ALLOCATE(rugos(klon), albedo(klon), stat=ierr)
225          IF (ierr /= 0) CALL abort_physic(modname, 'PB in allocating rugos and albedo',1)
226       END IF
227
228       first_call=.FALSE.
229    ENDIF
230 
231!****************************************************************************************
232! 1) Open the file limit.nc if it is the right moment to read, once a day.
233!    The file is read only by the master thread of the master mpi process(is_mpi_root)
234!
235!****************************************************************************************
236
237    is_modified = .FALSE.
238!ym    IF (MOD(itime-1, lmt_pas) == 0 .OR. jour_lu /= jour ) THEN   ! time to read
239!  not REALLY PERIODIC
240    IF (MOD(itime-1, lmt_pas) == 0) THEN   ! time to read
241       jour_lu = jour
242       is_modified = .TRUE.
243
244      IF (grid_type==unstructured) THEN
245
246!$OMP MASTER  ! Only master thread
247
248
249        IF ( type_ocean /= 'couple') THEN
250
251           CALL xios_recv_field("foce_limin",pct_mpi(:,is_oce))
252           CALL xios_recv_field("fsic_limin",pct_mpi(:,is_sic))
253  !         IF (read_continents .OR. itime == 1) THEN
254             CALL xios_recv_field("fter_limin",pct_mpi(:,is_ter))
255             CALL xios_recv_field("flic_limin",pct_mpi(:,is_lic))
256  !         ENDIF
257         ENDIF! type_ocean /= couple
258         
259         IF ( type_ocean /= 'couple') THEN                   
260             CALL xios_recv_field("sst_limin",sst_mpi)
261         ENDIF
262       
263         IF (.NOT. ok_veget) THEN
264           CALL xios_recv_field("alb_limin",alb_mpi)
265           CALL xios_recv_field("rug_limin",rug_mpi)
266         ENDIF
267
268       IF ( type_ocean /= 'couple') THEN
269          CALL Scatter_omp(sst_mpi,sst)
270          CALL Scatter_omp(pct_mpi(:,is_oce),pctsrf(:,is_oce))
271          CALL Scatter_omp(pct_mpi(:,is_sic),pctsrf(:,is_sic))
272!          IF (read_continents .OR. itime == 1) THEN
273             CALL Scatter_omp(pct_mpi(:,is_ter),pctsrf(:,is_ter))
274             CALL Scatter_omp(pct_mpi(:,is_lic),pctsrf(:,is_lic))
275!          END IF
276       END IF
277
278       IF (.NOT. ok_veget) THEN
279          CALL Scatter_omp(alb_mpi, albedo)
280          CALL Scatter_omp(rug_mpi, rugos)
281       END IF
282
283!$OMP END MASTER
284
285 
286     ELSE      ! grid_type==regular
287
288!$OMP MASTER  ! Only master thread
289
290       IF (is_mpi_root) THEN ! Only master processus!
291
292            ierr = NF90_OPEN ('limit.nc', NF90_NOWRITE, nid)
293            IF (ierr /= NF90_NOERR) CALL abort_physic(modname,&
294                 'Pb d''ouverture du fichier de conditions aux limites',1)
295           
296            ! La tranche de donnees a lire:
297            start(1) = 1
298            start(2) = jour
299            epais(1) = klon_glo
300            epais(2) = 1
301
302
303  !****************************************************************************************
304  ! 2) Read fraction if not type_ocean=couple
305  !
306  !****************************************************************************************
307
308            IF ( type_ocean /= 'couple') THEN
309  !
310  ! Ocean fraction
311               ierr = NF90_INQ_VARID(nid, 'FOCE', nvarid)
312               IF (ierr /= NF90_NOERR) CALL abort_physic(modname, 'Le champ <FOCE> est absent',1)
313               
314               ierr = NF90_GET_VAR(nid,nvarid,pct_glo(:,is_oce),start,epais)
315               IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Lecture echouee pour <FOCE>' ,1)
316  !
317  ! Sea-ice fraction
318               ierr = NF90_INQ_VARID(nid, 'FSIC', nvarid)
319               IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Le champ <FSIC> est absent',1)
320
321               ierr = NF90_GET_VAR(nid,nvarid,pct_glo(:,is_sic),start,epais)
322               IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Lecture echouee pour <FSIC>' ,1)
323
324
325  ! Read land and continentals fraction only if asked for
326               IF (read_continents .OR. itime == 1) THEN
327  !
328  ! Land fraction
329                  ierr = NF90_INQ_VARID(nid, 'FTER', nvarid)
330                  IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Le champ <FTER> est absent',1)
331                 
332                  ierr = NF90_GET_VAR(nid,nvarid,pct_glo(:,is_ter),start,epais)
333                  IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Lecture echouee pour <FTER>',1)
334  !
335  ! Continentale ice fraction
336                  ierr = NF90_INQ_VARID(nid, 'FLIC', nvarid)
337                  IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Le champ <FLIC> est absent',1)
338
339                  ierr = NF90_GET_VAR(nid,nvarid,pct_glo(:,is_lic),start,epais)
340                  IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Lecture echouee pour <FLIC>',1)
341               END IF
342
343
344          END IF ! type_ocean /= couple
345
346!****************************************************************************************
347! 3) Read sea-surface temperature, if not coupled ocean
348!
349!****************************************************************************************
350          IF ( type_ocean /= 'couple') THEN
351
352             ierr = NF90_INQ_VARID(nid, 'SST', nvarid)
353             IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Le champ <SST> est absent',1)
354
355             ierr = NF90_GET_VAR(nid,nvarid,sst_glo,start,epais)
356             IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Lecture echouee pour <SST>',1)
357         
358          END IF
359
360!****************************************************************************************
361! 4) Read albedo and rugosity for land surface, only in case of no vegetation model
362!
363!****************************************************************************************
364
365             IF (.NOT. ok_veget) THEN
366   !
367   ! Read albedo
368                ierr = NF90_INQ_VARID(nid, 'ALB', nvarid)
369                IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Le champ <ALB> est absent',1)
370
371                ierr = NF90_GET_VAR(nid,nvarid,alb_glo,start,epais)
372                IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Lecture echouee pour <ALB>',1)
373   !
374   ! Read rugosity
375                ierr = NF90_INQ_VARID(nid, 'RUG', nvarid)
376                IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Le champ <RUG> est absent',1)
377
378                ierr = NF90_GET_VAR(nid,nvarid,rug_glo,start,epais)
379                IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Lecture echouee pour <RUG>',1)
380
381             END IF
382
383!****************************************************************************************
384! 5) Close file and distribuate variables to all processus
385!
386!****************************************************************************************
387          ierr = NF90_CLOSE(nid)
388          IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Pb when closing file', 1)
389       ENDIF ! is_mpi_root
390   
391
392  !$OMP END MASTER
393  !$OMP BARRIER
394
395       IF ( type_ocean /= 'couple') THEN
396          CALL Scatter(sst_glo,sst)
397          CALL Scatter(pct_glo(:,is_oce),pctsrf(:,is_oce))
398          CALL Scatter(pct_glo(:,is_sic),pctsrf(:,is_sic))
399          IF (read_continents .OR. itime == 1) THEN
400             CALL Scatter(pct_glo(:,is_ter),pctsrf(:,is_ter))
401             CALL Scatter(pct_glo(:,is_lic),pctsrf(:,is_lic))
402          END IF
403       END IF
404
405       IF (.NOT. ok_veget) THEN
406          CALL Scatter(alb_glo, albedo)
407          CALL Scatter(rug_glo, rugos)
408       END IF
409     
410     ENDIF    ! grid_type
411
412    ENDIF ! time to read
413
414  END SUBROUTINE limit_read_tot
415
416
417END MODULE limit_read_mod
Note: See TracBrowser for help on using the repository browser.