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

Last change on this file since 3511 was 3505, checked in by yann meurdesoif, 5 years ago

Solve some wrong discrepency problem when comparing longitude from a restartphy file. The current discrepency test is not detecting that 360°==0°, so in some case it may stop the run for a wrong reason.

YM

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