source: LMDZ6/branches/Ocean_skin/libf/phylmd/phyetat0.F90 @ 3428

Last change on this file since 3428 was 3422, checked in by oboucher, 6 years ago

Changing initialisation of co2_send if not found in restart.
  • 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.4 KB
Line 
1! $Id: phyetat0.F90 3422 2018-12-08 18:52:56Z lguez $
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 (iflag_cycle_diurne.GE.1) tab_cntrl( 7) = iflag_cycle_diurne
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  ENDIF
438
439!--OB now this is for co2i
440  IF (type_trac == 'co2i') THEN
441     IF (carbon_cycle_cpl) THEN
442        ALLOCATE(co2_send(klon), stat=ierr)
443        IF (ierr /= 0) CALL abort_physic('phyetat0', 'pb allocation co2_send', 1)
444        found=phyetat0_get(1,co2_send,"co2_send","co2 send",co2_ppm)
445     END IF
446  END IF
447
448!===========================================
449!  ondes de gravite / relief
450!===========================================
451
452!  ondes de gravite non orographiques
453  if (ok_gwd_rando) found = &
454       phyetat0_get(klev,du_gwd_rando,"du_gwd_rando","du_gwd_rando",0.)
455  IF (.not. ok_hines .and. ok_gwd_rando) found &
456       = phyetat0_get(klev,du_gwd_front,"du_gwd_front","du_gwd_front",0.)
457
458!  prise en compte du relief sous-maille
459  found=phyetat0_get(1,zmea,"ZMEA","sub grid orography",0.)
460  found=phyetat0_get(1,zstd,"ZSTD","sub grid orography",0.)
461  found=phyetat0_get(1,zsig,"ZSIG","sub grid orography",0.)
462  found=phyetat0_get(1,zgam,"ZGAM","sub grid orography",0.)
463  found=phyetat0_get(1,zthe,"ZTHE","sub grid orography",0.)
464  found=phyetat0_get(1,zpic,"ZPIC","sub grid orography",0.)
465  found=phyetat0_get(1,zval,"ZVAL","sub grid orography",0.)
466  found=phyetat0_get(1,zmea,"ZMEA","sub grid orography",0.)
467  found=phyetat0_get(1,rugoro,"RUGSREL","sub grid orography",0.)
468
469!===========================================
470! Initialize ocean
471!===========================================
472
473  IF ( type_ocean == 'slab' ) THEN
474      CALL ocean_slab_init(dtime, pctsrf)
475      IF (nslay.EQ.1) THEN
476        found=phyetat0_get(1,tslab,"tslab01","tslab",0.)
477        IF (.NOT. found) THEN
478            found=phyetat0_get(1,tslab,"tslab","tslab",0.)
479        END IF
480      ELSE
481          DO i=1,nslay
482            WRITE(str2,'(i2.2)') i
483            found=phyetat0_get(1,tslab(:,i),"tslab"//str2,"tslab",0.) 
484          END DO
485      END IF
486      IF (.NOT. found) THEN
487          PRINT*, "phyetat0: Le champ <tslab> est absent"
488          PRINT*, "Initialisation a tsol_oce"
489          DO i=1,nslay
490              tslab(:,i)=MAX(ftsol(:,is_oce),271.35)
491          END DO
492      END IF
493
494      ! Sea ice variables
495      IF (version_ocean == 'sicINT') THEN
496          found=phyetat0_get(1,tice,"slab_tice","slab_tice",0.)
497          IF (.NOT. found) THEN
498              PRINT*, "phyetat0: Le champ <tice> est absent"
499              PRINT*, "Initialisation a tsol_sic"
500                  tice(:)=ftsol(:,is_sic)
501          END IF
502          found=phyetat0_get(1,seaice,"seaice","seaice",0.)
503          IF (.NOT. found) THEN
504              PRINT*, "phyetat0: Le champ <seaice> est absent"
505              PRINT*, "Initialisation a 0/1m suivant fraction glace"
506              seaice(:)=0.
507              WHERE (pctsrf(:,is_sic).GT.EPSFRA)
508                  seaice=917.
509              END WHERE
510          END IF
511      END IF !sea ice INT
512  END IF ! Slab       
513
514  ! on ferme le fichier
515  CALL close_startphy
516
517  ! Initialize module pbl_surface_mod
518
519  CALL pbl_surface_init(fder, snow, qsurf, tsoil)
520
521  ! Initialize module ocean_cpl_mod for the case of coupled ocean
522  IF ( type_ocean == 'couple' ) THEN
523     CALL ocean_cpl_init(dtime, longitude_deg, latitude_deg)
524  ENDIF
525
526  CALL init_iophy_new(latitude_deg, longitude_deg)
527
528  ! Initilialize module fonte_neige_mod     
529  CALL fonte_neige_init(run_off_lic_0)
530
531END SUBROUTINE phyetat0
532
533!===================================================================
534FUNCTION phyetat0_get(nlev,field,name,descr,default)
535!===================================================================
536! Lecture d'un champ avec contrôle
537! Function logique dont le resultat indique si la lecture
538! s'est bien passée
539! On donne une valeur par defaut dans le cas contraire
540!===================================================================
541
542USE iostart, ONLY : get_field
543USE dimphy, only: klon
544USE print_control_mod, ONLY: lunout
545
546IMPLICIT NONE
547
548LOGICAL phyetat0_get
549
550! arguments
551INTEGER,INTENT(IN) :: nlev
552CHARACTER*(*),INTENT(IN) :: name,descr
553REAL,INTENT(IN) :: default
554REAL,DIMENSION(klon,nlev),INTENT(INOUT) :: field
555
556! Local variables
557LOGICAL found
558
559   CALL get_field(name, field, found)
560   IF (.NOT. found) THEN
561     WRITE(lunout,*) "phyetat0: Le champ <",name,"> est absent"
562     WRITE(lunout,*) "Depart legerement fausse. Mais je continue"
563     field(:,:)=default
564   ENDIF
565   WRITE(lunout,*) name, descr, MINval(field),MAXval(field)
566   phyetat0_get=found
567
568RETURN
569END FUNCTION phyetat0_get
570
571!================================================================
572FUNCTION phyetat0_srf(nlev,field,name,descr,default)
573!===================================================================
574! Lecture d'un champ par sous-surface avec contrôle
575! Function logique dont le resultat indique si la lecture
576! s'est bien passée
577! On donne une valeur par defaut dans le cas contraire
578!===================================================================
579
580USE iostart, ONLY : get_field
581USE dimphy, only: klon
582USE indice_sol_mod, only: nbsrf
583USE print_control_mod, ONLY: lunout
584
585IMPLICIT NONE
586
587LOGICAL phyetat0_srf
588! arguments
589INTEGER,INTENT(IN) :: nlev
590CHARACTER*(*),INTENT(IN) :: name,descr
591REAL,INTENT(IN) :: default
592REAL,DIMENSION(klon,nlev,nbsrf),INTENT(INOUT) :: field
593
594! Local variables
595LOGICAL found,phyetat0_get
596INTEGER nsrf
597CHARACTER*2 str2
598 
599     IF (nbsrf.GT.99) THEN
600        WRITE(lunout,*) "Trop de sous-mailles"
601        call abort_physic("phyetat0", "", 1)
602     ENDIF
603
604     DO nsrf = 1, nbsrf
605        WRITE(str2, '(i2.2)') nsrf
606        found= phyetat0_get(nlev,field(:,:, nsrf), &
607        name//str2,descr//" srf:"//str2,default)
608     ENDDO
609
610     phyetat0_srf=found
611
612RETURN
613END FUNCTION phyetat0_srf
614
Note: See TracBrowser for help on using the repository browser.