source: trunk/LMDZ.GENERIC/libf/phystd/phyetat0_mod.F90 @ 1704

Last change on this file since 1704 was 1669, checked in by emillour, 8 years ago

Generic GCM:
Added possibility to run without a startfi.nc file (mainly usefull for
tests with coupling with dynamico dynamical core):

  • added flag "startphy_file" flag (.false. if doing an "academic" start on the physics side).
  • turned phyetat0.F90 into module phyetat0_mod.F90
  • turned tabfi.F into module tabfi_mod.F90 and added handling of startphy_file==.false. case
  • extra initializations in physiq_mod for startphy_file==.false. case.

EM
File size: 11.1 KB
Line 
1module phyetat0_mod
2
3implicit none
4
5contains
6
7subroutine phyetat0 (startphy_file, &
8                     ngrid,nlayer,fichnom,tab0,Lmodif,nsoil,nq, &
9                     day_ini,time,tsurf,tsoil, &
10                     emis,q2,qsurf,cloudfrac,totcloudfrac,hice, &
11                     rnat,pctsrf_sic,tslab,tsea_ice,sea_ice)
12
13  use tabfi_mod, only: tabfi
14  USE tracer_h, ONLY: noms
15  USE surfdat_h, only: phisfi, albedodat, zmea, zstd, zsig, zgam, zthe
16  use iostart, only: nid_start, open_startphy, close_startphy, &
17                     get_field, get_var, inquire_field, &
18                     inquire_dimension, inquire_dimension_length
19  use slab_ice_h, only: noceanmx
20  use callkeys_mod, only: CLFvarying
21
22  implicit none
23
24!======================================================================
25!  Arguments:
26!  ---------
27!  inputs:
28  logical,intent(in) :: startphy_file ! .true. if reading start file
29  integer,intent(in) :: ngrid
30  integer,intent(in) :: nlayer
31  character*(*),intent(in) :: fichnom ! "startfi.nc" file
32  integer,intent(in) :: tab0
33  integer,intent(in) :: Lmodif
34  integer,intent(in) :: nsoil ! # of soil layers
35  integer,intent(in) :: nq
36  integer,intent(out) :: day_ini
37  real,intent(out) :: time
38
39!  outputs:
40  real,intent(out) :: tsurf(ngrid) ! surface temperature
41  real,intent(out) :: tsoil(ngrid,nsoil) ! soil temperature
42  real,intent(out) :: emis(ngrid) ! surface emissivity
43  real,intent(out) :: q2(ngrid,nlayer+1) !
44  real,intent(out) :: qsurf(ngrid,nq) ! tracers on surface
45! real co2ice(ngrid) ! co2 ice cover
46  real,intent(out) :: cloudfrac(ngrid,nlayer)
47  real,intent(out) :: hice(ngrid), totcloudfrac(ngrid)
48  real,intent(out) :: pctsrf_sic(ngrid),tslab(ngrid,noceanmx) 
49  real,intent(out) :: tsea_ice(ngrid),sea_ice(ngrid)
50  real,intent(out) :: rnat(ngrid)
51
52!======================================================================
53!  Local variables:
54
55!      INTEGER radpas
56!      REAL co2_ppm
57!      REAL solaire
58
59      real xmin,xmax ! to display min and max of a field
60!
61      INTEGER ig,iq,lmax
62      INTEGER nid, nvarid
63      INTEGER ierr, i, nsrf
64!      integer isoil
65!      INTEGER length
66!      PARAMETER (length=100)
67      CHARACTER*7 str7
68      CHARACTER*2 str2
69      CHARACTER*1 yes
70!
71      REAL p_rad,p_omeg,p_g,p_cpp,p_mugaz,p_daysec
72      INTEGER nqold
73
74! flag which identifies if 'startfi.nc' file is using old names (qsurf01,...)
75!      logical :: oldtracernames=.false.
76      integer :: count
77      character(len=30) :: txt ! to store some text
78     
79      INTEGER :: indextime=1 ! index of selected time, default value=1
80      logical :: found
81     
82      character(len=8) :: modname="phyetat0"
83
84!
85! ALLOCATE ARRAYS IN surfdat_h
86!
87IF (.not. ALLOCATED(albedodat)) ALLOCATE(albedodat(ngrid))
88IF (.not. ALLOCATED(phisfi)) ALLOCATE(phisfi(ngrid))
89IF (.not. ALLOCATED(zmea)) ALLOCATE(zmea(ngrid))
90IF (.not. ALLOCATED(zstd)) ALLOCATE(zstd(ngrid))
91IF (.not. ALLOCATED(zsig)) ALLOCATE(zsig(ngrid))
92IF (.not. ALLOCATED(zgam)) ALLOCATE(zgam(ngrid))
93IF (.not. ALLOCATED(zthe)) ALLOCATE(zthe(ngrid))
94
95if (startphy_file) then
96  ! open physics initial state file:
97  call open_startphy(fichnom)
98
99  ! possibility to modify tab_cntrl in tabfi
100  write(*,*)
101  write(*,*) 'TABFI in phyeta0: Lmodif=',Lmodif," tab0=",tab0
102  call tabfi (ngrid,nid_start,Lmodif,tab0,day_ini,lmax,p_rad, &
103                   p_omeg,p_g,p_cpp,p_mugaz,p_daysec,time)
104
105else ! "academic" initialization of planetary parameters via tabfi
106  call tabfi (ngrid,0,0,0,day_ini,lmax,p_rad, &
107                   p_omeg,p_g,p_cpp,p_mugaz,p_daysec,time)
108endif ! of if (startphy_file)
109
110if (startphy_file) then
111  ! Load surface geopotential:
112  call get_field("phisfi",phisfi,found)
113  if (.not.found) then
114    call abort_physic(modname,"Failed loading <phisfi>",1)
115  endif
116else
117  phisfi(:)=0
118endif ! of if (startphy_file)
119write(*,*) "phyetat0: surface geopotential <phisfi> range:", &
120               minval(phisfi), maxval(phisfi)
121
122if (startphy_file) then
123  ! Load bare ground albedo:
124  call get_field("albedodat",albedodat,found)
125  if (.not.found) then
126    call abort_physic(modname,"Failed loading <albedodat>",1)
127  endif
128else
129  albedodat(:)=0.5 ! would be better to read value from def file...
130endif ! of if (startphy_file)
131write(*,*) "phyetat0: Bare ground albedo <albedodat> range:", &
132             minval(albedodat), maxval(albedodat)
133
134! ZMEA
135if (startphy_file) then
136  call get_field("ZMEA",zmea,found)
137  if (.not.found) then
138    call abort_physic(modname,"Failed loading <ZMEA>",1)
139  endif
140else
141  zmea(:)=0
142endif ! of if (startphy_file)
143write(*,*) "phyetat0: <ZMEA> range:", &
144             minval(zmea), maxval(zmea)
145
146! ZSTD
147if (startphy_file) then
148  call get_field("ZSTD",zstd,found)
149  if (.not.found) then
150    call abort_physic(modname,"Failed loading <ZSTD>",1)
151  endif
152else
153  zstd(:)=0
154endif ! of if (startphy_file)
155write(*,*) "phyetat0: <ZSTD> range:", &
156             minval(zstd), maxval(zstd)
157
158! ZSIG
159if (startphy_file) then
160  call get_field("ZSIG",zsig,found)
161  if (.not.found) then
162    call abort_physic(modname,"Failed loading <ZSIG>",1)
163  endif
164else
165  zsig(:)=0
166endif ! of if (startphy_file)
167write(*,*) "phyetat0: <ZSIG> range:", &
168             minval(zsig), maxval(zsig)
169
170! ZGAM
171if (startphy_file) then
172  call get_field("ZGAM",zgam,found)
173  if (.not.found) then
174    call abort_physic(modname,"Failed loading <ZGAM>",1)
175  endif
176else
177  zgam(:)=0
178endif ! of if (startphy_file)
179write(*,*) "phyetat0: <ZGAM> range:", &
180             minval(zgam), maxval(zgam)
181
182! ZTHE
183if (startphy_file) then
184  call get_field("ZTHE",zthe,found)
185  if (.not.found) then
186    call abort_physic(modname,"Failed loading <ZTHE>",1)
187  endif
188else
189  zthe(:)=0
190endif ! of if (startphy_file)
191write(*,*) "phyetat0: <ZTHE> range:", &
192             minval(zthe), maxval(zthe)
193
194! Surface temperature :
195if (startphy_file) then
196  call get_field("tsurf",tsurf,found,indextime)
197  if (.not.found) then
198    call abort_physic(modname,"Failed loading <tsurf>",1)
199  endif
200else
201  tsurf(:)=0 ! will be updated afterwards in physiq !
202endif ! of if (startphy_file)
203write(*,*) "phyetat0: Surface temperature <tsurf> range:", &
204             minval(tsurf), maxval(tsurf)
205
206! Surface emissivity
207if (startphy_file) then
208  call get_field("emis",emis,found,indextime)
209  if (.not.found) then
210    call abort_physic(modname,"Failed loading <emis>",1)
211  endif
212else
213  emis(:)=1 ! would be better to read value from def file...
214endif ! of if (startphy_file)
215write(*,*) "phyetat0: Surface emissivity <emis> range:", &
216             minval(emis), maxval(emis)
217
218! Cloud fraction (added by BC 2010)
219if (CLFvarying) then
220  if (startphy_file) then
221    call get_field("cloudfrac",cloudfrac,found,indextime)
222    if (.not.found) then
223      call abort_physic(modname,"Failed loading <cloudfrac>",1)
224    endif
225  else
226    cloudfrac(:,:)=0
227  endif ! of if (startphy_file)
228  write(*,*) "phyetat0: Cloud fraction <cloudfrac> range:", &
229               minval(cloudfrac), maxval(cloudfrac)
230else
231  cloudfrac(:,:)=0.0
232endif ! of if (CLFvarying)
233
234! Total cloud fraction (added by BC 2010)
235if (CLFvarying) then
236  if (startphy_file) then
237    call get_field("totcloudfrac",totcloudfrac,found,indextime)
238    if (.not.found) then
239      call abort_physic(modname,"Failed loading <totcloudfrac>",1)
240    endif
241  else
242    totcloudfrac(:)=0
243  endif ! of if (startphy_file)
244  write(*,*) "phyetat0: Total cloud fraction <totcloudfrac> range:", &
245             minval(totcloudfrac), maxval(totcloudfrac)
246else
247  totcloudfrac(:)=0.0
248endif ! of if (CLFvarying)
249
250! Height of oceanic ice (added by BC 2010)
251if (startphy_file) then
252  call get_field("hice",hice,found,indextime)
253  if (.not.found) then
254    write(*,*) "phyetat0: Failed loading <hice>"
255    !  call abort
256    hice(:)=0
257  endif
258else
259  hice(:)=0
260endif ! of if (startphy_file)
261write(*,*) "phyetat0: Height of oceanic ice <hice> range:", &
262             minval(hice), maxval(hice)
263
264! SLAB OCEAN (added by BC 2014)
265if (startphy_file) then
266  ! nature of the surface
267  call get_field("rnat",rnat,found,indextime)
268  if (.not.found) then
269    write(*,*) "phyetat0: Failed loading <rnat>"
270        rnat(1:ngrid)=1.
271  else
272      do ig=1,ngrid
273        if((nint(rnat(ig)).eq.2).or.(nint(rnat(ig)).eq.0))then
274          rnat(ig)=0.
275        else
276          rnat(ig)=1.
277        endif     
278      enddo
279  endif ! of if (.not.found)
280else
281  rnat(:)=1
282endif ! of if (startphy_file)
283write(*,*) "phyetat0: Nature of surface <rnat> range:", &
284             minval(rnat), maxval(rnat)
285
286if (startphy_file) then
287  ! Pourcentage of sea ice cover
288  call get_field("pctsrf_sic",pctsrf_sic,found,indextime)
289  if (.not.found) then
290    write(*,*) "phyetat0: Failed loading <pctsrf_sic>"
291    pctsrf_sic(1:ngrid)=0.
292  endif
293else
294  pctsrf_sic(:)=0
295endif ! of if (startphy_file)
296write(*,*) "phyetat0: Pourcentage of sea ice cover <pctsrf_sic> range:", &
297             minval(pctsrf_sic), maxval(pctsrf_sic)
298
299if (startphy_file) then
300  ! Slab ocean temperature (2 layers)
301  call get_field("tslab",tslab,found,indextime)
302  if (.not.found) then
303    write(*,*) "phyetat0: Failed loading <tslab>"
304    do iq=1,noceanmx
305      tslab(1:ngrid,iq)=tsurf(1:ngrid)
306    enddo
307  endif
308else
309  do iq=1,noceanmx
310    tslab(1:ngrid,iq)=tsurf(1:ngrid)
311  enddo
312endif ! of if (startphy_file)
313write(*,*) "phyetat0: Slab ocean temperature <tslab> range:", &
314             minval(tslab), maxval(tslab)
315
316if (startphy_file) then
317  ! Oceanic ice temperature
318  call get_field("tsea_ice",tsea_ice,found,indextime)
319  if (.not.found) then
320    write(*,*) "phyetat0: Failed loading <tsea_ice>"
321    tsea_ice(1:ngrid)=273.15-1.8
322  endif
323else
324  tsea_ice(1:ngrid)=273.15-1.8
325endif ! of if (startphy_file)
326write(*,*) "phyetat0: Oceanic ice temperature <tsea_ice> range:", &
327             minval(tsea_ice), maxval(tsea_ice)
328
329if (startphy_file) then
330  !  Oceanic ice quantity (kg/m^2)
331  call get_field("sea_ice",sea_ice,found,indextime)
332  if (.not.found) then
333    write(*,*) "phyetat0: Failed loading <sea_ice>"
334    tsea_ice(1:ngrid)=0.
335  endif
336else
337  tsea_ice(1:ngrid)=0
338endif ! of if (startphy_file)
339write(*,*) "phyetat0: Oceanic ice quantity <sea_ice> range:", &
340             minval(sea_ice), maxval(sea_ice)
341
342
343! pbl wind variance
344if (startphy_file) then
345  call get_field("q2",q2,found,indextime)
346  if (.not.found) then
347    call abort_physic(modname,"Failed loading <q2>",1)
348  endif
349else
350  q2(:,:)=0
351endif ! of if (startphy_file)
352write(*,*) "phyetat0: PBL wind variance <q2> range:", &
353             minval(q2), maxval(q2)
354
355! tracer on surface
356if (nq.ge.1) then
357  do iq=1,nq
358    txt=noms(iq)
359    if (startphy_file) then
360      call get_field(txt,qsurf(:,iq),found,indextime)
361      if (.not.found) then
362        write(*,*) "phyetat0: Failed loading <",trim(txt),">"
363        write(*,*) "         ",trim(txt)," is set to zero"
364        qsurf(:,iq) = 0.
365      endif
366    else
367      qsurf(:,iq)=0
368    endif ! of if (startphy_file)
369    write(*,*) "phyetat0: Surface tracer <",trim(txt),"> range:", &
370                 minval(qsurf(:,iq)), maxval(qsurf(:,iq))
371  enddo! of do iq=1,nq
372endif ! of if (nq.ge.1)
373
374
375if (startphy_file) then
376  ! Call to soil_settings, in order to read soil temperatures,
377  ! as well as thermal inertia and volumetric heat capacity
378  call soil_settings(nid_start,ngrid,nsoil,tsurf,tsoil,indextime)
379endif ! of if (startphy_file)
380!
381! close file:
382!
383if (startphy_file) call close_startphy
384
385end subroutine phyetat0
386
387end module phyetat0_mod
Note: See TracBrowser for help on using the repository browser.