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