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

Last change on this file since 3628 was 3628, checked in by lguez, 4 years ago

If the ocean skin parameterization is working actively
(activate_ocean_skin == 2) and we are coupled to the ocean then send
ocean-air interface salinity to the ocean. New dummy argument s_int
of procedures ocean_cpl_noice and cpl_send_ocean_fields. We can
only send interface salinity from the previous time-step since
communication with the ocean is before the call to bulk_flux. So make
s_int a state variable: move s_int from phys_output_var_mod to
phys_state_var_mod. Still, we only read s_int from startphy,
define it before the call to surf_ocean and write it to restartphy
if activate_ocean_skin == 2 and type_ocean == 'couple'. In
procedure pbl_surface, for clarity, move the definition of output
variables t_int, dter, dser, tkt, tks, rf, taur to missing_val to
after the call to surf_ocean, with the definition of s_int,
ds_ns, dt_ns to missing_val. This does not change anything for
t_int, dter, dser, tkt, tks, rf, taur. In pbl_surface_newfrac, we
choose to set s_int to 35 for an appearing ocean point, this is
questionable. In surf_ocean, change the intent of s_int from out
to inout.

  • 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: 22.3 KB
RevLine 
[1403]1! $Id: phyetat0.F90 3628 2020-02-04 21:36:32Z lguez $
[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
[3605]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, &
[2497]15       ftsol, pbl_tke, pctsrf, q_ancien, ql_ancien, qs_ancien, radpas, radsol, rain_fall, ratqs, &
[2399]16       rnebcon, rugoro, sig1, snow_fall, solaire_etat0, sollw, sollwdown, &
[1938]17       solsw, t_ancien, u_ancien, v_ancien, w01, wake_cstar, wake_deltaq, &
[2159]18       wake_deltat, wake_delta_pbl_TKE, delta_tsurf, wake_fip, wake_pe, &
[2635]19       wake_s, wake_dens, zgam, zmax0, zmea, zpic, zsig, &
[3080]20       zstd, zthe, zval, ale_bl, ale_bl_trig, alp_bl, u10m, v10m, treedrg, &
[3628]21       ale_wake, ale_bl_stat, s_int, ds_ns, dt_ns
[2952]22!FC
[2399]23  USE geometry_mod, ONLY : longitude_deg, latitude_deg
[1938]24  USE iostart, ONLY : close_startphy, get_field, get_var, open_startphy
[2320]25  USE infotrac_phy, only: nbtr, nqo, type_trac, tname, niadv
[1827]26  USE traclmdz_mod,    ONLY : traclmdz_from_restart
27  USE carbon_cycle_mod, ONLY : carbon_cycle_tr, carbon_cycle_cpl, co2_send
[1938]28  USE indice_sol_mod, only: nbsrf, is_ter, epsfra, is_lic, is_oce, is_sic
[2656]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
[3626]31#ifdef CPP_XIOS
32  USE wxios, ONLY: missing_val
33#else
34  use netcdf, only: missing_val => nf90_fill_real
35#endif
[3458]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)
[2237]77  INTEGER it, iiq, 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)
[3605]117     RCO2    = co2_ppm * 1.0e-06 * RMCO2 / RMD
[1827]118     ! ELSE : keep value from .def
[3605]119  ENDIF
[1279]120
[3605]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
[3605]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
[3605]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
[3605]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
[3605]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
[3605]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.)
314  found=phyetat0_get(1,sollw,"sollw","net LW radiation surf",0.)
315  found=phyetat0_get(1,sollwdown,"sollwdown","down LW radiation surf",0.)
[1827]316  IF (.NOT. found) THEN
[2252]317     sollwdown = 0. ;  zts=0.
[2188]318     do nsrf=1,nbsrf
319        zts(:)=zts(:)+ftsol(:,nsrf)*pctsrf(:,nsrf)
320     enddo
321     sollwdown(:)=sollw(:)+RSIGMA*zts(:)**4
322  ENDIF
323
[2252]324  found=phyetat0_get(1,radsol,"RADS","Solar radiation",0.)
325  found=phyetat0_get(1,fder,"fder","Flux derivative",0.)
[2188]326
[1827]327
328  ! Lecture de la longueur de rugosite
[2243]329  found=phyetat0_srf(1,z0m,"RUG","Z0m ancien",0.001)
330  IF (found) THEN
331     z0h(:,1:nbsrf)=z0m(:,1:nbsrf)
[1827]332  ELSE
[2243]333     found=phyetat0_srf(1,z0m,"Z0m","Roughness length, momentum ",0.001)
334     found=phyetat0_srf(1,z0h,"Z0h","Roughness length, enthalpy ",0.001)
[1827]335  ENDIF
[2952]336!FC
[2979]337  IF (ifl_pbltree>0) THEN
[2952]338!CALL get_field("FTER", pctsrf(:, is_ter), found)
[2979]339    treedrg(:,1:klev,1:nbsrf)= 0.0
340    CALL get_field("treedrg_ter", drg_ter(:,:), found)
[2952]341!  found=phyetat0_srf(1,treedrg,"treedrg","drag from vegetation" , 0.)
[2979]342    !lecture du profile de freinage des arbres
343    IF (.not. found ) THEN
344      treedrg(:,1:klev,1:nbsrf)= 0.0
345    ELSE
346      treedrg(:,1:klev,is_ter)= drg_ter(:,:)
[2952]347!     found=phyetat0_srf(klev,treedrg,"treedrg","freinage arbres",0.)
[2979]348    ENDIF
349  ELSE
350    ! initialize treedrg(), because it will be written in restartphy.nc
351    treedrg(:,:,:) = 0.0
352  ENDIF
[1827]353
354  ! Lecture de l'age de la neige:
[2252]355  found=phyetat0_srf(1,agesno,"AGESNO","SNOW AGE",0.001)
[1827]356
[2252]357  ancien_ok=.true.
358  ancien_ok=ancien_ok.AND.phyetat0_get(klev,t_ancien,"TANCIEN","TANCIEN",0.)
359  ancien_ok=ancien_ok.AND.phyetat0_get(klev,q_ancien,"QANCIEN","QANCIEN",0.)
[2497]360  ancien_ok=ancien_ok.AND.phyetat0_get(klev,ql_ancien,"QLANCIEN","QLANCIEN",0.)
361  ancien_ok=ancien_ok.AND.phyetat0_get(klev,qs_ancien,"QSANCIEN","QSANCIEN",0.)
[2252]362  ancien_ok=ancien_ok.AND.phyetat0_get(klev,u_ancien,"UANCIEN","UANCIEN",0.)
363  ancien_ok=ancien_ok.AND.phyetat0_get(klev,v_ancien,"VANCIEN","VANCIEN",0.)
[2499]364  ancien_ok=ancien_ok.AND.phyetat0_get(1,prw_ancien,"PRWANCIEN","PRWANCIEN",0.)
365  ancien_ok=ancien_ok.AND.phyetat0_get(1,prlw_ancien,"PRLWANCIEN","PRLWANCIEN",0.)
366  ancien_ok=ancien_ok.AND.phyetat0_get(1,prsw_ancien,"PRSWANCIEN","PRSWANCIEN",0.)
[1827]367
[2494]368  ! Ehouarn: addtional tests to check if t_ancien, q_ancien contain
369  !          dummy values (as is the case when generated by ce0l,
370  !          or by iniaqua)
[3605]371  IF ( (maxval(q_ancien).EQ.minval(q_ancien))       .OR. &
372       (maxval(ql_ancien).EQ.minval(ql_ancien))     .OR. &
373       (maxval(qs_ancien).EQ.minval(qs_ancien))     .OR. &
374       (maxval(prw_ancien).EQ.minval(prw_ancien))   .OR. &
375       (maxval(prlw_ancien).EQ.minval(prlw_ancien)) .OR. &
376       (maxval(prsw_ancien).EQ.minval(prsw_ancien)) .OR. &
377       (maxval(t_ancien).EQ.minval(t_ancien)) ) THEN
[2494]378    ancien_ok=.false.
[3605]379  ENDIF
[2494]380
[2252]381  found=phyetat0_get(klev,clwcon,"CLWCON","CLWCON",0.)
382  found=phyetat0_get(klev,rnebcon,"RNEBCON","RNEBCON",0.)
383  found=phyetat0_get(klev,ratqs,"RATQS","RATQS",0.)
[1827]384
[2252]385  found=phyetat0_get(1,run_off_lic_0,"RUNOFFLIC0","RUNOFFLIC0",0.)
[1827]386
[2252]387!==================================
388!  TKE
389!==================================
390!
[1827]391  IF (iflag_pbl>1) then
[2252]392     found=phyetat0_srf(klev+1,pbl_tke,"TKE","Turb. Kinetic. Energ. ",1.e-8)
[1827]393  ENDIF
[1403]394
[2252]395  IF (iflag_pbl>1 .AND. iflag_wake>=1  .AND. iflag_pbl_split >=1 ) then
396    found=phyetat0_srf(klev+1,wake_delta_pbl_tke,"DELTATKE","Del TKE wk/env",0.)
397    found=phyetat0_srf(1,delta_tsurf,"DELTA_TSURF","Delta Ts wk/env ",0.)
[2159]398  ENDIF   !(iflag_pbl>1 .AND. iflag_wake>=1 .AND. iflag_pbl_split >=1 )
399
[2251]400!==================================
401!  thermiques, poches, convection
402!==================================
[1403]403
[2252]404! Emanuel
405  found=phyetat0_get(klev,sig1,"sig1","sig1",0.)
406  found=phyetat0_get(klev,w01,"w01","w01",0.)
[1403]407
[2252]408! Wake
[2251]409  found=phyetat0_get(klev,wake_deltat,"WAKE_DELTAT","Delta T wake/env",0.)
[2243]410  found=phyetat0_get(klev,wake_deltaq,"WAKE_DELTAQ","Delta hum. wake/env",0.)
[2635]411  found=phyetat0_get(1,wake_s,"WAKE_S","Wake frac. area",0.)
[3000]412!jyg<
413!  Set wake_dens to -1000. when there is no restart so that the actual
414!  initialization is made in calwake.
415!!  found=phyetat0_get(1,wake_dens,"WAKE_DENS","Wake num. /unit area",0.)
416  found=phyetat0_get(1,wake_dens,"WAKE_DENS","Wake num. /unit area",-1000.)
417!>jyg
[2252]418  found=phyetat0_get(1,wake_cstar,"WAKE_CSTAR","WAKE_CSTAR",0.)
419  found=phyetat0_get(1,wake_pe,"WAKE_PE","WAKE_PE",0.)
420  found=phyetat0_get(1,wake_fip,"WAKE_FIP","WAKE_FIP",0.)
[879]421
[2252]422! Thermiques
[2251]423  found=phyetat0_get(1,zmax0,"ZMAX0","ZMAX0",40.)
424  found=phyetat0_get(1,f0,"F0","F0",1.e-5)
[2252]425  found=phyetat0_get(klev+1,fm_therm,"FM_THERM","Thermals mass flux",0.)
[2251]426  found=phyetat0_get(klev,entr_therm,"ENTR_THERM","Thermals Entrain.",0.)
427  found=phyetat0_get(klev,detr_therm,"DETR_THERM","Thermals Detrain.",0.)
[782]428
[2252]429! ALE/ALP
[2251]430  found=phyetat0_get(1,ale_bl,"ALE_BL","ALE BL",0.)
431  found=phyetat0_get(1,ale_bl_trig,"ALE_BL_TRIG","ALE BL_TRIG",0.)
432  found=phyetat0_get(1,alp_bl,"ALP_BL","ALP BL",0.)
[3080]433  found=phyetat0_get(1,ale_wake,"ALE_WAKE","ALE_WAKE",0.)
434  found=phyetat0_get(1,ale_bl_stat,"ALE_BL_STAT","ALE_BL_STAT",0.)
[1279]435
[2251]436!===========================================
[1827]437  ! Read and send field trs to traclmdz
[2251]438!===========================================
[1827]439
440  IF (type_trac == 'lmdz') THEN
[2265]441     DO it=1, nbtr                                                                 
442!!        iiq=niadv(it+2)                                                           ! jyg
443        iiq=niadv(it+nqo)                                                           ! jyg
[2252]444        found=phyetat0_get(1,trs(:,it),"trs_"//tname(iiq), &
445              "Surf trac"//tname(iiq),0.)
[3605]446     ENDDO
[1827]447     CALL traclmdz_from_restart(trs)
[3421]448  ENDIF
[1827]449
[3421]450!--OB now this is for co2i
451  IF (type_trac == 'co2i') THEN
[1827]452     IF (carbon_cycle_cpl) THEN
[2252]453        ALLOCATE(co2_send(klon), stat=ierr)
[2311]454        IF (ierr /= 0) CALL abort_physic('phyetat0', 'pb allocation co2_send', 1)
[3422]455        found=phyetat0_get(1,co2_send,"co2_send","co2 send",co2_ppm)
[3605]456     ENDIF
457  ENDIF
[1827]458
[2251]459!===========================================
[2252]460!  ondes de gravite / relief
[2251]461!===========================================
462
[2252]463!  ondes de gravite non orographiques
[3605]464  IF (ok_gwd_rando) found = &
[2333]465       phyetat0_get(klev,du_gwd_rando,"du_gwd_rando","du_gwd_rando",0.)
[3605]466  IF (.NOT. ok_hines .AND. ok_gwd_rando) found &
[2333]467       = phyetat0_get(klev,du_gwd_front,"du_gwd_front","du_gwd_front",0.)
[1938]468
[2252]469!  prise en compte du relief sous-maille
470  found=phyetat0_get(1,zmea,"ZMEA","sub grid orography",0.)
471  found=phyetat0_get(1,zstd,"ZSTD","sub grid orography",0.)
472  found=phyetat0_get(1,zsig,"ZSIG","sub grid orography",0.)
473  found=phyetat0_get(1,zgam,"ZGAM","sub grid orography",0.)
474  found=phyetat0_get(1,zthe,"ZTHE","sub grid orography",0.)
475  found=phyetat0_get(1,zpic,"ZPIC","sub grid orography",0.)
476  found=phyetat0_get(1,zval,"ZVAL","sub grid orography",0.)
477  found=phyetat0_get(1,zmea,"ZMEA","sub grid orography",0.)
478  found=phyetat0_get(1,rugoro,"RUGSREL","sub grid orography",0.)
479
[2251]480!===========================================
481! Initialize ocean
482!===========================================
483
[2057]484  IF ( type_ocean == 'slab' ) THEN
[3605]485      CALL ocean_slab_init(phys_tstep, pctsrf)
[2656]486      IF (nslay.EQ.1) THEN
487        found=phyetat0_get(1,tslab,"tslab01","tslab",0.)
488        IF (.NOT. found) THEN
489            found=phyetat0_get(1,tslab,"tslab","tslab",0.)
[3605]490        ENDIF
[2656]491      ELSE
492          DO i=1,nslay
493            WRITE(str2,'(i2.2)') i
494            found=phyetat0_get(1,tslab(:,i),"tslab"//str2,"tslab",0.) 
[3605]495          ENDDO
496      ENDIF
[2057]497      IF (.NOT. found) THEN
498          PRINT*, "phyetat0: Le champ <tslab> est absent"
499          PRINT*, "Initialisation a tsol_oce"
500          DO i=1,nslay
[2209]501              tslab(:,i)=MAX(ftsol(:,is_oce),271.35)
[3605]502          ENDDO
503      ENDIF
[2251]504
[2209]505      ! Sea ice variables
506      IF (version_ocean == 'sicINT') THEN
[2656]507          found=phyetat0_get(1,tice,"slab_tice","slab_tice",0.)
[2209]508          IF (.NOT. found) THEN
509              PRINT*, "phyetat0: Le champ <tice> est absent"
510              PRINT*, "Initialisation a tsol_sic"
511                  tice(:)=ftsol(:,is_sic)
[3605]512          ENDIF
[2656]513          found=phyetat0_get(1,seaice,"seaice","seaice",0.)
[2209]514          IF (.NOT. found) THEN
515              PRINT*, "phyetat0: Le champ <seaice> est absent"
516              PRINT*, "Initialisation a 0/1m suivant fraction glace"
517              seaice(:)=0.
518              WHERE (pctsrf(:,is_sic).GT.EPSFRA)
519                  seaice=917.
[3605]520              ENDWHERE
521          ENDIF
522      ENDIF !sea ice INT
523  ENDIF ! Slab       
[2057]524
[1827]525  ! on ferme le fichier
526  CALL close_startphy
527
528  ! Initialize module pbl_surface_mod
529
[2243]530  CALL pbl_surface_init(fder, snow, qsurf, tsoil)
[1827]531
532  ! Initialize module ocean_cpl_mod for the case of coupled ocean
533  IF ( type_ocean == 'couple' ) THEN
[3605]534     CALL ocean_cpl_init(phys_tstep, longitude_deg, latitude_deg)
[1827]535  ENDIF
536
[3605]537!  CALL init_iophy_new(latitude_deg, longitude_deg)
[2054]538
[1827]539  ! Initilialize module fonte_neige_mod     
540  CALL fonte_neige_init(run_off_lic_0)
541
[3458]542  if (activate_ocean_skin >= 1) then
[3628]543     if (activate_ocean_skin == 2 .and. type_ocean == 'couple') found &
544          = phyetat0_get(1, s_int, "S_int", "salinity at air-ocean interface", &
545          35.)
546     found = phyetat0_get(1, ds_ns, "dS_ns", "delta salinity near surface", 0.)
[3458]547     found = phyetat0_get(1, dt_ns, "dT_ns", "delta temperature near surface", &
548          0.)
[3429]549
[3458]550     where (pctsrf(:, is_oce) == 0.)
[3628]551        s_int = missing_val
[3626]552        ds_ns = missing_val
553        dt_ns = missing_val
[3458]554     end where
555  end if
[3429]556
[1827]557END SUBROUTINE phyetat0
[2243]558
559!===================================================================
560FUNCTION phyetat0_get(nlev,field,name,descr,default)
561!===================================================================
562! Lecture d'un champ avec contrôle
563! Function logique dont le resultat indique si la lecture
564! s'est bien passée
565! On donne une valeur par defaut dans le cas contraire
566!===================================================================
567
568USE iostart, ONLY : get_field
569USE dimphy, only: klon
[2311]570USE print_control_mod, ONLY: lunout
[2243]571
572IMPLICIT NONE
573
574LOGICAL phyetat0_get
575
576! arguments
577INTEGER,INTENT(IN) :: nlev
578CHARACTER*(*),INTENT(IN) :: name,descr
579REAL,INTENT(IN) :: default
580REAL,DIMENSION(klon,nlev),INTENT(INOUT) :: field
581
582! Local variables
583LOGICAL found
584
585   CALL get_field(name, field, found)
586   IF (.NOT. found) THEN
587     WRITE(lunout,*) "phyetat0: Le champ <",name,"> est absent"
588     WRITE(lunout,*) "Depart legerement fausse. Mais je continue"
589     field(:,:)=default
590   ENDIF
591   WRITE(lunout,*) name, descr, MINval(field),MAXval(field)
592   phyetat0_get=found
593
594RETURN
595END FUNCTION phyetat0_get
596
597!================================================================
598FUNCTION phyetat0_srf(nlev,field,name,descr,default)
599!===================================================================
600! Lecture d'un champ par sous-surface avec contrôle
601! Function logique dont le resultat indique si la lecture
602! s'est bien passée
603! On donne une valeur par defaut dans le cas contraire
604!===================================================================
605
606USE iostart, ONLY : get_field
607USE dimphy, only: klon
608USE indice_sol_mod, only: nbsrf
[2311]609USE print_control_mod, ONLY: lunout
[2243]610
611IMPLICIT NONE
612
613LOGICAL phyetat0_srf
614! arguments
615INTEGER,INTENT(IN) :: nlev
616CHARACTER*(*),INTENT(IN) :: name,descr
617REAL,INTENT(IN) :: default
618REAL,DIMENSION(klon,nlev,nbsrf),INTENT(INOUT) :: field
619
620! Local variables
621LOGICAL found,phyetat0_get
622INTEGER nsrf
623CHARACTER*2 str2
624 
625     IF (nbsrf.GT.99) THEN
626        WRITE(lunout,*) "Trop de sous-mailles"
[2311]627        call abort_physic("phyetat0", "", 1)
[2243]628     ENDIF
629
630     DO nsrf = 1, nbsrf
631        WRITE(str2, '(i2.2)') nsrf
632        found= phyetat0_get(nlev,field(:,:, nsrf), &
633        name//str2,descr//" srf:"//str2,default)
634     ENDDO
635
636     phyetat0_srf=found
637
638RETURN
639END FUNCTION phyetat0_srf
640
Note: See TracBrowser for help on using the repository browser.