source: LMDZ5/branches/testing/libf/phylmd/phyetat0.F90 @ 5425

Last change on this file since 5425 was 2669, checked in by Laurent Fairhead, 8 years ago

Merged trunk changes r2640:2664 into testing branch
  • 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: 20.3 KB
Line 
1! $Id: phyetat0.F90 2669 2016-10-14 12:57:28Z jyg $
2
3SUBROUTINE phyetat0 (fichnom, clesphy0, tabcntr0)
4
5  USE dimphy, only: klon, zmasq, klev
6  USE iophy, ONLY : init_iophy_new
7  USE ocean_cpl_mod,    ONLY : ocean_cpl_init
8  USE fonte_neige_mod,  ONLY : fonte_neige_init
9  USE pbl_surface_mod,  ONLY : pbl_surface_init
10  USE surface_data,     ONLY : type_ocean, version_ocean
11  USE phys_state_var_mod, ONLY : ancien_ok, clwcon, detr_therm, dtime, &
12       qsol, fevap, z0m, z0h, agesno, &
13       du_gwd_rando, du_gwd_front, entr_therm, f0, fm_therm, &
14       falb_dir, falb_dif, prw_ancien, prlw_ancien, prsw_ancien, &
15       ftsol, pbl_tke, pctsrf, q_ancien, ql_ancien, qs_ancien, radpas, radsol, rain_fall, ratqs, &
16       rnebcon, rugoro, sig1, snow_fall, solaire_etat0, sollw, sollwdown, &
17       solsw, t_ancien, u_ancien, v_ancien, w01, wake_cstar, wake_deltaq, &
18       wake_deltat, wake_delta_pbl_TKE, delta_tsurf, wake_fip, wake_pe, &
19       wake_s, wake_dens, zgam, zmax0, zmea, zpic, zsig, &
20       zstd, zthe, zval, ale_bl, ale_bl_trig, alp_bl, u10m, v10m
21  USE geometry_mod, ONLY : longitude_deg, latitude_deg
22  USE iostart, ONLY : close_startphy, get_field, get_var, open_startphy
23  USE infotrac_phy, only: nbtr, nqo, type_trac, tname, niadv
24  USE traclmdz_mod,    ONLY : traclmdz_from_restart
25  USE carbon_cycle_mod, ONLY : carbon_cycle_tr, carbon_cycle_cpl, co2_send
26  USE indice_sol_mod, only: nbsrf, is_ter, epsfra, is_lic, is_oce, is_sic
27  USE ocean_slab_mod, ONLY: nslay, tslab, seaice, tice, ocean_slab_init
28  USE time_phylmdz_mod, ONLY: init_iteration, pdtphys, itau_phy
29
30  IMPLICIT none
31  !======================================================================
32  ! Auteur(s) Z.X. Li (LMD/CNRS) date: 19930818
33  ! Objet: Lecture de l'etat initial pour la physique
34  !======================================================================
35  include "netcdf.inc"
36  include "dimsoil.h"
37  include "clesphys.h"
38  include "thermcell.h"
39  include "compbl.h"
40  include "YOMCST.h"
41  !======================================================================
42  CHARACTER*(*) fichnom
43
44  ! les variables globales lues dans le fichier restart
45
46  REAL tsoil(klon, nsoilmx, nbsrf)
47  REAL qsurf(klon, nbsrf)
48  REAL snow(klon, nbsrf)
49  real fder(klon)
50  REAL run_off_lic_0(klon)
51  REAL fractint(klon)
52  REAL trs(klon, nbtr)
53  REAL zts(klon)
54
55  CHARACTER*6 ocean_in
56  LOGICAL ok_veget_in
57
58  INTEGER        longcles
59  PARAMETER    ( longcles = 20 )
60  REAL clesphy0( longcles )
61
62  REAL xmin, xmax
63
64  INTEGER nid, nvarid
65  INTEGER ierr, i, nsrf, isoil , k
66  INTEGER length
67  PARAMETER (length=100)
68  INTEGER it, iiq, isw
69  REAL tab_cntrl(length), tabcntr0(length)
70  CHARACTER*7 str7
71  CHARACTER*2 str2
72  LOGICAL :: found,phyetat0_get,phyetat0_srf
73  REAL :: lon_startphy(klon), lat_startphy(klon)
74
75  ! FH1D
76  !     real iolat(jjm+1)
77  !real iolat(jjm+1-1/(iim*jjm))
78
79  ! Ouvrir le fichier contenant l'etat initial:
80
81  CALL open_startphy(fichnom)
82
83  ! Lecture des parametres de controle:
84
85  CALL get_var("controle", tab_cntrl)
86
87!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
88  ! FH 2008/05/09 On elimine toutes les clefs physiques dans la dynamique
89  ! Les constantes de la physiques sont lues dans la physique seulement.
90  ! Les egalites du type
91  !             tab_cntrl( 5 )=clesphy0(1)
92  ! sont remplacees par
93  !             clesphy0(1)=tab_cntrl( 5 )
94  ! On inverse aussi la logique.
95  ! On remplit les tab_cntrl avec les parametres lus dans les .def
96!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
97
98  DO i = 1, length
99     tabcntr0( i ) = tab_cntrl( i )
100  ENDDO
101
102  tab_cntrl(1)=pdtphys
103  tab_cntrl(2)=radpas
104
105  ! co2_ppm : value from the previous time step
106  IF (carbon_cycle_tr .OR. carbon_cycle_cpl) THEN
107     co2_ppm = tab_cntrl(3)
108     RCO2    = co2_ppm * 1.0e-06  * 44.011/28.97
109     ! ELSE : keep value from .def
110  END IF
111
112  ! co2_ppm0 : initial value of atmospheric CO2 (from create_etat0_limit.e .def)
113  co2_ppm0   = tab_cntrl(16)
114
115  solaire_etat0      = tab_cntrl(4)
116  tab_cntrl(5)=iflag_con
117  tab_cntrl(6)=nbapp_rad
118
119  if (cycle_diurne) tab_cntrl( 7) =1.
120  if (soil_model) tab_cntrl( 8) =1.
121  if (new_oliq) tab_cntrl( 9) =1.
122  if (ok_orodr) tab_cntrl(10) =1.
123  if (ok_orolf) tab_cntrl(11) =1.
124  if (ok_limitvrai) tab_cntrl(12) =1.
125
126  itau_phy = tab_cntrl(15)
127
128  clesphy0(1)=tab_cntrl( 5 )
129  clesphy0(2)=tab_cntrl( 6 )
130  clesphy0(3)=tab_cntrl( 7 )
131  clesphy0(4)=tab_cntrl( 8 )
132  clesphy0(5)=tab_cntrl( 9 )
133  clesphy0(6)=tab_cntrl( 10 )
134  clesphy0(7)=tab_cntrl( 11 )
135  clesphy0(8)=tab_cntrl( 12 )
136
137  ! set time iteration
138   CALL init_iteration(itau_phy)
139
140  ! read latitudes and make a sanity check (because already known from dyn)
141  CALL get_field("latitude",lat_startphy)
142  DO i=1,klon
143    IF (ABS(lat_startphy(i)-latitude_deg(i))>=1) THEN
144      WRITE(*,*) "phyetat0: Error! Latitude discrepancy wrt startphy file:",&
145                 " i=",i," lat_startphy(i)=",lat_startphy(i),&
146                 " latitude_deg(i)=",latitude_deg(i)
147      ! This is presumably serious enough to abort run
148      CALL abort_physic("phyetat0","discrepancy in latitudes!",1)
149    ENDIF
150    IF (ABS(lat_startphy(i)-latitude_deg(i))>=0.0001) THEN
151      WRITE(*,*) "phyetat0: Warning! Latitude discrepancy wrt startphy file:",&
152                 " i=",i," lat_startphy(i)=",lat_startphy(i),&
153                 " latitude_deg(i)=",latitude_deg(i)
154    ENDIF
155  ENDDO
156
157  ! read longitudes and make a sanity check (because already known from dyn)
158  CALL get_field("longitude",lon_startphy)
159  DO i=1,klon
160    IF (ABS(lon_startphy(i)-longitude_deg(i))>=1) THEN
161      WRITE(*,*) "phyetat0: Error! Longitude discrepancy wrt startphy file:",&
162                 " i=",i," lon_startphy(i)=",lon_startphy(i),&
163                 " longitude_deg(i)=",longitude_deg(i)
164      ! This is presumably serious enough to abort run
165      CALL abort_physic("phyetat0","discrepancy in longitudes!",1)
166    ENDIF
167    IF (ABS(lon_startphy(i)-longitude_deg(i))>=0.0001) THEN
168      WRITE(*,*) "phyetat0: Warning! Longitude discrepancy wrt startphy file:",&
169                 " i=",i," lon_startphy(i)=",lon_startphy(i),&
170                 " longitude_deg(i)=",longitude_deg(i)
171    ENDIF
172  ENDDO
173
174  ! Lecture du masque terre mer
175
176  CALL get_field("masque", zmasq, found)
177  IF (.NOT. found) THEN
178     PRINT*, 'phyetat0: Le champ <masque> est absent'
179     PRINT *, 'fichier startphy non compatible avec phyetat0'
180  ENDIF
181
182  ! Lecture des fractions pour chaque sous-surface
183
184  ! initialisation des sous-surfaces
185
186  pctsrf = 0.
187
188  ! fraction de terre
189
190  CALL get_field("FTER", pctsrf(:, is_ter), found)
191  IF (.NOT. found) PRINT*, 'phyetat0: Le champ <FTER> est absent'
192
193  ! fraction de glace de terre
194
195  CALL get_field("FLIC", pctsrf(:, is_lic), found)
196  IF (.NOT. found) PRINT*, 'phyetat0: Le champ <FLIC> est absent'
197
198  ! fraction d'ocean
199
200  CALL get_field("FOCE", pctsrf(:, is_oce), found)
201  IF (.NOT. found) PRINT*, 'phyetat0: Le champ <FOCE> est absent'
202
203  ! fraction glace de mer
204
205  CALL get_field("FSIC", pctsrf(:, is_sic), found)
206  IF (.NOT. found) PRINT*, 'phyetat0: Le champ <FSIC> est absent'
207
208  !  Verification de l'adequation entre le masque et les sous-surfaces
209
210  fractint( 1 : klon) = pctsrf(1 : klon, is_ter)  &
211       + pctsrf(1 : klon, is_lic)
212  DO i = 1 , klon
213     IF ( abs(fractint(i) - zmasq(i) ) .GT. EPSFRA ) THEN
214        WRITE(*, *) 'phyetat0: attention fraction terre pas ',  &
215             'coherente ', i, zmasq(i), pctsrf(i, is_ter) &
216             , pctsrf(i, is_lic)
217        WRITE(*, *) 'Je force la coherence zmasq=fractint'
218        zmasq(i) = fractint(i)
219     ENDIF
220  END DO
221  fractint (1 : klon) =  pctsrf(1 : klon, is_oce)  &
222       + pctsrf(1 : klon, is_sic)
223  DO i = 1 , klon
224     IF ( abs( fractint(i) - (1. - zmasq(i))) .GT. EPSFRA ) THEN
225        WRITE(*, *) 'phyetat0 attention fraction ocean pas ',  &
226             'coherente ', i, zmasq(i) , pctsrf(i, is_oce) &
227             , pctsrf(i, is_sic)
228        WRITE(*, *) 'Je force la coherence zmasq=1.-fractint'
229        zmasq(i) = 1. - fractint(i)
230     ENDIF
231  END DO
232
233!===================================================================
234! Lecture des temperatures du sol:
235!===================================================================
236
237  found=phyetat0_get(1,ftsol(:,1),"TS","Surface temperature",283.)
238  IF (found) THEN
239     DO nsrf=2,nbsrf
240        ftsol(:,nsrf)=ftsol(:,1)
241     ENDDO
242  ELSE
243     found=phyetat0_srf(1,ftsol,"TS","Surface temperature",283.)
244  ENDIF
245
246!===================================================================
247  ! Lecture des albedo difus et direct
248!===================================================================
249
250  DO nsrf = 1, nbsrf
251     DO isw=1, nsw
252        IF (isw.GT.99) THEN
253           PRINT*, "Trop de bandes SW"
254           call abort_physic("phyetat0", "", 1)
255        ENDIF
256        WRITE(str2, '(i2.2)') isw
257        found=phyetat0_srf(1,falb_dir(:, isw,:),"A_dir_SW"//str2//"srf","Direct Albedo",0.2)
258        found=phyetat0_srf(1,falb_dif(:, isw,:),"A_dif_SW"//str2//"srf","Direct Albedo",0.2)
259     ENDDO
260  ENDDO
261
262  found=phyetat0_srf(1,u10m,"U10M","u a 10m",0.)
263  found=phyetat0_srf(1,v10m,"V10M","v a 10m",0.)
264
265!===================================================================
266  ! Lecture des temperatures du sol profond:
267!===================================================================
268
269   DO isoil=1, nsoilmx
270        IF (isoil.GT.99) THEN
271           PRINT*, "Trop de couches "
272           call abort_physic("phyetat0", "", 1)
273        ENDIF
274        WRITE(str2,'(i2.2)') isoil
275        found=phyetat0_srf(1,tsoil(:, isoil,:),"Tsoil"//str2//"srf","Temp soil",0.)
276        IF (.NOT. found) THEN
277           PRINT*, "phyetat0: Le champ <Tsoil"//str7//"> est absent"
278           PRINT*, "          Il prend donc la valeur de surface"
279           tsoil(:, isoil, :)=ftsol(:, :)
280        ENDIF
281   ENDDO
282
283!=======================================================================
284! Lecture precipitation/evaporation
285!=======================================================================
286
287  found=phyetat0_srf(1,qsurf,"QS","Near surface hmidity",0.)
288  found=phyetat0_get(1,qsol,"QSOL","Surface hmidity / bucket",0.)
289  found=phyetat0_srf(1,snow,"SNOW","Surface snow",0.)
290  found=phyetat0_srf(1,fevap,"EVAP","evaporation",0.)
291  found=phyetat0_get(1,snow_fall,"snow_f","snow fall",0.)
292  found=phyetat0_get(1,rain_fall,"rain_f","rain fall",0.)
293
294!=======================================================================
295! Radiation
296!=======================================================================
297
298  found=phyetat0_get(1,solsw,"solsw","net SW radiation surf",0.)
299  found=phyetat0_get(1,sollw,"sollw","net LW radiation surf",0.)
300  found=phyetat0_get(1,sollwdown,"sollwdown","down LW radiation surf",0.)
301  IF (.NOT. found) THEN
302     sollwdown = 0. ;  zts=0.
303     do nsrf=1,nbsrf
304        zts(:)=zts(:)+ftsol(:,nsrf)*pctsrf(:,nsrf)
305     enddo
306     sollwdown(:)=sollw(:)+RSIGMA*zts(:)**4
307  ENDIF
308
309  found=phyetat0_get(1,radsol,"RADS","Solar radiation",0.)
310  found=phyetat0_get(1,fder,"fder","Flux derivative",0.)
311
312
313  ! Lecture de la longueur de rugosite
314  found=phyetat0_srf(1,z0m,"RUG","Z0m ancien",0.001)
315  IF (found) THEN
316     z0h(:,1:nbsrf)=z0m(:,1:nbsrf)
317  ELSE
318     found=phyetat0_srf(1,z0m,"Z0m","Roughness length, momentum ",0.001)
319     found=phyetat0_srf(1,z0h,"Z0h","Roughness length, enthalpy ",0.001)
320  ENDIF
321
322  ! Lecture de l'age de la neige:
323  found=phyetat0_srf(1,agesno,"AGESNO","SNOW AGE",0.001)
324
325  ancien_ok=.true.
326  ancien_ok=ancien_ok.AND.phyetat0_get(klev,t_ancien,"TANCIEN","TANCIEN",0.)
327  ancien_ok=ancien_ok.AND.phyetat0_get(klev,q_ancien,"QANCIEN","QANCIEN",0.)
328  ancien_ok=ancien_ok.AND.phyetat0_get(klev,ql_ancien,"QLANCIEN","QLANCIEN",0.)
329  ancien_ok=ancien_ok.AND.phyetat0_get(klev,qs_ancien,"QSANCIEN","QSANCIEN",0.)
330  ancien_ok=ancien_ok.AND.phyetat0_get(klev,u_ancien,"UANCIEN","UANCIEN",0.)
331  ancien_ok=ancien_ok.AND.phyetat0_get(klev,v_ancien,"VANCIEN","VANCIEN",0.)
332  ancien_ok=ancien_ok.AND.phyetat0_get(1,prw_ancien,"PRWANCIEN","PRWANCIEN",0.)
333  ancien_ok=ancien_ok.AND.phyetat0_get(1,prlw_ancien,"PRLWANCIEN","PRLWANCIEN",0.)
334  ancien_ok=ancien_ok.AND.phyetat0_get(1,prsw_ancien,"PRSWANCIEN","PRSWANCIEN",0.)
335
336  ! Ehouarn: addtional tests to check if t_ancien, q_ancien contain
337  !          dummy values (as is the case when generated by ce0l,
338  !          or by iniaqua)
339  if ( (maxval(q_ancien).eq.minval(q_ancien))       .or. &
340       (maxval(ql_ancien).eq.minval(ql_ancien))     .or. &
341       (maxval(qs_ancien).eq.minval(qs_ancien))     .or. &
342       (maxval(prw_ancien).eq.minval(prw_ancien))   .or. &
343       (maxval(prlw_ancien).eq.minval(prlw_ancien)) .or. &
344       (maxval(prsw_ancien).eq.minval(prsw_ancien)) .or. &
345       (maxval(t_ancien).eq.minval(t_ancien)) ) then
346    ancien_ok=.false.
347  endif
348
349  found=phyetat0_get(klev,clwcon,"CLWCON","CLWCON",0.)
350  found=phyetat0_get(klev,rnebcon,"RNEBCON","RNEBCON",0.)
351  found=phyetat0_get(klev,ratqs,"RATQS","RATQS",0.)
352
353
354  found=phyetat0_get(1,run_off_lic_0,"RUNOFFLIC0","RUNOFFLIC0",0.)
355
356!==================================
357!  TKE
358!==================================
359!
360  IF (iflag_pbl>1) then
361     found=phyetat0_srf(klev+1,pbl_tke,"TKE","Turb. Kinetic. Energ. ",1.e-8)
362  ENDIF
363
364  IF (iflag_pbl>1 .AND. iflag_wake>=1  .AND. iflag_pbl_split >=1 ) then
365    found=phyetat0_srf(klev+1,wake_delta_pbl_tke,"DELTATKE","Del TKE wk/env",0.)
366    found=phyetat0_srf(1,delta_tsurf,"DELTA_TSURF","Delta Ts wk/env ",0.)
367  ENDIF   !(iflag_pbl>1 .AND. iflag_wake>=1 .AND. iflag_pbl_split >=1 )
368
369!==================================
370!  thermiques, poches, convection
371!==================================
372
373! Emanuel
374  found=phyetat0_get(klev,sig1,"sig1","sig1",0.)
375  found=phyetat0_get(klev,w01,"w01","w01",0.)
376
377! Wake
378  found=phyetat0_get(klev,wake_deltat,"WAKE_DELTAT","Delta T wake/env",0.)
379  found=phyetat0_get(klev,wake_deltaq,"WAKE_DELTAQ","Delta hum. wake/env",0.)
380  found=phyetat0_get(1,wake_s,"WAKE_S","Wake frac. area",0.)
381  found=phyetat0_get(1,wake_dens,"WAKE_DENS","Wake num. /unit area",0.)
382  found=phyetat0_get(1,wake_cstar,"WAKE_CSTAR","WAKE_CSTAR",0.)
383  found=phyetat0_get(1,wake_pe,"WAKE_PE","WAKE_PE",0.)
384  found=phyetat0_get(1,wake_fip,"WAKE_FIP","WAKE_FIP",0.)
385
386! Thermiques
387  found=phyetat0_get(1,zmax0,"ZMAX0","ZMAX0",40.)
388  found=phyetat0_get(1,f0,"F0","F0",1.e-5)
389  found=phyetat0_get(klev+1,fm_therm,"FM_THERM","Thermals mass flux",0.)
390  found=phyetat0_get(klev,entr_therm,"ENTR_THERM","Thermals Entrain.",0.)
391  found=phyetat0_get(klev,detr_therm,"DETR_THERM","Thermals Detrain.",0.)
392
393! ALE/ALP
394  found=phyetat0_get(1,ale_bl,"ALE_BL","ALE BL",0.)
395  found=phyetat0_get(1,ale_bl_trig,"ALE_BL_TRIG","ALE BL_TRIG",0.)
396  found=phyetat0_get(1,alp_bl,"ALP_BL","ALP BL",0.)
397
398!===========================================
399  ! Read and send field trs to traclmdz
400!===========================================
401
402  IF (type_trac == 'lmdz') THEN
403     DO it=1, nbtr                                                                 
404!!        iiq=niadv(it+2)                                                           ! jyg
405        iiq=niadv(it+nqo)                                                           ! jyg
406        found=phyetat0_get(1,trs(:,it),"trs_"//tname(iiq), &
407              "Surf trac"//tname(iiq),0.)
408     END DO
409     CALL traclmdz_from_restart(trs)
410
411     IF (carbon_cycle_cpl) THEN
412        ALLOCATE(co2_send(klon), stat=ierr)
413        IF (ierr /= 0) CALL abort_physic('phyetat0', 'pb allocation co2_send', 1)
414        found=phyetat0_get(1,co2_send,"co2_send","co2 send",0.)
415     END IF
416  END IF
417
418!===========================================
419!  ondes de gravite / relief
420!===========================================
421
422!  ondes de gravite non orographiques
423  if (ok_gwd_rando) found = &
424       phyetat0_get(klev,du_gwd_rando,"du_gwd_rando","du_gwd_rando",0.)
425  IF (.not. ok_hines .and. ok_gwd_rando) found &
426       = phyetat0_get(klev,du_gwd_front,"du_gwd_front","du_gwd_front",0.)
427
428!  prise en compte du relief sous-maille
429  found=phyetat0_get(1,zmea,"ZMEA","sub grid orography",0.)
430  found=phyetat0_get(1,zstd,"ZSTD","sub grid orography",0.)
431  found=phyetat0_get(1,zsig,"ZSIG","sub grid orography",0.)
432  found=phyetat0_get(1,zgam,"ZGAM","sub grid orography",0.)
433  found=phyetat0_get(1,zthe,"ZTHE","sub grid orography",0.)
434  found=phyetat0_get(1,zpic,"ZPIC","sub grid orography",0.)
435  found=phyetat0_get(1,zval,"ZVAL","sub grid orography",0.)
436  found=phyetat0_get(1,zmea,"ZMEA","sub grid orography",0.)
437  found=phyetat0_get(1,rugoro,"RUGSREL","sub grid orography",0.)
438
439!===========================================
440! Initialize ocean
441!===========================================
442
443  IF ( type_ocean == 'slab' ) THEN
444      CALL ocean_slab_init(dtime, pctsrf)
445      IF (nslay.EQ.1) THEN
446        found=phyetat0_get(1,tslab,"tslab01","tslab",0.)
447        IF (.NOT. found) THEN
448            found=phyetat0_get(1,tslab,"tslab","tslab",0.)
449        END IF
450      ELSE
451          DO i=1,nslay
452            WRITE(str2,'(i2.2)') i
453            found=phyetat0_get(1,tslab(:,i),"tslab"//str2,"tslab",0.) 
454          END DO
455      END IF
456      IF (.NOT. found) THEN
457          PRINT*, "phyetat0: Le champ <tslab> est absent"
458          PRINT*, "Initialisation a tsol_oce"
459          DO i=1,nslay
460              tslab(:,i)=MAX(ftsol(:,is_oce),271.35)
461          END DO
462      END IF
463
464      ! Sea ice variables
465      IF (version_ocean == 'sicINT') THEN
466          found=phyetat0_get(1,tice,"slab_tice","slab_tice",0.)
467          IF (.NOT. found) THEN
468              PRINT*, "phyetat0: Le champ <tice> est absent"
469              PRINT*, "Initialisation a tsol_sic"
470                  tice(:)=ftsol(:,is_sic)
471          END IF
472          found=phyetat0_get(1,seaice,"seaice","seaice",0.)
473          IF (.NOT. found) THEN
474              PRINT*, "phyetat0: Le champ <seaice> est absent"
475              PRINT*, "Initialisation a 0/1m suivant fraction glace"
476              seaice(:)=0.
477              WHERE (pctsrf(:,is_sic).GT.EPSFRA)
478                  seaice=917.
479              END WHERE
480          END IF
481      END IF !sea ice INT
482  END IF ! Slab       
483
484  ! on ferme le fichier
485  CALL close_startphy
486
487  ! Initialize module pbl_surface_mod
488
489  CALL pbl_surface_init(fder, snow, qsurf, tsoil)
490
491  ! Initialize module ocean_cpl_mod for the case of coupled ocean
492  IF ( type_ocean == 'couple' ) THEN
493     CALL ocean_cpl_init(dtime, longitude_deg, latitude_deg)
494  ENDIF
495
496  CALL init_iophy_new(latitude_deg, longitude_deg)
497
498  ! Initilialize module fonte_neige_mod     
499  CALL fonte_neige_init(run_off_lic_0)
500
501END SUBROUTINE phyetat0
502
503!===================================================================
504FUNCTION phyetat0_get(nlev,field,name,descr,default)
505!===================================================================
506! Lecture d'un champ avec contrôle
507! Function logique dont le resultat indique si la lecture
508! s'est bien passée
509! On donne une valeur par defaut dans le cas contraire
510!===================================================================
511
512USE iostart, ONLY : get_field
513USE dimphy, only: klon
514USE print_control_mod, ONLY: lunout
515
516IMPLICIT NONE
517
518LOGICAL phyetat0_get
519
520! arguments
521INTEGER,INTENT(IN) :: nlev
522CHARACTER*(*),INTENT(IN) :: name,descr
523REAL,INTENT(IN) :: default
524REAL,DIMENSION(klon,nlev),INTENT(INOUT) :: field
525
526! Local variables
527LOGICAL found
528
529   CALL get_field(name, field, found)
530   IF (.NOT. found) THEN
531     WRITE(lunout,*) "phyetat0: Le champ <",name,"> est absent"
532     WRITE(lunout,*) "Depart legerement fausse. Mais je continue"
533     field(:,:)=default
534   ENDIF
535   WRITE(lunout,*) name, descr, MINval(field),MAXval(field)
536   phyetat0_get=found
537
538RETURN
539END FUNCTION phyetat0_get
540
541!================================================================
542FUNCTION phyetat0_srf(nlev,field,name,descr,default)
543!===================================================================
544! Lecture d'un champ par sous-surface avec contrôle
545! Function logique dont le resultat indique si la lecture
546! s'est bien passée
547! On donne une valeur par defaut dans le cas contraire
548!===================================================================
549
550USE iostart, ONLY : get_field
551USE dimphy, only: klon
552USE indice_sol_mod, only: nbsrf
553USE print_control_mod, ONLY: lunout
554
555IMPLICIT NONE
556
557LOGICAL phyetat0_srf
558! arguments
559INTEGER,INTENT(IN) :: nlev
560CHARACTER*(*),INTENT(IN) :: name,descr
561REAL,INTENT(IN) :: default
562REAL,DIMENSION(klon,nlev,nbsrf),INTENT(INOUT) :: field
563
564! Local variables
565LOGICAL found,phyetat0_get
566INTEGER nsrf
567CHARACTER*2 str2
568 
569     IF (nbsrf.GT.99) THEN
570        WRITE(lunout,*) "Trop de sous-mailles"
571        call abort_physic("phyetat0", "", 1)
572     ENDIF
573
574     DO nsrf = 1, nbsrf
575        WRITE(str2, '(i2.2)') nsrf
576        found= phyetat0_get(nlev,field(:,:, nsrf), &
577        name//str2,descr//" srf:"//str2,default)
578     ENDDO
579
580     phyetat0_srf=found
581
582RETURN
583END FUNCTION phyetat0_srf
584
Note: See TracBrowser for help on using the repository browser.