source: trunk/LMDZ.MARS/libf/phymars/phyredem.F90 @ 3509

Last change on this file since 3509 was 3509, checked in by jbclement, 13 days ago

Dynamic + Mars PCM:
Addition of the description for the 'controle' array in the "start.nc" and "startfi.nc" files. It is given by the variable 'controle_descriptor' whose the element 'controle_descriptor(i)' explains 'controle(i)'.
JBC

File size: 15.4 KB
Line 
1module phyredem
2
3implicit none
4
5!=======================================================================
6contains
7!=======================================================================
8
9subroutine physdem0(filename,lonfi,latfi,nsoil,ngrid,nlay,nq, &
10                         phystep,day_ini,time,airefi, &
11                         alb,ith,def_slope, &
12                         subslope_dist)
13
14! create physics restart file and write time-independent variables
15  use comsoil_h, only: inertiedat, volcapa, mlayer
16  use geometry_mod, only: cell_area
17  use surfdat_h, only: zmea, zstd, zsig, zgam, zthe, &
18                       z0_default, albedice, emisice, emissiv, &
19                       iceradius, dtemisice, phisfi, z0,   &
20                       hmons,summit,base,watercaptag
21  use dimradmars_mod, only: tauvis
22  use iostart, only : open_restartphy, close_restartphy, &
23                      put_var, put_field, length, ldscrpt, ndscrpt
24  use mod_grid_phy_lmdz, only : klon_glo
25  use planete_h, only: aphelie, emin_turb, lmixmin, obliquit, &
26                       peri_day, periheli, year_day
27  use comcstfi_h, only: g, mugaz, omeg, rad, rcp
28  use time_phylmdz_mod, only: daysec
29  use comslope_mod, ONLY: nslope
30  USE paleoclimate_mod, ONLY: paleoclimate, h2o_ice_depth, lag_co2_ice, d_coef
31
32  implicit none
33 
34  character(len=*), intent(in) :: filename
35  real,intent(in) :: lonfi(ngrid)
36  real,intent(in) :: latfi(ngrid)
37  integer,intent(in) :: nsoil
38  integer,intent(in) :: ngrid
39  integer,intent(in) :: nlay
40  integer,intent(in) :: nq
41  real,intent(in) :: phystep
42  real,intent(in) :: day_ini
43  real,intent(in) :: time
44  real,intent(in) :: airefi(ngrid)
45  real,intent(in) :: alb(ngrid)
46  real,intent(in) :: ith(ngrid,nsoil) ! thermal inertia for present day
47  real, intent(in) :: def_slope(nslope+1) !boundaries for bining of the slopes
48  real, intent(in) :: subslope_dist(ngrid,nslope) !undermesh statistics
49
50  character(ndscrpt), dimension(ldscrpt), parameter :: dscrpt_tab_cntrl = (/ &
51      "(1)  Total number of nodes on physics grid        ", &
52      "(2)  Number of atmospheric layers                 ", &
53      "(3)  Initial day                                  ", &
54      "(4)  Initial time of day                          ", &
55      "(5)  Radius of Mars (m) ~3397200                  ", &
56      "(6)  Rotation rate (rad.s-1)                      ", &
57      "(7)  Gravity (m.s-2) ~3.72                        ", &
58      "(8)  Molar mass of the atmosphere (g.mol-1) ~43.49", &
59      "(9)  = r/Cp ~0.256793 (=kappa in the dynamics)    ", &
60      "(10) Length of a sol (s) ~88775                   ", &
61      "(11) Time step in the physics                     ", &
62      "(12) -                                            ", &
63      "(13) -                                            ", &
64      "(14) Length of year (sols) ~668.6                 ", &
65      "(15) Min. Sun-Mars distance (Mkm) ~206.66         ", &
66      "(16) Max. Sun-Mars distance (Mkm) ~249.22         ", &
67      "(17) Date of perihelion (sols since N. spring)    ", &
68      "(18) Obliquity of the planet (deg) ~23.98         ", &
69      "(19) Default surface roughness (m) ~0.01          ", &
70      "(20) Mixing length ~100                           ", &
71      "(21) Minimal energy ~1.e-8                        ", &
72      "(22) Albedo of northern cap ~0.5                  ", &
73      "(23) Albedo of southern cap ~0.5                  ", &
74      "(24) Emissivity of northern cap ~0.95             ", &
75      "(25) Emissivity of southern cap ~0.95             ", &
76      "(26) Emissivity of martian soil ~.95              ", &
77      "(27) Mean visible optical depth                   ", &
78      "(28) -                                            ", &
79      "(29) -                                            ", &
80      "(30) -                                            ", &
81      "(31) Mean scat radius of CO2 snow (north)         ", &
82      "(32) Mean scat radius of CO2 snow (south)         ", &
83      "(33) Time scale for snow metamorphism (north)     ", &
84      "(34) Time scale for snow metamorphism (south)     ", &
85      "(35) Soil volumetric heat capacity                "/)
86  real, dimension(length) :: tab_cntrl ! nb "length=100" defined in iostart module
87  integer                 :: ig
88  real, dimension(ngrid)  :: watercaptag_tmp
89
90  ! Create physics start file
91  call open_restartphy(filename)
92
93  ! Build tab_cntrl(:) array
94  tab_cntrl(:)=0.0
95  !cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
96  ! Fill control array tab_cntrl(:) with parameters for this run
97  !cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
98  ! Informations on the physics grid
99  tab_cntrl(1) = float(klon_glo)     ! Total number of nodes on physics grid
100  tab_cntrl(2) = float(nlay)         ! Number of atmospheric layers
101  tab_cntrl(3) = day_ini + int(time) ! Initial day
102  tab_cntrl(4) = time -int(time)     ! Initial time of day
103
104  ! Informations about Mars, used by dynamics and physics
105  tab_cntrl(5) = rad      ! Radius of Mars (m) ~3397200
106  tab_cntrl(6) = omeg     ! Rotation rate (rad.s-1)
107  tab_cntrl(7) = g        ! Gravity (m.s-2) ~3.72
108  tab_cntrl(8) = mugaz    ! Molar mass of the atmosphere (g.mol-1) ~43.49
109  tab_cntrl(9) = rcp      ! = r/cp  ~0.256793 (=kappa dans dynamique)
110  tab_cntrl(10) = daysec  ! Length of a sol (s) ~88775
111  tab_cntrl(11) = phystep ! Time step in the physics
112
113  ! Informations about Mars, only for physics
114  tab_cntrl(14) = year_day  ! Length of year (sols) ~668.6
115  tab_cntrl(15) = periheli  ! Min. Sun-Mars distance (Mkm) ~206.66
116  tab_cntrl(16) = aphelie   ! Max. SUn-Mars distance (Mkm) ~249.22
117  tab_cntrl(17) = peri_day  ! Date of perihelion (sols since N. spring)
118  tab_cntrl(18) = obliquit  ! Obliquity of the planet (deg) ~23.98
119
120  ! Boundary layer and turbulence
121  tab_cntrl(19) = z0_default ! Default surface roughness (m) ~0.01
122  tab_cntrl(20) = lmixmin    ! Mixing length ~100
123  tab_cntrl(21) = emin_turb  ! Minimal energy ~1.e-8
124
125  ! Optical properties of polar caps and ground emissivity
126  tab_cntrl(22) = albedice(1)  ! Albedo of northern cap ~0.5
127  tab_cntrl(23) = albedice(2)  ! Albedo of southern cap ~0.5
128  tab_cntrl(24) = emisice(1)   ! Emissivity of northern cap ~0.95
129  tab_cntrl(25) = emisice(2)   ! Emissivity of southern cap ~0.95
130  tab_cntrl(26) = emissiv      ! Emissivity of martian soil ~.95
131  tab_cntrl(31) = iceradius(1) ! Mean scat radius of CO2 snow (north)
132  tab_cntrl(32) = iceradius(2) ! Mean scat radius of CO2 snow (south)
133  tab_cntrl(33) = dtemisice(1) ! Time scale for snow metamorphism (north)
134  tab_cntrl(34) = dtemisice(2) ! Time scale for snow metamorphism (south)
135
136  ! dust aerosol properties
137  tab_cntrl(27) = tauvis ! Mean visible optical depth
138
139  ! Soil properties:
140  tab_cntrl(35) = volcapa ! Soil volumetric heat capacity
141
142  ! Write the controle array
143  call put_var("controle","Control parameters",tab_cntrl)
144
145  ! Write the controle array descriptor
146  call put_var("controle_descriptor","Description of control parameters",dscrpt_tab_cntrl)
147 
148  ! Write the mid-layer depths
149  call put_var("soildepth","Soil mid-layer depth",mlayer)
150 
151  ! Write longitudes
152  call put_field("longitude","Longitudes of physics grid",lonfi)
153 
154  ! Write latitudes
155  call put_field("latitude","Latitudes of physics grid",latfi)
156 
157  ! Write mesh areas
158  call put_field("area","Mesh area",cell_area)
159 
160  ! Write surface geopotential
161  call put_field("phisfi","Geopotential at the surface",phisfi)
162 
163  ! Write surface albedo
164  call put_field("albedodat","Albedo of bare ground",alb)
165 
166  ! Subgrid topogaphy variables
167  call put_field("ZMEA","Relief: mean relief",zmea)
168  call put_field("ZSTD","Relief: standard deviation",zstd)
169  call put_field("ZSIG","Relief: sigma parameter",zsig)
170  call put_field("ZGAM","Relief: gamma parameter",zgam)
171  call put_field("ZTHE","Relief: theta parameter",zthe)
172  call put_field("hmons","Relief: hmons parameter (summit - base)",hmons)
173  call put_field("summit","Relief: altitude from the aeroid of the summit of the highest subgrid topography",summit)
174  call put_field("base","Relief: altitude from the aeroid of the base of the highest subgrid topography",base)
175     
176  ! Soil thermal inertia
177  call put_field("inertiedat","Soil thermal inertia - present day",ith)
178 
179  ! Surface roughness length
180  call put_field("z0","Surface roughness length",z0)
181
182  ! Water cap
183  watercaptag_tmp = 0
184  where (watercaptag) watercaptag_tmp = 1
185  call put_field("watercaptag","Infinite water reservoir",watercaptag_tmp)
186
187  ! Sub grid scale slope parametrization
188  call put_var("def_slope","slope criterium stages",def_slope)
189  call put_field("subslope_dist","under mesh slope distribution",subslope_dist)
190
191  ! Paleoclimate outputs
192  if (paleoclimate) then
193    call put_field("h2o_ice_depth","Depth to the fisrt H2O ice",h2o_ice_depth)
194    call put_field("lag_co2_ice","Depth of the CO2 lags",lag_co2_ice)
195    call put_field("d_coef","Diffusion coefficent",d_coef)
196  endif
197
198  ! Close file
199  call close_restartphy
200 
201end subroutine physdem0
202
203!=======================================================================
204
205subroutine physdem1(filename,nsoil,ngrid,nlay,nq,nqsoil, &
206                    phystep,time,tsurf,tsoil,inertiesoil, &
207                    albedo,emis,q2,qsurf,qsoil,&
208                    tauscaling,totcloudfrac,wstar, &
209                    watercap,perennial_co2ice)
210  ! write time-dependent variable to restart file
211  use iostart, only : open_restartphy, close_restartphy, &
212                      put_var, put_field
213  use tracer_mod, only: noms ! tracer names
214  use nonoro_gwd_ran_mod, only: du_nonoro_gwd, dv_nonoro_gwd
215  use compute_dtau_mod, only: dtau
216  use dust_rad_adjust_mod, only: dust_rad_adjust_prev,dust_rad_adjust_next
217  use dust_param_mod, only: dustscaling_mode
218  use comsoil_h,only: flux_geo,adsorption_soil,igcm_h2o_vap_soil, &
219                      igcm_h2o_ice_soil,igcm_h2o_vap_ads
220  use comslope_mod, only: nslope
221  use paleoclimate_mod, only: paleoclimate
222  implicit none
223 
224  include "callkeys.h"
225 
226  character(len=*),intent(in) :: filename
227  integer,intent(in) :: nsoil
228  integer,intent(in) :: ngrid
229  integer,intent(in) :: nlay
230  integer,intent(in) :: nq
231  integer,intent(in) :: nqsoil
232  real,intent(in) :: phystep
233  real,intent(in) :: time
234  real,intent(in) :: tsurf(ngrid,nslope)
235  real,intent(in) :: tsoil(ngrid,nsoil,nslope)
236  real,intent(in) :: inertiesoil(ngrid,nsoil,nslope)
237  real,intent(in) :: albedo(ngrid,2,nslope)
238  real,intent(in) :: emis(ngrid,nslope)
239  real,intent(in) :: q2(ngrid,nlay+1)
240  real,intent(in) :: qsurf(ngrid,nq,nslope)
241  real,intent(in) :: qsoil(ngrid,nsoil,nqsoil,nslope)
242  real,intent(in) :: tauscaling(ngrid)
243  real,intent(in) :: totcloudfrac(ngrid)
244  real,intent(in) :: wstar(ngrid)
245  real,intent(in) :: watercap(ngrid,nslope)
246  real,intent(in) :: perennial_co2ice(ngrid,nslope)
247
248  integer :: iq
249  character(len=30) :: txt ! to store some text
250  ! indexes of water vapour & water ice tracers (if any):
251  integer :: i_h2o_vap=0
252  integer :: i_h2o_ice=0
253  integer :: i_hdo_vap=0
254  integer :: i_hdo_ice=0
255
256 
257  ! Open file
258  call open_restartphy(filename)
259 
260  ! First variable to write must be "Time", in order to correctly
261  ! set time counter in file
262  call put_var("Time","Temps de simulation",time)
263
264  ! Water ice layer
265  call put_field("watercap","H2O ice cover",watercap,time)
266 
267  ! Perennial CO2 ice layer
268  if (paleoclimate) call put_field("perennial_co2ice","CO2 ice cover",perennial_co2ice,time)
269
270  ! Surface temperature
271  call put_field("tsurf","Surface temperature",tsurf,time)
272
273  ! Soil temperature
274  call put_field("inertiesoil","Soil thermal inertia",inertiesoil,time)
275 
276  ! Soil temperature
277  call put_field("tsoil","Soil temperature",tsoil,time)
278 
279  ! Albedo of the surface
280  call put_field("albedo","Surface albedo",albedo(:,1,:),time)
281 
282  ! Emissivity of the surface
283  call put_field("emis","Surface emissivity",emis,time)
284 
285  ! Planetary Boundary Layer
286  call put_field("q2","pbl wind variance",q2,time)
287
288  ! Sub-grid cloud fraction
289  call put_field("totcloudfrac","Total cloud fraction",totcloudfrac,time)
290 
291  ! Dust conversion factor
292  ! Only to be read by newstart to convert to actual dust values
293  ! Or by any user who wants to reconstruct dust, opacity from the start files.
294  call put_field("tauscaling","dust conversion factor",tauscaling,time)
295
296  ! Radiative scaling coefficients
297  if (dustscaling_mode==2) then
298    call put_field("dust_rad_adjust_prev", &
299                   "radiative dust adjustement factor prev. sol", &
300                   dust_rad_adjust_prev,time)
301    call put_field("dust_rad_adjust_next", &
302                   "radiative dust adjustement factor next sol", &
303                   dust_rad_adjust_next,time)
304  endif
305 
306  if (dustinjection.gt.0) then
307    call put_field("dtau","dust opacity difference between GCM and scenario",&
308                   dtau,time)
309  endif
310
311  if (calltherm) then
312    call put_field("wstar","Max vertical velocity in thermals",wstar,time)
313  endif
314
315  ! Tracers on the surface
316  ! preliminary stuff: look for water vapour & water ice tracers (if any)
317  do iq=1,nq
318    if (noms(iq).eq."h2o_vap") then
319      i_h2o_vap=iq
320    endif
321    if (noms(iq).eq."h2o_ice") then
322      i_h2o_ice=iq
323    endif
324
325  ! look for HDO vapour & HDO ice tracers (if any)
326    if (noms(iq).eq."hdo_vap") then
327      i_hdo_vap=iq
328    endif
329    if (noms(iq).eq."hdo_ice") then
330      i_hdo_ice=iq
331    endif
332  enddo
333
334 
335  if (nq.gt.0) then
336    do iq=1,nq
337      txt=noms(iq)
338      ! Exception: there is no water vapour surface tracer
339      if (txt.eq."h2o_vap") then
340        write(*,*)"physdem1: skipping water vapour tracer"
341        if (i_h2o_ice.eq.i_h2o_vap) then
342          ! then there is no "water ice" tracer; but still
343          ! there is some water ice on the surface
344          write(*,*)"          writing water ice instead"
345          txt="h2o_ice"
346        else
347          ! there is a "water ice" tracer which has been / will be
348          ! delt with in due time
349          cycle
350        endif ! of if (igcm_h2o_ice.eq.igcm_h2o_vap)
351      endif ! of if (txt.eq."h2o_vap")
352
353      if (txt.eq."hdo_vap") then
354        write(*,*)"physdem1: skipping HDO vapour tracer"
355        if (i_hdo_ice.eq.i_hdo_vap) then
356          ! then there is no "water ice" tracer; but still
357          ! there is some water ice on the surface
358          write(*,*)"          writing HDO ice instead"
359          txt="hdo_ice"
360        else
361          ! there is a "water ice" tracer which has been / will be
362          ! delt with in due time
363          cycle
364        endif ! of if (igcm_hdo_ice.eq.igcm_hdo_vap)
365      endif ! of if (txt.eq."hdo_vap")
366
367      ! co2_ice has been added to qsurf(:,igcm_co2) in co2condens4micro
368      if (txt.eq."co2_ice") then
369        write(*,*)"physdem1: skipping co2_ice tracer"
370        cycle
371      end if     
372
373      call put_field(trim(txt),"tracer on surface",qsurf(:,iq,:),time)
374    enddo
375  endif
376
377  ! Non-orographic gavity waves
378  if (calllott_nonoro) then
379     call put_field("du_nonoro_gwd","Zonal wind tendency due to GW",du_nonoro_gwd,time)
380     call put_field("dv_nonoro_gwd","Meridional wind tendency due to GW",dv_nonoro_gwd,time)
381  endif
382
383  ! Geothermal Flux
384     call put_field('flux_geo','Geothermal flux',flux_geo,time)
385
386  ! Adsorption
387  if (adsorption_soil) then
388     call put_field("h2o_vap_soil","subsurface water vapour", &
389                        qsoil(:,:,igcm_h2o_vap_soil,:), time)
390     call put_field("h2o_ice_soil","subsurface water ice", &
391                        qsoil(:,:,igcm_h2o_ice_soil,:), time)
392     call put_field("h2o_vap_ads", "adsorbed water", &
393                        qsoil(:,:,igcm_h2o_vap_ads,:), time)
394  endif
395
396  ! Close file
397  call close_restartphy
398 
399end subroutine physdem1
400
401end module phyredem
Note: See TracBrowser for help on using the repository browser.