source: LMDZ5/branches/LMDZ_tree_FC/libf/phylmd/phyetat0.F90 @ 2938

Last change on this file since 2938 was 2938, checked in by fcheruy, 7 years ago

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