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

Last change on this file since 3605 was 3605, checked in by lguez, 5 years ago

Merge revisions 3427:3600 of trunk into branch Ocean_skin

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