source: LMDZ5/trunk/libf/phylmd/limit_read_mod.F90 @ 2679

Last change on this file since 2679 was 2667, checked in by fhourdin, 8 years ago

Suite de la correction precedente
  • 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: 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  SUBROUTINE limit_read_frac(itime, dtime, jour, pctsrf_new, is_modified)
34!
35! This subroutine is called from "change_srf_frac" for case of
36! ocean=force or from ocean_slab_frac for ocean=slab.
37! The fraction for all sub-surfaces at actual time step is returned.
38
39    USE dimphy
40    USE indice_sol_mod
41
42! Input arguments
43!****************************************************************************************
44    INTEGER, INTENT(IN) :: itime   ! time step
45    INTEGER, INTENT(IN) :: jour    ! current day
46    REAL   , INTENT(IN) :: dtime   ! length of time step
47 
48! Output arguments
49!****************************************************************************************
50    REAL, DIMENSION(klon,nbsrf), INTENT(OUT) :: pctsrf_new  ! sub surface fractions
51    LOGICAL, INTENT(OUT)                     :: is_modified ! true if pctsrf is modified at this time step
52
53! End declaration
54!****************************************************************************************
55
56! 1) Read file limit.nc
57    CALL limit_read_tot(itime, dtime, jour, is_modified)
58
59! 2) Return the fraction read in limit_read_tot
60    pctsrf_new(:,:) = pctsrf(:,:)
61   
62  END SUBROUTINE limit_read_frac
63
64!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
65
66  SUBROUTINE limit_read_rug_alb(itime, dtime, jour, &
67       knon, knindex, &
68       rugos_out, alb_out)
69!
70! This subroutine is called from surf_land_bucket.
71! The flag "ok_veget" must can not be true. If coupled run, "ocean=couple"
72! then this routine will call limit_read_tot.
73!
74    USE dimphy
75    USE surface_data
76
77! Input arguments
78!****************************************************************************************
79    INTEGER, INTENT(IN) :: itime                     ! numero du pas de temps courant
80    INTEGER, INTENT(IN) :: jour                      ! jour a lire dans l'annee
81    REAL   , INTENT(IN) :: dtime                     ! pas de temps de la physique (en s)
82    INTEGER, INTENT(IN) :: knon                      ! nomber of points on compressed grid
83    INTEGER, DIMENSION(klon), INTENT(IN) :: knindex  ! grid point number for compressed grid
84! Output arguments
85!****************************************************************************************
86    REAL, DIMENSION(klon), INTENT(OUT) :: rugos_out
87    REAL, DIMENSION(klon), INTENT(OUT) :: alb_out
88   
89! Local variables
90!****************************************************************************************
91    INTEGER :: i
92    LOGICAL :: is_modified
93!****************************************************************************************
94
95IF (type_ocean == 'couple'.OR. &
96         (type_ocean == 'slab' .AND. version_ocean == 'sicINT')) THEN
97       ! limit.nc has not yet been read. Do it now!
98       CALL limit_read_tot(itime, dtime, jour, is_modified)
99    END IF
100
101    DO i=1,knon
102       rugos_out(i) = rugos(knindex(i))
103       alb_out(i)  = albedo(knindex(i))
104    END DO
105
106  END SUBROUTINE limit_read_rug_alb
107
108!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
109
110  SUBROUTINE limit_read_sst(knon, knindex, sst_out)
111!
112! This subroutine returns the sea surface temperature already read from limit.nc.
113!
114    USE dimphy, ONLY : klon
115
116    INTEGER, INTENT(IN)                  :: knon     ! nomber of points on compressed grid
117    INTEGER, DIMENSION(klon), INTENT(IN) :: knindex  ! grid point number for compressed grid
118    REAL, DIMENSION(klon), INTENT(OUT)   :: sst_out
119
120    INTEGER :: i
121
122    DO i = 1, knon
123       sst_out(i) = sst(knindex(i))
124    END DO
125
126  END SUBROUTINE limit_read_sst
127
128!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
129!!
130!! Private subroutine :
131!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
132
133  SUBROUTINE limit_read_tot(itime, dtime, jour, is_modified)
134!
135! Read everything needed from limit.nc
136!
137! 0) Initialize
138! 1) Open the file limit.nc, if it is time
139! 2) Read fraction, if not type_ocean=couple
140! 3) Read sea surface temperature, if not type_ocean=couple
141! 4) Read albedo and rugosity for land surface, only in case of no vegetation model
142! 5) Close file and distribuate variables to all processus
143
144    USE dimphy
145    USE mod_grid_phy_lmdz
146    USE mod_phys_lmdz_para
147    USE surface_data, ONLY : type_ocean, ok_veget
148    USE netcdf
149    USE indice_sol_mod
150    USE phys_cal_mod, ONLY : calend, year_len
151    USE print_control_mod, ONLY: lunout, prt_level
152
153    IMPLICIT NONE
154   
155! In- and ouput arguments
156!****************************************************************************************
157    INTEGER, INTENT(IN) :: itime   ! numero du pas de temps courant
158    INTEGER, INTENT(IN) :: jour    ! jour a lire dans l'annee
159    REAL   , INTENT(IN) :: dtime   ! pas de temps de la physique (en s)
160
161    LOGICAL, INTENT(OUT) :: is_modified  ! true if pctsrf is modified at this time step
162
163! Locals variables with attribute SAVE
164!****************************************************************************************
165! frequence de lecture des conditions limites (en pas de physique)
166    INTEGER,SAVE                              :: lmt_pas
167!$OMP THREADPRIVATE(lmt_pas)
168    LOGICAL, SAVE                             :: first_call=.TRUE.
169!$OMP THREADPRIVATE(first_call) 
170    INTEGER, SAVE                             :: jour_lu = -1
171!$OMP THREADPRIVATE(jour_lu) 
172! Locals variables
173!****************************************************************************************
174    INTEGER                                   :: nid, nvarid, ndimid, nn
175    INTEGER                                   :: ii, ierr
176    INTEGER, DIMENSION(2)                     :: start, epais
177    REAL, DIMENSION(klon_glo,nbsrf)           :: pct_glo  ! fraction at global grid
178    REAL, DIMENSION(klon_glo)                 :: sst_glo  ! sea-surface temperature at global grid
179    REAL, DIMENSION(klon_glo)                 :: rug_glo  ! rugosity at global grid
180    REAL, DIMENSION(klon_glo)                 :: alb_glo  ! albedo at global grid
181    CHARACTER(len=20)                         :: modname='limit_read_mod'     
182    CHARACTER(LEN=99)                         :: abort_message, calendar, str
183
184! End declaration
185!****************************************************************************************
186
187!****************************************************************************************
188! 0) Initialization
189!
190!****************************************************************************************
191    IF (first_call) THEN
192       ! calculate number of time steps for one day
193       lmt_pas = NINT(86400./dtime * 1.0)
194       
195       ! Allocate module save variables
196       IF ( type_ocean /= 'couple' ) THEN
197          ALLOCATE(pctsrf(klon,nbsrf), sst(klon), stat=ierr)
198          IF (ierr /= 0) CALL abort_physic(modname, 'PB in allocating pctsrf and sst',1)
199       END IF
200
201       IF ( .NOT. ok_veget ) THEN
202          ALLOCATE(rugos(klon), albedo(klon), stat=ierr)
203          IF (ierr /= 0) CALL abort_physic(modname, 'PB in allocating rugos and albedo',1)
204       END IF
205
206       first_call=.FALSE.
207    ENDIF
208 
209!****************************************************************************************
210! 1) Open the file limit.nc if it is the right moment to read, once a day.
211!    The file is read only by the master thread of the master mpi process(is_mpi_root)
212!    Check by the way if the number of records is correct.
213!
214!****************************************************************************************
215
216    is_modified = .FALSE.
217    IF (MOD(itime-1, lmt_pas) == 0 .OR. jour_lu /= jour ) THEN   ! time to read
218       jour_lu = jour
219       is_modified = .TRUE.
220!$OMP MASTER  ! Only master thread
221       IF (is_mpi_root) THEN ! Only master processus
222
223          ierr = NF90_OPEN ('limit.nc', NF90_NOWRITE, nid)
224          IF (ierr /= NF90_NOERR) CALL abort_physic(modname,&
225               'Pb d''ouverture du fichier de conditions aux limites',1)
226
227          !--- WARNING IF CALENDAR IS KNOWN AND DOES NOT MATCH THE ONE OF LMDZ
228          ierr=NF90_INQ_VARID(nid, 'TEMPS', nvarid)
229          ierr=NF90_GET_ATT(nid, nvarid, 'calendar', calendar)
230          IF(ierr==NF90_NOERR.AND.calendar/=calend.AND.prt_level>=1) THEN
231             WRITE(lunout,*)'BEWARE: gcm and limit.nc calendars differ: '
232             WRITE(lunout,*)'  '//TRIM(calend)//' for gcm'
233             WRITE(lunout,*)'  '//TRIM(calendar)//' for limit.nc file'
234          END IF
235
236          !--- ERROR IF FILE RECORDS NUMBER IS NOT EQUAL TO EXPECTED NUMBER OF DAYS
237          ierr=NF90_INQUIRE(nid, UnlimitedDimID=ndimid)
238          ierr=NF90_INQUIRE_DIMENSION(nid, ndimid, len=nn)
239          CALL num2str(nn,str)
240          abort_message='limit.nc records number ('//TRIM(str)//') does'//&
241            ' not match year length ('
242          CALL num2str(year_len,str)
243          abort_message=TRIM(abort_message)//TRIM(str)//')'
244          IF(nn/=year_len) CALL abort_physic(modname,abort_message,1)
245
246          !--- ERROR IF FILES AND LMDZ HORIZONTAL RESOLUTIONS DO NOT MATCH
247          ierr=NF90_INQ_DIMID(nid, 'points_physiques', ndimid)
248          ierr=NF90_INQUIRE_DIMENSION(nid, ndimid, len=nn)
249          CALL num2str(nn,str)
250          abort_message='limit.nc horizontal number of cells ('//TRIM(str)//') does'//&
251            ' not match LMDZ klon_glo ('
252          CALL num2str(klon_glo,str)
253          abort_message=TRIM(abort_message)//TRIM(str)//')'
254          IF(nn/=klon_glo) CALL abort_physic(modname,abort_message,1)
255
256          ! La tranche de donnees a lire:
257          start(1) = 1
258          start(2) = jour
259          epais(1) = klon_glo
260          epais(2) = 1
261
262
263!****************************************************************************************
264! 2) Read fraction if not type_ocean=couple
265!
266!****************************************************************************************
267
268          IF ( type_ocean /= 'couple') THEN
269!
270! Ocean fraction
271             ierr = NF90_INQ_VARID(nid, 'FOCE', nvarid)
272             IF (ierr /= NF90_NOERR) CALL abort_physic(modname, 'Le champ <FOCE> est absent',1)
273             
274             ierr = NF90_GET_VAR(nid,nvarid,pct_glo(:,is_oce),start,epais)
275             IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Lecture echouee pour <FOCE>' ,1)
276!
277! Sea-ice fraction
278             ierr = NF90_INQ_VARID(nid, 'FSIC', nvarid)
279             IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Le champ <FSIC> est absent',1)
280
281             ierr = NF90_GET_VAR(nid,nvarid,pct_glo(:,is_sic),start,epais)
282             IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Lecture echouee pour <FSIC>' ,1)
283
284
285! Read land and continentals fraction only if asked for
286             IF (read_continents .OR. itime == 1) THEN
287!
288! Land fraction
289                ierr = NF90_INQ_VARID(nid, 'FTER', nvarid)
290                IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Le champ <FTER> est absent',1)
291               
292                ierr = NF90_GET_VAR(nid,nvarid,pct_glo(:,is_ter),start,epais)
293                IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Lecture echouee pour <FTER>',1)
294!
295! Continentale ice fraction
296                ierr = NF90_INQ_VARID(nid, 'FLIC', nvarid)
297                IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Le champ <FLIC> est absent',1)
298
299                ierr = NF90_GET_VAR(nid,nvarid,pct_glo(:,is_lic),start,epais)
300                IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Lecture echouee pour <FLIC>',1)
301             END IF
302
303          END IF ! type_ocean /= couple
304
305!****************************************************************************************
306! 3) Read sea-surface temperature, if not coupled ocean
307!
308!****************************************************************************************
309          IF ( type_ocean /= 'couple') THEN
310
311             ierr = NF90_INQ_VARID(nid, 'SST', nvarid)
312             IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Le champ <SST> est absent',1)
313
314             ierr = NF90_GET_VAR(nid,nvarid,sst_glo,start,epais)
315             IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Lecture echouee pour <SST>',1)
316         
317          END IF
318
319!****************************************************************************************
320! 4) Read albedo and rugosity for land surface, only in case of no vegetation model
321!
322!****************************************************************************************
323
324          IF (.NOT. ok_veget) THEN
325!
326! Read albedo
327             ierr = NF90_INQ_VARID(nid, 'ALB', nvarid)
328             IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Le champ <ALB> est absent',1)
329
330             ierr = NF90_GET_VAR(nid,nvarid,alb_glo,start,epais)
331             IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Lecture echouee pour <ALB>',1)
332!
333! Read rugosity
334             ierr = NF90_INQ_VARID(nid, 'RUG', nvarid)
335             IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Le champ <RUG> est absent',1)
336
337             ierr = NF90_GET_VAR(nid,nvarid,rug_glo,start,epais)
338             IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Lecture echouee pour <RUG>',1)
339
340          END IF
341
342!****************************************************************************************
343! 5) Close file and distribuate variables to all processus
344!
345!****************************************************************************************
346          ierr = NF90_CLOSE(nid)
347          IF (ierr /= NF90_NOERR) CALL abort_physic(modname,'Pb when closing file', 1)
348       ENDIF ! is_mpi_root
349
350!$OMP END MASTER
351!$OMP BARRIER
352
353       IF ( type_ocean /= 'couple') THEN
354          CALL Scatter(sst_glo,sst)
355          CALL Scatter(pct_glo(:,is_oce),pctsrf(:,is_oce))
356          CALL Scatter(pct_glo(:,is_sic),pctsrf(:,is_sic))
357          IF (read_continents .OR. itime == 1) THEN
358             CALL Scatter(pct_glo(:,is_ter),pctsrf(:,is_ter))
359             CALL Scatter(pct_glo(:,is_lic),pctsrf(:,is_lic))
360          END IF
361       END IF
362
363       IF (.NOT. ok_veget) THEN
364          CALL Scatter(alb_glo, albedo)
365          CALL Scatter(rug_glo, rugos)
366       END IF
367
368    ENDIF ! time to read
369
370  END SUBROUTINE limit_read_tot
371
372  !--------------------------------------------------------------------------------------
373  SUBROUTINE num2str(n,str)
374  !--------------------------------------------------------------------------------------
375  ! Arguments:
376    INTEGER,           INTENT(IN)  :: n
377    CHARACTER(LEN=99), INTENT(OUT) :: str
378  !--------------------------------------------------------------------------------------
379  ! Local variables:
380    INTEGER :: nn
381  !--------------------------------------------------------------------------------------
382    nn=n; str=''; DO WHILE(nn>0); str=CHAR(nn-10*(nn/10)-48)//TRIM(str); nn=nn/10; END DO
383  END SUBROUTINE num2str
384  !--------------------------------------------------------------------------------------
385
386END MODULE limit_read_mod
Note: See TracBrowser for help on using the repository browser.