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

Last change on this file since 4059 was 4059, checked in by oboucher, 2 years ago

Audran Borella's parametrisation for ice supersaturation
activated with flag_ice_sursat (FALSE by default)

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