source: LMDZ6/trunk/libf/phylmd/phyetat0.F90 @ 3766

Last change on this file since 3766 was 3756, checked in by oboucher, 4 years ago

The fraction of SW down radiation that is diffuse is extracted from radlwsw.F90 and passed on to pbl_surface for inclusion in fields_out for transfer to ORCHIDEE if it requested in input file coupling_fields.def. Hence there is nothing automatic here. The fraction of SW down radiation is also added as a diagnostic. To satisfy the case when this new quantity is used in ORCHIDEE, then it is added in the restart (and read in case it is there).

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