source: LMDZ4/trunk/libf/phylmd/phyetat0.F @ 5434

Last change on this file since 5434 was 1403, checked in by Laurent Fairhead, 14 years ago

Merged LMDZ4V5.0-dev branch changes r1292:r1399 to trunk.

Validation:
Validation consisted in compiling the HEAD revision of the trunk,
LMDZ4V5.0-dev branch and the merged sources and running different
configurations on local and SX8 machines comparing results.

Local machine: bench configuration, 32x24x11, gfortran

  • IPSLCM5A configuration (comparison between trunk and merged sources):
    • numerical convergence on dynamical fields over 3 days
    • start files are equivalent (except for RN and PB fields)
    • daily history files equivalent
  • MH07 configuration, new physics package (comparison between LMDZ4V5.0-dev branch and merged sources):
    • numerical convergence on dynamical fields over 3 days
    • start files are equivalent (except for RN and PB fields)
    • daily history files equivalent

SX8 machine (brodie), 96x95x39 on 4 processors:

  • IPSLCM5A configuration:
    • start files are equivalent (except for RN and PB fields)
    • monthly history files equivalent
  • MH07 configuration:
    • start files are equivalent (except for RN and PB fields)
    • monthly history files equivalent

Changes to the makegcm and create_make_gcm scripts to take into account
main programs in F90 files

Fusion de la branche LMDZ4V5.0-dev (r1292:r1399) au tronc principal

Validation:
La validation a consisté à compiler la HEAD de le trunk et de la banche
LMDZ4V5.0-dev et les sources fusionnées et de faire tourner le modéle selon
différentes configurations en local et sur SX8 et de comparer les résultats

En local: 32x24x11, config bench/gfortran

  • pour une config IPSLCM5A (comparaison tronc/fusion):
    • convergence numérique sur les champs dynamiques après 3 jours
    • restart et restartphy égaux (à part sur RN et Pb)
    • fichiers histoire égaux
  • pour une config nlle physique (MH07) (comparaison LMDZ4v5.0-dev/fusion):
    • convergence numérique sur les champs dynamiques après 3 jours
    • restart et restartphy égaux
    • fichiers histoire équivalents

Sur brodie, 96x95x39 sur 4 proc:

  • pour une config IPSLCM5A:
    • restart et restartphy égaux (à part sur RN et PB)
    • pas de différence dans les fichiers histmth.nc
  • pour une config MH07
    • restart et restartphy égaux (à part sur RN et PB)
    • pas de différence dans les fichiers histmth.nc

Changement sur makegcm et create_make-gcm pour pouvoir prendre en compte des
programmes principaux en *F90
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 30.6 KB
Line 
1!
2! $Id: phyetat0.F 1403 2010-07-01 09:02:53Z fhourdin $
3!
4c
5c
6      SUBROUTINE phyetat0 (fichnom,
7     .           clesphy0,
8     .           tabcntr0)
9
10      USE dimphy
11      USE mod_grid_phy_lmdz
12      USE mod_phys_lmdz_para
13      USE iophy
14      USE ocean_cpl_mod,    ONLY : ocean_cpl_init
15      USE fonte_neige_mod,  ONLY : fonte_neige_init
16      USE pbl_surface_mod,  ONLY : pbl_surface_init
17      USE surface_data,     ONLY : type_ocean
18      USE phys_state_var_mod
19      USE iostart
20      USE write_field_phy
21      USE infotrac
22      USE traclmdz_mod,    ONLY : traclmdz_from_restart
23      USE carbon_cycle_mod,ONLY : carbon_cycle_tr, carbon_cycle_cpl
24
25      IMPLICIT none
26c======================================================================
27c Auteur(s) Z.X. Li (LMD/CNRS) date: 19930818
28c Objet: Lecture de l'etat initial pour la physique
29c======================================================================
30#include "dimensions.h"
31#include "netcdf.inc"
32#include "indicesol.h"
33#include "dimsoil.h"
34#include "clesphys.h"
35#include "temps.h"
36#include "thermcell.h"
37#include "compbl.h"
38c======================================================================
39      CHARACTER*(*) fichnom
40
41c les variables globales lues dans le fichier restart
42
43      REAL tsoil(klon,nsoilmx,nbsrf)
44      REAL tslab(klon), seaice(klon)
45      REAL qsurf(klon,nbsrf)
46      REAL qsol(klon)
47      REAL snow(klon,nbsrf)
48      REAL evap(klon,nbsrf)
49      real fder(klon)
50      REAL frugs(klon,nbsrf)
51      REAL agesno(klon,nbsrf)
52      REAL run_off_lic_0(klon)
53      REAL fractint(klon)
54      REAL trs(klon,nbtr)
55
56      CHARACTER*6 ocean_in
57      LOGICAL ok_veget_in
58
59      INTEGER        longcles
60      PARAMETER    ( longcles = 20 )
61      REAL clesphy0( longcles )
62c
63      REAL xmin, xmax
64c
65      INTEGER nid, nvarid
66      INTEGER ierr, i, nsrf, isoil ,k
67      INTEGER length
68      PARAMETER (length=100)
69      INTEGER it, iiq
70      REAL tab_cntrl(length), tabcntr0(length)
71      CHARACTER*7 str7
72      CHARACTER*2 str2
73      LOGICAL :: found
74
75c FH1D
76c     real iolat(jjm+1)
77      real iolat(jjm+1-1/iim)
78c
79c Ouvrir le fichier contenant l'etat initial:
80c
81
82     
83      CALL open_startphy(fichnom)
84     
85
86c
87c Lecture des parametres de controle:
88c
89      CALL get_var("controle",tab_cntrl)
90       
91c
92!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
93! FH 2008/05/09 On elimine toutes les clefs physiques dans la dynamique
94! Les constantes de la physiques sont lues dans la physique seulement.
95! Les egalites du type
96!             tab_cntrl( 5 )=clesphy0(1)
97! sont remplacees par
98!             clesphy0(1)=tab_cntrl( 5 )
99! On inverse aussi la logique.
100! On remplit les tab_cntrl avec les parametres lus dans les .def
101!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
102
103         DO i = 1, length
104           tabcntr0( i ) = tab_cntrl( i )
105         ENDDO
106c
107         tab_cntrl(1)=dtime
108         tab_cntrl(2)=radpas
109
110c co2_ppm : value from the previous time step
111         IF (carbon_cycle_tr .OR. carbon_cycle_cpl) THEN
112            co2_ppm = tab_cntrl(3)
113            RCO2    = co2_ppm * 1.0e-06  * 44.011/28.97
114c ELSE : keep value from .def
115         END IF
116
117c co2_ppm0 : initial value of atmospheric CO2 (from create_etat0_limit.e .def)
118         co2_ppm0   = tab_cntrl(16)
119
120         solaire_etat0      = tab_cntrl(4)
121         tab_cntrl(5)=iflag_con
122         tab_cntrl(6)=nbapp_rad
123
124         if (cycle_diurne) tab_cntrl( 7) =1.
125         if (soil_model) tab_cntrl( 8) =1.
126         if (new_oliq) tab_cntrl( 9) =1.
127         if (ok_orodr) tab_cntrl(10) =1.
128         if (ok_orolf) tab_cntrl(11) =1.
129         if (ok_limitvrai) tab_cntrl(12) =1.
130
131
132      itau_phy = tab_cntrl(15)
133
134       
135
136         IF( clesphy0(1).NE.tab_cntrl( 5 ) )  THEN
137             clesphy0(1)=tab_cntrl( 5 )
138         ENDIF
139
140         IF( clesphy0(2).NE.tab_cntrl( 6 ) )  THEN
141             clesphy0(2)=tab_cntrl( 6 )
142         ENDIF
143
144         IF( clesphy0(3).NE.tab_cntrl( 7 ) )  THEN
145             clesphy0(3)=tab_cntrl( 7 )
146         ENDIF
147
148         IF( clesphy0(4).NE.tab_cntrl( 8 ) )  THEN
149             clesphy0(4)=tab_cntrl( 8 )
150         ENDIF
151
152         IF( clesphy0(5).NE.tab_cntrl( 9 ) )  THEN
153             clesphy0(5)=tab_cntrl( 9 )
154         ENDIF
155
156         IF( clesphy0(6).NE.tab_cntrl( 10 ) )  THEN
157             clesphy0(6)=tab_cntrl( 10 )
158         ENDIF
159
160         IF( clesphy0(7).NE.tab_cntrl( 11 ) )  THEN
161             clesphy0(7)=tab_cntrl( 11 )
162         ENDIF
163
164         IF( clesphy0(8).NE.tab_cntrl( 12 ) )  THEN
165             clesphy0(8)=tab_cntrl( 12 )
166         ENDIF
167
168
169c
170c Lecture des latitudes (coordonnees):
171c
172      CALL get_field("latitude",rlat)
173
174c
175c Lecture des longitudes (coordonnees):
176c
177      CALL get_field("longitude",rlon)
178
179C
180C
181C Lecture du masque terre mer
182C
183      CALL get_field("masque",zmasq,found)
184      IF (.NOT. found) THEN
185        PRINT*, 'phyetat0: Le champ <masque> est absent'
186        PRINT *, 'fichier startphy non compatible avec phyetat0'
187      ENDIF
188
189       
190C Lecture des fractions pour chaque sous-surface
191C
192C initialisation des sous-surfaces
193C
194      pctsrf = 0.
195C
196C fraction de terre
197C
198
199      CALL get_field("FTER",pctsrf(:,is_ter),found)
200      IF (.NOT. found) PRINT*, 'phyetat0: Le champ <FTER> est absent'
201
202C
203C fraction de glace de terre
204C
205      CALL get_field("FLIC",pctsrf(:,is_lic),found)
206      IF (.NOT. found) PRINT*, 'phyetat0: Le champ <FLIC> est absent'
207
208C
209C fraction d'ocean
210C
211      CALL get_field("FOCE",pctsrf(:,is_oce),found)
212      IF (.NOT. found) PRINT*, 'phyetat0: Le champ <FOCE> est absent'
213
214C
215C fraction glace de mer
216C
217      CALL get_field("FSIC",pctsrf(:,is_sic),found)
218      IF (.NOT. found) PRINT*, 'phyetat0: Le champ <FSIC> est absent'
219
220C
221C  Verification de l'adequation entre le masque et les sous-surfaces
222C
223      fractint( 1 : klon) = pctsrf(1 : klon, is_ter)
224     $    + pctsrf(1 : klon, is_lic)
225      DO i = 1 , klon
226        IF ( abs(fractint(i) - zmasq(i) ) .GT. EPSFRA ) THEN
227            WRITE(*,*) 'phyetat0: attention fraction terre pas ',
228     $          'coherente ', i, zmasq(i), pctsrf(i, is_ter)
229     $          ,pctsrf(i, is_lic)
230            WRITE(*,*) 'Je force la coherence zmasq=fractint'
231            zmasq(i) = fractint(i)
232        ENDIF
233      END DO
234      fractint (1 : klon) =  pctsrf(1 : klon, is_oce)
235     $    + pctsrf(1 : klon, is_sic)
236      DO i = 1 , klon
237        IF ( abs( fractint(i) - (1. - zmasq(i))) .GT. EPSFRA ) THEN
238            WRITE(*,*) 'phyetat0 attention fraction ocean pas ',
239     $          'coherente ', i, zmasq(i) , pctsrf(i, is_oce)
240     $          ,pctsrf(i, is_sic)
241            WRITE(*,*) 'Je force la coherence zmasq=fractint'
242            zmasq(i) = fractint(i)
243        ENDIF
244      END DO
245
246C
247c Lecture des temperatures du sol:
248c
249
250       CALL get_field("TS",ftsol(:,1),found)
251       IF (.NOT. found) THEN
252         PRINT*, 'phyetat0: Le champ <TS> est absent'
253         PRINT*, '          Mais je vais essayer de lire TS**'
254         DO nsrf = 1, nbsrf
255           IF (nsrf.GT.99) THEN
256             PRINT*, "Trop de sous-mailles"
257             CALL abort
258           ENDIF
259           WRITE(str2,'(i2.2)') nsrf
260           CALL get_field("TS"//str2,ftsol(:,nsrf))
261
262           xmin = 1.0E+20
263           xmax = -1.0E+20
264           DO i = 1, klon
265              xmin = MIN(ftsol(i,nsrf),xmin)
266              xmax = MAX(ftsol(i,nsrf),xmax)
267           ENDDO
268           PRINT*,'Temperature du sol TS**:', nsrf, xmin, xmax
269         ENDDO
270      ELSE
271         PRINT*, 'phyetat0: Le champ <TS> est present'
272         PRINT*, '          J ignore donc les autres temperatures TS**'
273         xmin = 1.0E+20
274         xmax = -1.0E+20
275         DO i = 1, klon
276            xmin = MIN(ftsol(i,1),xmin)
277            xmax = MAX(ftsol(i,1),xmax)
278         ENDDO
279         PRINT*,'Temperature du sol <TS>', xmin, xmax
280         DO nsrf = 2, nbsrf
281         DO i = 1, klon
282            ftsol(i,nsrf) = ftsol(i,1)
283         ENDDO
284         ENDDO
285      ENDIF
286
287c
288c Lecture des temperatures du sol profond:
289c
290      DO nsrf = 1, nbsrf
291        DO isoil=1, nsoilmx
292          IF (isoil.GT.99 .AND. nsrf.GT.99) THEN
293            PRINT*, "Trop de couches ou sous-mailles"
294            CALL abort
295          ENDIF
296          WRITE(str7,'(i2.2,"srf",i2.2)') isoil, nsrf
297         
298          CALL get_field('Tsoil'//str7,tsoil(:,isoil,nsrf),found)
299          IF (.NOT. found) THEN
300            PRINT*, "phyetat0: Le champ <Tsoil"//str7//"> est absent"
301            PRINT*, "          Il prend donc la valeur de surface"
302            DO i=1, klon
303               tsoil(i,isoil,nsrf)=ftsol(i,nsrf)
304            ENDDO
305          ENDIF
306        ENDDO
307      ENDDO
308c
309c Lecture de l'humidite de l'air juste au dessus du sol:
310c
311
312      CALL get_field("QS",qsurf(:,1),found)
313      IF (.NOT. found) THEN
314         PRINT*, 'phyetat0: Le champ <QS> est absent'
315         PRINT*, '          Mais je vais essayer de lire QS**'
316         DO nsrf = 1, nbsrf
317           IF (nsrf.GT.99) THEN
318             PRINT*, "Trop de sous-mailles"
319             CALL abort
320           ENDIF
321           WRITE(str2,'(i2.2)') nsrf
322           CALL get_field("QS"//str2,qsurf(:,nsrf))
323           xmin = 1.0E+20
324           xmax = -1.0E+20
325           DO i = 1, klon
326              xmin = MIN(qsurf(i,nsrf),xmin)
327              xmax = MAX(qsurf(i,nsrf),xmax)
328           ENDDO
329           PRINT*,'Humidite pres du sol QS**:', nsrf, xmin, xmax
330         ENDDO
331      ELSE
332         PRINT*, 'phyetat0: Le champ <QS> est present'
333         PRINT*, '          J ignore donc les autres humidites QS**'
334         xmin = 1.0E+20
335         xmax = -1.0E+20
336         DO i = 1, klon
337            xmin = MIN(qsurf(i,1),xmin)
338            xmax = MAX(qsurf(i,1),xmax)
339         ENDDO
340         PRINT*,'Humidite pres du sol <QS>', xmin, xmax
341         DO nsrf = 2, nbsrf
342           DO i = 1, klon
343             qsurf(i,nsrf) = qsurf(i,1)
344           ENDDO
345         ENDDO
346      ENDIF
347
348C
349C Eau dans le sol (pour le modele de sol "bucket")
350C
351      CALL get_field("QSOL",qsol,found)
352      IF (.NOT. found) THEN
353        PRINT*, 'phyetat0: Le champ <QSOL> est absent'
354        PRINT*, '          Valeur par defaut nulle'
355          qsol(:)=0.
356      ENDIF
357
358      xmin = 1.0E+20
359      xmax = -1.0E+20
360      DO i = 1, klon
361        xmin = MIN(qsol(i),xmin)
362        xmax = MAX(qsol(i),xmax)
363      ENDDO
364      PRINT*,'Eau dans le sol (mm) <QSOL>', xmin, xmax
365
366c
367c Lecture de neige au sol:
368c
369
370      CALL get_field("SNOW",snow(:,1),found)
371      IF (.NOT. found) THEN
372        PRINT*, 'phyetat0: Le champ <SNOW> est absent'
373        PRINT*, '          Mais je vais essayer de lire SNOW**'
374        DO nsrf = 1, nbsrf
375          IF (nsrf.GT.99) THEN
376            PRINT*, "Trop de sous-mailles"
377            CALL abort
378          ENDIF
379          WRITE(str2,'(i2.2)') nsrf
380          CALL get_field( "SNOW"//str2,snow(:,nsrf))
381          xmin = 1.0E+20
382          xmax = -1.0E+20
383          DO i = 1, klon
384            xmin = MIN(snow(i,nsrf),xmin)
385            xmax = MAX(snow(i,nsrf),xmax)
386          ENDDO
387          PRINT*,'Neige du sol SNOW**:', nsrf, xmin, xmax
388        ENDDO
389      ELSE
390         PRINT*, 'phyetat0: Le champ <SNOW> est present'
391         PRINT*, '          J ignore donc les autres neiges SNOW**'
392         xmin = 1.0E+20
393         xmax = -1.0E+20
394         DO i = 1, klon
395            xmin = MIN(snow(i,1),xmin)
396            xmax = MAX(snow(i,1),xmax)
397         ENDDO
398         PRINT*,'Neige du sol <SNOW>', xmin, xmax
399         DO nsrf = 2, nbsrf
400         DO i = 1, klon
401            snow(i,nsrf) = snow(i,1)
402         ENDDO
403         ENDDO
404      ENDIF
405c
406c Lecture de albedo de l'interval visible au sol:
407c
408      CALL get_field("ALBE",falb1(:,1),found)
409      IF (.NOT. found) THEN
410         PRINT*, 'phyetat0: Le champ <ALBE> est absent'
411         PRINT*, '          Mais je vais essayer de lire ALBE**'
412         DO nsrf = 1, nbsrf
413           IF (nsrf.GT.99) THEN
414             PRINT*, "Trop de sous-mailles"
415             CALL abort
416           ENDIF
417           WRITE(str2,'(i2.2)') nsrf
418           CALL get_field("ALBE"//str2,falb1(:,nsrf))
419           xmin = 1.0E+20
420           xmax = -1.0E+20
421           DO i = 1, klon
422              xmin = MIN(falb1(i,nsrf),xmin)
423              xmax = MAX(falb1(i,nsrf),xmax)
424           ENDDO
425           PRINT*,'Albedo du sol ALBE**:', nsrf, xmin, xmax
426         ENDDO
427      ELSE
428         PRINT*, 'phyetat0: Le champ <ALBE> est present'
429         PRINT*, '          J ignore donc les autres ALBE**'
430         xmin = 1.0E+20
431         xmax = -1.0E+20
432         DO i = 1, klon
433            xmin = MIN(falb1(i,1),xmin)
434            xmax = MAX(falb1(i,1),xmax)
435         ENDDO
436         PRINT*,'Neige du sol <ALBE>', xmin, xmax
437         DO nsrf = 2, nbsrf
438           DO i = 1, klon
439            falb1(i,nsrf) = falb1(i,1)
440           ENDDO
441         ENDDO
442      ENDIF
443
444c
445c Lecture de albedo au sol dans l'interval proche infra-rouge:
446c
447      CALL get_field("ALBLW",falb2(:,1),found)
448      IF (.NOT. found) THEN
449         PRINT*, 'phyetat0: Le champ <ALBLW> est absent'
450         PRINT*, '          Mais je vais prendre ALBE**'
451         DO nsrf = 1, nbsrf
452           DO i = 1, klon
453             falb2(i,nsrf) = falb1(i,nsrf)
454           ENDDO
455         ENDDO
456      ELSE
457         PRINT*, 'phyetat0: Le champ <ALBLW> est present'
458         PRINT*, '          J ignore donc les autres ALBLW**'
459         xmin = 1.0E+20
460         xmax = -1.0E+20
461         DO i = 1, klon
462            xmin = MIN(falb2(i,1),xmin)
463            xmax = MAX(falb2(i,1),xmax)
464         ENDDO
465         PRINT*,'Neige du sol <ALBLW>', xmin, xmax
466         DO nsrf = 2, nbsrf
467           DO i = 1, klon
468             falb2(i,nsrf) = falb2(i,1)
469           ENDDO
470         ENDDO
471      ENDIF
472c
473c Lecture de evaporation: 
474c
475      CALL get_field("EVAP",evap(:,1),found)
476      IF (.NOT. found) THEN
477         PRINT*, 'phyetat0: Le champ <EVAP> est absent'
478         PRINT*, '          Mais je vais essayer de lire EVAP**'
479         DO nsrf = 1, nbsrf
480           IF (nsrf.GT.99) THEN
481             PRINT*, "Trop de sous-mailles"
482             CALL abort
483           ENDIF
484           WRITE(str2,'(i2.2)') nsrf
485           CALL get_field("EVAP"//str2, evap(:,nsrf))
486           xmin = 1.0E+20
487           xmax = -1.0E+20
488           DO i = 1, klon
489              xmin = MIN(evap(i,nsrf),xmin)
490              xmax = MAX(evap(i,nsrf),xmax)
491           ENDDO
492           PRINT*,'evap du sol EVAP**:', nsrf, xmin, xmax
493         ENDDO
494      ELSE
495         PRINT*, 'phyetat0: Le champ <EVAP> est present'
496         PRINT*, '          J ignore donc les autres EVAP**'
497         xmin = 1.0E+20
498         xmax = -1.0E+20
499         DO i = 1, klon
500            xmin = MIN(evap(i,1),xmin)
501            xmax = MAX(evap(i,1),xmax)
502         ENDDO
503         PRINT*,'Evap du sol <EVAP>', xmin, xmax
504         DO nsrf = 2, nbsrf
505         DO i = 1, klon
506            evap(i,nsrf) = evap(i,1)
507         ENDDO
508         ENDDO
509      ENDIF
510c
511c Lecture precipitation liquide:
512c
513      CALL get_field("rain_f",rain_fall)
514      xmin = 1.0E+20
515      xmax = -1.0E+20
516      DO i = 1, klon
517         xmin = MIN(rain_fall(i),xmin)
518         xmax = MAX(rain_fall(i),xmax)
519      ENDDO
520      PRINT*,'Precipitation liquide rain_f:', xmin, xmax
521c
522c Lecture precipitation solide:
523c
524      CALL get_field("snow_f",snow_fall)
525      xmin = 1.0E+20
526      xmax = -1.0E+20
527      DO i = 1, klon
528         xmin = MIN(snow_fall(i),xmin)
529         xmax = MAX(snow_fall(i),xmax)
530      ENDDO
531      PRINT*,'Precipitation solide snow_f:', xmin, xmax
532c
533c Lecture rayonnement solaire au sol:
534c
535      CALL get_field("solsw",solsw,found)
536      IF (.NOT. found) THEN
537         PRINT*, 'phyetat0: Le champ <solsw> est absent'
538         PRINT*, 'mis a zero'
539         solsw(:) = 0.
540      ENDIF
541      xmin = 1.0E+20
542      xmax = -1.0E+20
543      DO i = 1, klon
544         xmin = MIN(solsw(i),xmin)
545         xmax = MAX(solsw(i),xmax)
546      ENDDO
547      PRINT*,'Rayonnement solaire au sol solsw:', xmin, xmax
548c
549c Lecture rayonnement IF au sol:
550c
551      CALL get_field("sollw",sollw,found)
552      IF (.NOT. found) THEN
553         PRINT*, 'phyetat0: Le champ <sollw> est absent'
554         PRINT*, 'mis a zero'
555         sollw = 0.
556      ENDIF
557      xmin = 1.0E+20
558      xmax = -1.0E+20
559      DO i = 1, klon
560         xmin = MIN(sollw(i),xmin)
561         xmax = MAX(sollw(i),xmax)
562      ENDDO
563      PRINT*,'Rayonnement IF au sol sollw:', xmin, xmax
564     
565c
566c Lecture derive des flux:
567c
568      CALL get_field("fder",fder,found)
569      IF (.NOT. found) THEN
570         PRINT*, 'phyetat0: Le champ <fder> est absent'
571         PRINT*, 'mis a zero'
572         fder = 0.
573      ENDIF
574      xmin = 1.0E+20
575      xmax = -1.0E+20
576      DO i = 1, klon
577         xmin = MIN(fder(i),xmin)
578         xmax = MAX(fder(i),xmax)
579      ENDDO
580      PRINT*,'Derive des flux fder:', xmin, xmax
581
582c
583c Lecture du rayonnement net au sol:
584c
585      CALL get_field("RADS",radsol)
586      xmin = 1.0E+20
587      xmax = -1.0E+20
588      DO i = 1, klon
589         xmin = MIN(radsol(i),xmin)
590         xmax = MAX(radsol(i),xmax)
591      ENDDO
592      PRINT*,'Rayonnement net au sol radsol:', xmin, xmax
593c
594c Lecture de la longueur de rugosite
595c
596c
597      CALL get_field("RUG",frugs(:,1),found)
598      IF (.NOT. found) THEN
599         PRINT*, 'phyetat0: Le champ <RUG> est absent'
600         PRINT*, '          Mais je vais essayer de lire RUG**'
601         DO nsrf = 1, nbsrf
602           IF (nsrf.GT.99) THEN
603             PRINT*, "Trop de sous-mailles"
604             CALL abort
605           ENDIF
606           WRITE(str2,'(i2.2)') nsrf
607           CALL get_field("RUG"//str2,frugs(:,nsrf))
608           xmin = 1.0E+20
609           xmax = -1.0E+20
610           DO i = 1, klon
611              xmin = MIN(frugs(i,nsrf),xmin)
612              xmax = MAX(frugs(i,nsrf),xmax)
613           ENDDO
614           PRINT*,'rugosite du sol RUG**:', nsrf, xmin, xmax
615         ENDDO
616      ELSE
617         PRINT*, 'phyetat0: Le champ <RUG> est present'
618         PRINT*, '          J ignore donc les autres RUG**'
619         xmin = 1.0E+20
620         xmax = -1.0E+20
621         DO i = 1, klon
622            xmin = MIN(frugs(i,1),xmin)
623            xmax = MAX(frugs(i,1),xmax)
624         ENDDO
625         PRINT*,'rugosite <RUG>', xmin, xmax
626         DO nsrf = 2, nbsrf
627         DO i = 1, klon
628            frugs(i,nsrf) = frugs(i,1)
629         ENDDO
630         ENDDO
631      ENDIF
632
633c
634c Lecture de l'age de la neige:
635c
636      CALL get_field("AGESNO",agesno(:,1),found)
637      IF (.NOT. found) THEN
638         PRINT*, 'phyetat0: Le champ <AGESNO> est absent'
639         PRINT*, '          Mais je vais essayer de lire AGESNO**'
640         DO nsrf = 1, nbsrf
641           IF (nsrf.GT.99) THEN
642             PRINT*, "Trop de sous-mailles"
643             CALL abort
644           ENDIF
645           WRITE(str2,'(i2.2)') nsrf
646           CALL get_field("AGESNO"//str2,agesno(:,nsrf),found)
647           IF (.NOT. found) THEN
648              PRINT*, "phyetat0: Le champ <AGESNO"//str2//"> est absent"
649              agesno = 50.0
650           ENDIF
651           xmin = 1.0E+20
652           xmax = -1.0E+20
653           DO i = 1, klon
654              xmin = MIN(agesno(i,nsrf),xmin)
655              xmax = MAX(agesno(i,nsrf),xmax)
656           ENDDO
657           PRINT*,'Age de la neige AGESNO**:', nsrf, xmin, xmax
658         ENDDO
659      ELSE
660         PRINT*, 'phyetat0: Le champ <AGESNO> est present'
661         PRINT*, '          J ignore donc les autres AGESNO**'
662         xmin = 1.0E+20
663         xmax = -1.0E+20
664         DO i = 1, klon
665            xmin = MIN(agesno(i,1),xmin)
666            xmax = MAX(agesno(i,1),xmax)
667         ENDDO
668         PRINT*,'Age de la neige <AGESNO>', xmin, xmax
669         DO nsrf = 2, nbsrf
670         DO i = 1, klon
671            agesno(i,nsrf) = agesno(i,1)
672         ENDDO
673         ENDDO
674      ENDIF
675
676c
677      CALL get_field("ZMEA", zmea)
678      xmin = 1.0E+20
679      xmax = -1.0E+20
680      DO i = 1, klon
681         xmin = MIN(zmea(i),xmin)
682         xmax = MAX(zmea(i),xmax)
683      ENDDO
684      PRINT*,'OROGRAPHIE SOUS-MAILLE zmea:', xmin, xmax
685c
686c
687      CALL get_field("ZSTD",zstd)
688      xmin = 1.0E+20
689      xmax = -1.0E+20
690      DO i = 1, klon
691         xmin = MIN(zstd(i),xmin)
692         xmax = MAX(zstd(i),xmax)
693      ENDDO
694      PRINT*,'OROGRAPHIE SOUS-MAILLE zstd:', xmin, xmax
695c
696c
697      CALL get_field("ZSIG",zsig)
698      xmin = 1.0E+20
699      xmax = -1.0E+20
700      DO i = 1, klon
701         xmin = MIN(zsig(i),xmin)
702         xmax = MAX(zsig(i),xmax)
703      ENDDO
704      PRINT*,'OROGRAPHIE SOUS-MAILLE zsig:', xmin, xmax
705c
706c
707      CALL get_field("ZGAM",zgam)
708      xmin = 1.0E+20
709      xmax = -1.0E+20
710      DO i = 1, klon
711         xmin = MIN(zgam(i),xmin)
712         xmax = MAX(zgam(i),xmax)
713      ENDDO
714      PRINT*,'OROGRAPHIE SOUS-MAILLE zgam:', xmin, xmax
715c
716c
717      CALL get_field("ZTHE",zthe)
718      xmin = 1.0E+20
719      xmax = -1.0E+20
720      DO i = 1, klon
721         xmin = MIN(zthe(i),xmin)
722         xmax = MAX(zthe(i),xmax)
723      ENDDO
724      PRINT*,'OROGRAPHIE SOUS-MAILLE zthe:', xmin, xmax
725c
726c
727      CALL get_field("ZPIC",zpic)
728      xmin = 1.0E+20
729      xmax = -1.0E+20
730      DO i = 1, klon
731         xmin = MIN(zpic(i),xmin)
732         xmax = MAX(zpic(i),xmax)
733      ENDDO
734      PRINT*,'OROGRAPHIE SOUS-MAILLE zpic:', xmin, xmax
735c
736      CALL get_field("ZVAL",zval)
737      xmin = 1.0E+20
738      xmax = -1.0E+20
739      DO i = 1, klon
740         xmin = MIN(zval(i),xmin)
741         xmax = MAX(zval(i),xmax)
742      ENDDO
743      PRINT*,'OROGRAPHIE SOUS-MAILLE zval:', xmin, xmax
744c
745c
746      CALL get_field("RUGSREL",rugoro)
747      xmin = 1.0E+20
748      xmax = -1.0E+20
749      DO i = 1, klon
750         xmin = MIN(rugoro(i),xmin)
751         xmax = MAX(rugoro(i),xmax)
752      ENDDO
753      PRINT*,'Rugosite relief (ecart-type) rugsrel:', xmin, xmax
754c
755c
756     
757c
758      ancien_ok = .TRUE.
759
760      CALL get_field("TANCIEN",t_ancien,found)
761      IF (.NOT. found) THEN
762         PRINT*, "phyetat0: Le champ <TANCIEN> est absent"
763         PRINT*, "Depart legerement fausse. Mais je continue"
764         ancien_ok = .FALSE.
765      ENDIF
766
767
768      CALL get_field("QANCIEN",q_ancien,found)
769      IF (.NOT. found) THEN
770         PRINT*, "phyetat0: Le champ <QANCIEN> est absent"
771         PRINT*, "Depart legerement fausse. Mais je continue"
772         ancien_ok = .FALSE.
773      ENDIF
774
775      u_ancien = 0.0   !AXC: We don't have u_ancien and v_ancien in the start
776      v_ancien = 0.0   !AXC: files, therefore they have to be initialized.
777c
778
779      clwcon=0.
780      CALL get_field("CLWCON",clwcon(:,1),found)
781      IF (.NOT. found) THEN
782         PRINT*, "phyetat0: Le champ CLWCON est absent"
783         PRINT*, "Depart legerement fausse. Mais je continue"
784      ENDIF
785      xmin = 1.0E+20
786      xmax = -1.0E+20
787      xmin = MINval(clwcon)
788      xmax = MAXval(clwcon)
789      PRINT*,'Eau liquide convective (ecart-type) clwcon:', xmin, xmax
790c
791      rnebcon = 0.
792      CALL get_field("RNEBCON",rnebcon(:,1),found)
793      IF (.NOT. found) THEN
794         PRINT*, "phyetat0: Le champ RNEBCON est absent"
795         PRINT*, "Depart legerement fausse. Mais je continue"
796      ENDIF
797      xmin = 1.0E+20
798      xmax = -1.0E+20
799      xmin = MINval(rnebcon)
800      xmax = MAXval(rnebcon)
801      PRINT*,'Nebulosite convective (ecart-type) rnebcon:', xmin, xmax
802
803c
804c Lecture ratqs
805c
806      ratqs=0.
807      CALL get_field("RATQS",ratqs(:,1),found)
808      IF (.NOT. found) THEN
809         PRINT*, "phyetat0: Le champ <RATQS> est absent"
810         PRINT*, "Depart legerement fausse. Mais je continue"
811      ENDIF
812      xmin = 1.0E+20
813      xmax = -1.0E+20
814      xmin = MINval(ratqs)
815      xmax = MAXval(ratqs)
816      PRINT*,'(ecart-type) ratqs:', xmin, xmax
817c
818c Lecture run_off_lic_0
819c
820      CALL get_field("RUNOFFLIC0",run_off_lic_0,found)
821      IF (.NOT. found) THEN
822         PRINT*, "phyetat0: Le champ <RUNOFFLIC0> est absent"
823         PRINT*, "Depart legerement fausse. Mais je continue"
824         run_off_lic_0 = 0.
825      ENDIF
826      xmin = 1.0E+20
827      xmax = -1.0E+20
828      xmin = MINval(run_off_lic_0)
829      xmax = MAXval(run_off_lic_0)
830      PRINT*,'(ecart-type) run_off_lic_0:', xmin, xmax
831
832
833c Lecture de l'energie cinetique turbulente
834c
835
836      IF (iflag_pbl>1) then
837        DO nsrf = 1, nbsrf
838          IF (nsrf.GT.99) THEN
839            PRINT*, "Trop de sous-mailles"
840            CALL abort
841          ENDIF
842          WRITE(str2,'(i2.2)') nsrf
843          CALL get_field("TKE"//str2,pbl_tke(:,1:klev,nsrf),found)
844          IF (.NOT. found) THEN
845            PRINT*, "phyetat0: <TKE"//str2//"> est absent"
846            pbl_tke(:,:,nsrf)=1.e-8
847          ENDIF
848          xmin = 1.0E+20
849          xmax = -1.0E+20
850          DO k = 1, klev
851            DO i = 1, klon
852              xmin = MIN(pbl_tke(i,k,nsrf),xmin)
853              xmax = MAX(pbl_tke(i,k,nsrf),xmax)
854            ENDDO
855          ENDDO
856          PRINT*,'Temperature du sol TKE**:', nsrf, xmin, xmax
857        ENDDO
858      ENDIF
859c
860c zmax0
861      CALL get_field("ZMAX0",zmax0,found)
862      IF (.NOT. found) THEN
863        PRINT*, "phyetat0: Le champ <ZMAX0> est absent"
864        PRINT*, "Depart legerement fausse. Mais je continue"
865        zmax0=40.
866      ENDIF
867      xmin = 1.0E+20
868      xmax = -1.0E+20
869      xmin = MINval(zmax0)
870      xmax = MAXval(zmax0)
871      PRINT*,'(ecart-type) zmax0:', xmin, xmax
872c
873c           f0(ig)=1.e-5
874c f0
875      CALL get_field("F0",f0,found)
876      IF (.NOT. found) THEN
877         PRINT*, "phyetat0: Le champ <f0> est absent"
878         PRINT*, "Depart legerement fausse. Mais je continue"
879         f0=1.e-5
880      ENDIF
881      xmin = 1.0E+20
882      xmax = -1.0E+20
883      xmin = MINval(f0)
884      xmax = MAXval(f0)
885      PRINT*,'(ecart-type) f0:', xmin, xmax
886c
887c ema_work1
888c
889      CALL get_field("EMA_WORK1",ema_work1,found)
890      IF (.NOT. found) THEN
891        PRINT*, "phyetat0: Le champ <EMA_WORK1> est absent"
892        PRINT*, "Depart legerement fausse. Mais je continue"
893        ema_work1=0.
894      ELSE
895        xmin = 1.0E+20
896        xmax = -1.0E+20
897        DO k = 1, klev
898          DO i = 1, klon
899            xmin = MIN(ema_work1(i,k),xmin)
900            xmax = MAX(ema_work1(i,k),xmax)
901          ENDDO
902        ENDDO
903        PRINT*,'ema_work1:', xmin, xmax
904      ENDIF
905c
906c ema_work2
907c
908      CALL get_field("EMA_WORK2",ema_work2,found)
909      IF (.NOT. found) THEN
910        PRINT*, "phyetat0: Le champ <EMA_WORK2> est absent"
911        PRINT*, "Depart legerement fausse. Mais je continue"
912        ema_work2=0.
913      ELSE
914        xmin = 1.0E+20
915        xmax = -1.0E+20
916        DO k = 1, klev
917          DO i = 1, klon
918            xmin = MIN(ema_work2(i,k),xmin)
919            xmax = MAX(ema_work2(i,k),xmax)
920          ENDDO
921        ENDDO
922        PRINT*,'ema_work2:', xmin, xmax
923      ENDIF
924c
925c wake_deltat
926c
927      CALL get_field("WAKE_DELTAT",wake_deltat,found)
928      IF (.NOT. found) THEN
929        PRINT*, "phyetat0: Le champ <WAKE_DELTAT> est absent"
930        PRINT*, "Depart legerement fausse. Mais je continue"
931        wake_deltat=0.
932      ELSE
933        xmin = 1.0E+20
934        xmax = -1.0E+20
935        DO k = 1, klev
936          DO i = 1, klon
937            xmin = MIN(wake_deltat(i,k),xmin)
938            xmax = MAX(wake_deltat(i,k),xmax)
939          ENDDO
940        ENDDO
941        PRINT*,'wake_deltat:', xmin, xmax
942      ENDIF
943c
944c wake_deltaq
945c   
946      CALL get_field("WAKE_DELTAQ",wake_deltaq,found)
947      IF (.NOT. found) THEN
948        PRINT*, "phyetat0: Le champ <WAKE_DELTAQ> est absent"
949        PRINT*, "Depart legerement fausse. Mais je continue"
950        wake_deltaq=0.
951      ELSE
952        xmin = 1.0E+20
953        xmax = -1.0E+20
954        DO k = 1, klev
955          DO i = 1, klon
956            xmin = MIN(wake_deltaq(i,k),xmin)
957            xmax = MAX(wake_deltaq(i,k),xmax)
958          ENDDO
959        ENDDO
960        PRINT*,'wake_deltaq:', xmin, xmax
961      ENDIF
962c
963c wake_s
964c
965      CALL get_field("WAKE_S",wake_s,found)
966      IF (.NOT. found) THEN
967        PRINT*, "phyetat0: Le champ <WAKE_S> est absent"
968        PRINT*, "Depart legerement fausse. Mais je continue"
969        wake_s=0.
970      ENDIF
971      xmin = 1.0E+20
972      xmax = -1.0E+20
973      xmin = MINval(wake_s)
974      xmax = MAXval(wake_s)
975      PRINT*,'(ecart-type) wake_s:', xmin, xmax
976c
977c wake_cstar
978c
979      CALL get_field("WAKE_CSTAR",wake_cstar,found)
980      IF (.NOT. found) THEN
981         PRINT*, "phyetat0: Le champ <WAKE_CSTAR> est absent"
982         PRINT*, "Depart legerement fausse. Mais je continue"
983         wake_cstar=0.
984      ENDIF
985      xmin = 1.0E+20
986      xmax = -1.0E+20
987      xmin = MINval(wake_cstar)
988      xmax = MAXval(wake_cstar)
989      PRINT*,'(ecart-type) wake_cstar:', xmin, xmax
990c
991c wake_pe
992c
993      CALL get_field("WAKE_PE",wake_pe,found)
994      IF (.NOT. found) THEN
995         PRINT*, "phyetat0: Le champ <WAKE_PE> est absent"
996         PRINT*, "Depart legerement fausse. Mais je continue"
997         wake_pe=0.
998      ENDIF
999      xmin = 1.0E+20
1000      xmax = -1.0E+20
1001      xmin = MINval(wake_pe)
1002      xmax = MAXval(wake_pe)
1003      PRINT*,'(ecart-type) wake_pe:', xmin, xmax
1004c
1005c wake_fip
1006c
1007      CALL get_field("WAKE_FIP",wake_fip,found)
1008      IF (.NOT. found) THEN
1009         PRINT*, "phyetat0: Le champ <WAKE_FIP> est absent"
1010         PRINT*, "Depart legerement fausse. Mais je continue"
1011         wake_fip=0.
1012      ENDIF
1013      xmin = 1.0E+20
1014      xmax = -1.0E+20
1015      xmin = MINval(wake_fip)
1016      xmax = MAXval(wake_fip)
1017      PRINT*,'(ecart-type) wake_fip:', xmin, xmax
1018c
1019c  thermiques
1020c
1021
1022      CALL get_field("FM_THERM",fm_therm,found)
1023      IF (.NOT. found) THEN
1024         PRINT*, "phyetat0: Le champ <fm_therm> est absent"
1025         PRINT*, "Depart legerement fausse. Mais je continue"
1026         fm_therm=0.
1027      ENDIF
1028      xmin = 1.0E+20
1029      xmax = -1.0E+20
1030      xmin = MINval(fm_therm)
1031      xmax = MAXval(fm_therm)
1032      PRINT*,'(ecart-type) fm_therm:', xmin, xmax
1033
1034      CALL get_field("ENTR_THERM",entr_therm,found)
1035      IF (.NOT. found) THEN
1036         PRINT*, "phyetat0: Le champ <entr_therm> est absent"
1037         PRINT*, "Depart legerement fausse. Mais je continue"
1038         entr_therm=0.
1039      ENDIF
1040      xmin = 1.0E+20
1041      xmax = -1.0E+20
1042      xmin = MINval(entr_therm)
1043      xmax = MAXval(entr_therm)
1044      PRINT*,'(ecart-type) entr_therm:', xmin, xmax
1045
1046      CALL get_field("DETR_THERM",detr_therm,found)
1047      IF (.NOT. found) THEN
1048         PRINT*, "phyetat0: Le champ <detr_therm> est absent"
1049         PRINT*, "Depart legerement fausse. Mais je continue"
1050         detr_therm=0.
1051      ENDIF
1052      xmin = 1.0E+20
1053      xmax = -1.0E+20
1054      xmin = MINval(detr_therm)
1055      xmax = MAXval(detr_therm)
1056      PRINT*,'(ecart-type) detr_therm:', xmin, xmax
1057
1058
1059
1060c
1061c Read and send field trs to traclmdz
1062c
1063      IF (type_trac == 'lmdz') THEN
1064         DO it=1,nbtr
1065            iiq=niadv(it+2)
1066            CALL get_field("trs_"//tname(iiq),trs(:,it),found)
1067            IF (.NOT. found) THEN
1068               PRINT*,
1069     $           "phyetat0: Le champ <trs_"//tname(iiq)//"> est absent"
1070               PRINT*, "Depart legerement fausse. Mais je continue"
1071               trs(:,it) = 0.
1072            ENDIF
1073            xmin = 1.0E+20
1074            xmax = -1.0E+20
1075            xmin = MINval(trs(:,it))
1076            xmax = MAXval(trs(:,it))
1077            PRINT*,"(ecart-type) trs_"//tname(iiq)//" :", xmin, xmax
1078
1079         END DO
1080         
1081         CALL traclmdz_from_restart(trs)
1082      END IF
1083
1084
1085c on ferme le fichier
1086      CALL close_startphy
1087
1088      CALL init_iophy_new(rlat,rlon)
1089       
1090
1091c
1092c Initialize module pbl_surface_mod
1093c
1094      CALL pbl_surface_init(qsol, fder, snow, qsurf,
1095     $     evap, frugs, agesno, tsoil)
1096
1097c Initialize module ocean_cpl_mod for the case of coupled ocean
1098      IF ( type_ocean == 'couple' ) THEN
1099         CALL ocean_cpl_init(dtime, rlon, rlat)
1100      ENDIF
1101c
1102c Initilialize module fonte_neige_mod     
1103c
1104      CALL fonte_neige_init(run_off_lic_0)
1105
1106
1107      RETURN
1108      END
Note: See TracBrowser for help on using the repository browser.