source: LMDZ6/branches/IPSLCM6.0.13/libf/phylmd/phyetat0.F90 @ 5091

Last change on this file since 5091 was 3081, checked in by Laurent Fairhead, 7 years ago

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