source: LMDZ5/branches/testing/libf/phylmd/phyetat0.F @ 1707

Last change on this file since 1707 was 1707, checked in by Laurent Fairhead, 11 years ago

Version testing basée sur la r1706


Testing release based on r1706

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