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

Last change on this file since 3605 was 3605, checked in by lguez, 4 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
Line 
1! $Id: phyetat0.F90 3605 2019-11-21 15:43:45Z lguez $
2
3SUBROUTINE phyetat0 (fichnom, clesphy0, tabcntr0)
4
5  USE dimphy, only: klon, zmasq, klev
6  USE iophy, ONLY : init_iophy_new
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
10  USE surface_data,     ONLY : type_ocean, version_ocean
11  USE phys_state_var_mod, ONLY : ancien_ok, clwcon, detr_therm, phys_tstep, &
12       qsol, fevap, z0m, z0h, agesno, &
13       du_gwd_rando, du_gwd_front, entr_therm, f0, fm_therm, &
14       falb_dir, falb_dif, prw_ancien, prlw_ancien, prsw_ancien, &
15       ftsol, pbl_tke, pctsrf, q_ancien, ql_ancien, qs_ancien, radpas, radsol, rain_fall, ratqs, &
16       rnebcon, rugoro, sig1, snow_fall, solaire_etat0, sollw, sollwdown, &
17       solsw, t_ancien, u_ancien, v_ancien, w01, wake_cstar, wake_deltaq, &
18       wake_deltat, wake_delta_pbl_TKE, delta_tsurf, wake_fip, wake_pe, &
19       wake_s, wake_dens, zgam, zmax0, zmea, zpic, zsig, &
20       zstd, zthe, zval, ale_bl, ale_bl_trig, alp_bl, u10m, v10m, treedrg, &
21       ale_wake, ale_bl_stat, ds_ns, dt_ns
22!FC
23  USE geometry_mod, ONLY : longitude_deg, latitude_deg
24  USE iostart, ONLY : close_startphy, get_field, get_var, open_startphy
25  USE infotrac_phy, only: nbtr, nqo, type_trac, tname, niadv
26  USE traclmdz_mod,    ONLY : traclmdz_from_restart
27  USE carbon_cycle_mod, ONLY : carbon_cycle_tr, carbon_cycle_cpl, co2_send
28  USE indice_sol_mod, only: nbsrf, is_ter, epsfra, is_lic, is_oce, is_sic
29  USE ocean_slab_mod, ONLY: nslay, tslab, seaice, tice, ocean_slab_init
30  USE time_phylmdz_mod, ONLY: init_iteration, pdtphys, itau_phy
31  use netcdf, only: nf90_fill_real
32  use config_ocean_skin_m, only: activate_ocean_skin
33
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"
43  include "YOMCST.h"
44  !======================================================================
45  CHARACTER*(*) fichnom
46
47  ! les variables globales lues dans le fichier restart
48
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)
56  REAL zts(klon)
57  ! pour drag arbres FC
58  REAL drg_ter(klon,klev)
59
60  CHARACTER*6 ocean_in
61  LOGICAL ok_veget_in
62
63  INTEGER        longcles
64  PARAMETER    ( longcles = 20 )
65  REAL clesphy0( longcles )
66
67  REAL xmin, xmax
68
69  INTEGER nid, nvarid
70  INTEGER ierr, i, nsrf, isoil , k
71  INTEGER length
72  PARAMETER (length=100)
73  INTEGER it, iiq, isw
74  REAL tab_cntrl(length), tabcntr0(length)
75  CHARACTER*7 str7
76  CHARACTER*2 str2
77  LOGICAL :: found,phyetat0_get,phyetat0_srf
78  REAL :: lon_startphy(klon), lat_startphy(klon)
79
80  ! FH1D
81  !     real iolat(jjm+1)
82  !real iolat(jjm+1-1/(iim*jjm))
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
92!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
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
101!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
102
103  DO i = 1, length
104     tabcntr0( i ) = tab_cntrl( i )
105  ENDDO
106
107  tab_cntrl(1)=pdtphys
108  tab_cntrl(2)=radpas
109
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)
113     RCO2    = co2_ppm * 1.0e-06 * RMCO2 / RMD
114     ! ELSE : keep value from .def
115  ENDIF
116
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
121
122  solaire_etat0      = tab_cntrl(4)
123  tab_cntrl(5)=iflag_con
124  tab_cntrl(6)=nbapp_rad
125
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.
132
133  itau_phy = tab_cntrl(15)
134
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 )
143
144  ! set time iteration
145   CALL init_iteration(itau_phy)
146
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
163
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
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
175    ENDIF
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
182    ENDIF
183  ENDDO
184
185  ! Lecture du masque terre mer
186
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
192
193  ! Lecture des fractions pour chaque sous-surface
194
195  ! initialisation des sous-surfaces
196
197  pctsrf = 0.
198
199  ! fraction de terre
200
201  CALL get_field("FTER", pctsrf(:, is_ter), found)
202  IF (.NOT. found) PRINT*, 'phyetat0: Le champ <FTER> est absent'
203
204  ! fraction de glace de terre
205
206  CALL get_field("FLIC", pctsrf(:, is_lic), found)
207  IF (.NOT. found) PRINT*, 'phyetat0: Le champ <FLIC> est absent'
208
209  ! fraction d'ocean
210
211  CALL get_field("FOCE", pctsrf(:, is_oce), found)
212  IF (.NOT. found) PRINT*, 'phyetat0: Le champ <FOCE> est absent'
213
214  ! fraction glace de mer
215
216  CALL get_field("FSIC", pctsrf(:, is_sic), found)
217  IF (.NOT. found) PRINT*, 'phyetat0: Le champ <FSIC> est absent'
218
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
231  ENDDO
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)
239        WRITE(*, *) 'Je force la coherence zmasq=1.-fractint'
240        zmasq(i) = 1. - fractint(i)
241     ENDIF
242  ENDDO
243
244!===================================================================
245! Lecture des temperatures du sol:
246!===================================================================
247
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)
252     ENDDO
253  ELSE
254     found=phyetat0_srf(1,ftsol,"TS","Surface temperature",283.)
255  ENDIF
256
257!===================================================================
258  ! Lecture des albedo difus et direct
259!===================================================================
260
261  DO nsrf = 1, nbsrf
262     DO isw=1, nsw
263        IF (isw.GT.99) THEN
264           PRINT*, "Trop de bandes SW"
265           call abort_physic("phyetat0", "", 1)
266        ENDIF
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)
270     ENDDO
271  ENDDO
272
273  found=phyetat0_srf(1,u10m,"U10M","u a 10m",0.)
274  found=phyetat0_srf(1,v10m,"V10M","v a 10m",0.)
275
276!===================================================================
277  ! Lecture des temperatures du sol profond:
278!===================================================================
279
280   DO isoil=1, nsoilmx
281        IF (isoil.GT.99) THEN
282           PRINT*, "Trop de couches "
283           call abort_physic("phyetat0", "", 1)
284        ENDIF
285        WRITE(str2,'(i2.2)') isoil
286        found=phyetat0_srf(1,tsoil(:, isoil,:),"Tsoil"//str2//"srf","Temp soil",0.)
287        IF (.NOT. found) THEN
288           PRINT*, "phyetat0: Le champ <Tsoil"//str7//"> est absent"
289           PRINT*, "          Il prend donc la valeur de surface"
290           tsoil(:, isoil, :)=ftsol(:, :)
291        ENDIF
292   ENDDO
293
294!=======================================================================
295! Lecture precipitation/evaporation
296!=======================================================================
297
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.)
304
305!=======================================================================
306! Radiation
307!=======================================================================
308
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.)
312  IF (.NOT. found) THEN
313     sollwdown = 0. ;  zts=0.
314     do nsrf=1,nbsrf
315        zts(:)=zts(:)+ftsol(:,nsrf)*pctsrf(:,nsrf)
316     enddo
317     sollwdown(:)=sollw(:)+RSIGMA*zts(:)**4
318  ENDIF
319
320  found=phyetat0_get(1,radsol,"RADS","Solar radiation",0.)
321  found=phyetat0_get(1,fder,"fder","Flux derivative",0.)
322
323
324  ! Lecture de la longueur de rugosite
325  found=phyetat0_srf(1,z0m,"RUG","Z0m ancien",0.001)
326  IF (found) THEN
327     z0h(:,1:nbsrf)=z0m(:,1:nbsrf)
328  ELSE
329     found=phyetat0_srf(1,z0m,"Z0m","Roughness length, momentum ",0.001)
330     found=phyetat0_srf(1,z0h,"Z0h","Roughness length, enthalpy ",0.001)
331  ENDIF
332!FC
333  IF (ifl_pbltree>0) THEN
334!CALL get_field("FTER", pctsrf(:, is_ter), found)
335    treedrg(:,1:klev,1:nbsrf)= 0.0
336    CALL get_field("treedrg_ter", drg_ter(:,:), found)
337!  found=phyetat0_srf(1,treedrg,"treedrg","drag from vegetation" , 0.)
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(:,:)
343!     found=phyetat0_srf(klev,treedrg,"treedrg","freinage arbres",0.)
344    ENDIF
345  ELSE
346    ! initialize treedrg(), because it will be written in restartphy.nc
347    treedrg(:,:,:) = 0.0
348  ENDIF
349
350  ! Lecture de l'age de la neige:
351  found=phyetat0_srf(1,agesno,"AGESNO","SNOW AGE",0.001)
352
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.)
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.)
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.)
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.)
363
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)
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
374    ancien_ok=.false.
375  ENDIF
376
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.)
380
381  found=phyetat0_get(1,run_off_lic_0,"RUNOFFLIC0","RUNOFFLIC0",0.)
382
383!==================================
384!  TKE
385!==================================
386!
387  IF (iflag_pbl>1) then
388     found=phyetat0_srf(klev+1,pbl_tke,"TKE","Turb. Kinetic. Energ. ",1.e-8)
389  ENDIF
390
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.)
394  ENDIF   !(iflag_pbl>1 .AND. iflag_wake>=1 .AND. iflag_pbl_split >=1 )
395
396!==================================
397!  thermiques, poches, convection
398!==================================
399
400! Emanuel
401  found=phyetat0_get(klev,sig1,"sig1","sig1",0.)
402  found=phyetat0_get(klev,w01,"w01","w01",0.)
403
404! Wake
405  found=phyetat0_get(klev,wake_deltat,"WAKE_DELTAT","Delta T wake/env",0.)
406  found=phyetat0_get(klev,wake_deltaq,"WAKE_DELTAQ","Delta hum. wake/env",0.)
407  found=phyetat0_get(1,wake_s,"WAKE_S","Wake frac. area",0.)
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
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.)
417
418! Thermiques
419  found=phyetat0_get(1,zmax0,"ZMAX0","ZMAX0",40.)
420  found=phyetat0_get(1,f0,"F0","F0",1.e-5)
421  found=phyetat0_get(klev+1,fm_therm,"FM_THERM","Thermals mass flux",0.)
422  found=phyetat0_get(klev,entr_therm,"ENTR_THERM","Thermals Entrain.",0.)
423  found=phyetat0_get(klev,detr_therm,"DETR_THERM","Thermals Detrain.",0.)
424
425! ALE/ALP
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.)
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.)
431
432!===========================================
433  ! Read and send field trs to traclmdz
434!===========================================
435
436  IF (type_trac == 'lmdz') THEN
437     DO it=1, nbtr                                                                 
438!!        iiq=niadv(it+2)                                                           ! jyg
439        iiq=niadv(it+nqo)                                                           ! jyg
440        found=phyetat0_get(1,trs(:,it),"trs_"//tname(iiq), &
441              "Surf trac"//tname(iiq),0.)
442     ENDDO
443     CALL traclmdz_from_restart(trs)
444  ENDIF
445
446!--OB now this is for co2i
447  IF (type_trac == 'co2i') THEN
448     IF (carbon_cycle_cpl) THEN
449        ALLOCATE(co2_send(klon), stat=ierr)
450        IF (ierr /= 0) CALL abort_physic('phyetat0', 'pb allocation co2_send', 1)
451        found=phyetat0_get(1,co2_send,"co2_send","co2 send",co2_ppm)
452     ENDIF
453  ENDIF
454
455!===========================================
456!  ondes de gravite / relief
457!===========================================
458
459!  ondes de gravite non orographiques
460  IF (ok_gwd_rando) found = &
461       phyetat0_get(klev,du_gwd_rando,"du_gwd_rando","du_gwd_rando",0.)
462  IF (.NOT. ok_hines .AND. ok_gwd_rando) found &
463       = phyetat0_get(klev,du_gwd_front,"du_gwd_front","du_gwd_front",0.)
464
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
476!===========================================
477! Initialize ocean
478!===========================================
479
480  IF ( type_ocean == 'slab' ) THEN
481      CALL ocean_slab_init(phys_tstep, pctsrf)
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.)
486        ENDIF
487      ELSE
488          DO i=1,nslay
489            WRITE(str2,'(i2.2)') i
490            found=phyetat0_get(1,tslab(:,i),"tslab"//str2,"tslab",0.) 
491          ENDDO
492      ENDIF
493      IF (.NOT. found) THEN
494          PRINT*, "phyetat0: Le champ <tslab> est absent"
495          PRINT*, "Initialisation a tsol_oce"
496          DO i=1,nslay
497              tslab(:,i)=MAX(ftsol(:,is_oce),271.35)
498          ENDDO
499      ENDIF
500
501      ! Sea ice variables
502      IF (version_ocean == 'sicINT') THEN
503          found=phyetat0_get(1,tice,"slab_tice","slab_tice",0.)
504          IF (.NOT. found) THEN
505              PRINT*, "phyetat0: Le champ <tice> est absent"
506              PRINT*, "Initialisation a tsol_sic"
507                  tice(:)=ftsol(:,is_sic)
508          ENDIF
509          found=phyetat0_get(1,seaice,"seaice","seaice",0.)
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.
516              ENDWHERE
517          ENDIF
518      ENDIF !sea ice INT
519  ENDIF ! Slab       
520
521  ! on ferme le fichier
522  CALL close_startphy
523
524  ! Initialize module pbl_surface_mod
525
526  CALL pbl_surface_init(fder, snow, qsurf, tsoil)
527
528  ! Initialize module ocean_cpl_mod for the case of coupled ocean
529  IF ( type_ocean == 'couple' ) THEN
530     CALL ocean_cpl_init(phys_tstep, longitude_deg, latitude_deg)
531  ENDIF
532
533!  CALL init_iophy_new(latitude_deg, longitude_deg)
534
535  ! Initilialize module fonte_neige_mod     
536  CALL fonte_neige_init(run_off_lic_0)
537
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.)
542
543     where (pctsrf(:, is_oce) == 0.)
544        ds_ns = nf90_fill_real
545        dt_ns = nf90_fill_real
546     end where
547  end if
548
549END SUBROUTINE phyetat0
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
562USE print_control_mod, ONLY: lunout
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
601USE print_control_mod, ONLY: lunout
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"
619        call abort_physic("phyetat0", "", 1)
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.