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

Last change on this file since 4352 was 4298, checked in by pcadule, 2 years ago

modifications to ensure restartability of the model when CO2 tracer is passed to radiative code, and to add diagnostics variables

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