[4603] | 1 | #include "conf_gcm.F90" |
---|
[2018] | 2 | |
---|
[2017] | 3 | ! |
---|
[3540] | 4 | ! $Id: 1DUTILS.h 5084 2024-07-19 16:40:44Z idelkadi $ |
---|
[2017] | 5 | ! |
---|
[2019] | 6 | ! |
---|
| 7 | ! |
---|
[2017] | 8 | SUBROUTINE conf_unicol |
---|
[2019] | 9 | ! |
---|
[2017] | 10 | #ifdef CPP_IOIPSL |
---|
| 11 | use IOIPSL |
---|
| 12 | #else |
---|
| 13 | ! if not using IOIPSL, we still need to use (a local version of) getin |
---|
| 14 | use ioipsl_getincom |
---|
| 15 | #endif |
---|
[2311] | 16 | USE print_control_mod, ONLY: lunout |
---|
[2017] | 17 | IMPLICIT NONE |
---|
[2019] | 18 | !----------------------------------------------------------------------- |
---|
| 19 | ! Auteurs : A. Lahellec . |
---|
| 20 | ! |
---|
| 21 | ! Declarations : |
---|
| 22 | ! -------------- |
---|
[2017] | 23 | |
---|
[4603] | 24 | #include "compar1d.h" |
---|
| 25 | #include "flux_arp.h" |
---|
| 26 | #include "tsoilnudge.h" |
---|
| 27 | #include "fcg_gcssold.h" |
---|
| 28 | #include "fcg_racmo.h" |
---|
[2019] | 29 | ! |
---|
| 30 | ! |
---|
| 31 | ! local: |
---|
| 32 | ! ------ |
---|
[2017] | 33 | |
---|
[2019] | 34 | ! CHARACTER ch1*72,ch2*72,ch3*72,ch4*12 |
---|
[2017] | 35 | |
---|
[2019] | 36 | ! |
---|
| 37 | ! ------------------------------------------------------------------- |
---|
| 38 | ! |
---|
| 39 | ! ......... Initilisation parametres du lmdz1D .......... |
---|
| 40 | ! |
---|
| 41 | !--------------------------------------------------------------------- |
---|
| 42 | ! initialisations: |
---|
| 43 | ! ---------------- |
---|
[2017] | 44 | |
---|
| 45 | !Config Key = lunout |
---|
| 46 | !Config Desc = unite de fichier pour les impressions |
---|
| 47 | !Config Def = 6 |
---|
| 48 | !Config Help = unite de fichier pour les impressions |
---|
| 49 | !Config (defaut sortie standard = 6) |
---|
| 50 | lunout=6 |
---|
| 51 | ! CALL getin('lunout', lunout) |
---|
| 52 | IF (lunout /= 5 .and. lunout /= 6) THEN |
---|
| 53 | OPEN(lunout,FILE='lmdz.out') |
---|
| 54 | ENDIF |
---|
| 55 | |
---|
| 56 | !Config Key = prt_level |
---|
[2716] | 57 | !Config Desc = niveau d'impressions de debogage |
---|
[2017] | 58 | !Config Def = 0 |
---|
[2716] | 59 | !Config Help = Niveau d'impression pour le debogage |
---|
[2017] | 60 | !Config (0 = minimum d'impression) |
---|
[2019] | 61 | ! prt_level = 0 |
---|
| 62 | ! CALL getin('prt_level',prt_level) |
---|
[2017] | 63 | |
---|
[2019] | 64 | !----------------------------------------------------------------------- |
---|
| 65 | ! Parametres de controle du run: |
---|
| 66 | !----------------------------------------------------------------------- |
---|
[2017] | 67 | |
---|
| 68 | !Config Key = restart |
---|
| 69 | !Config Desc = on repart des startphy et start1dyn |
---|
| 70 | !Config Def = false |
---|
| 71 | !Config Help = les fichiers restart doivent etre renomme en start |
---|
[2019] | 72 | restart =.false. |
---|
[2017] | 73 | CALL getin('restart',restart) |
---|
| 74 | |
---|
| 75 | !Config Key = forcing_type |
---|
| 76 | !Config Desc = defines the way the SCM is forced: |
---|
| 77 | !Config Def = 0 |
---|
| 78 | !!Config Help = 0 ==> forcing_les = .true. |
---|
| 79 | ! initial profiles from file prof.inp.001 |
---|
| 80 | ! no forcing by LS convergence ; |
---|
| 81 | ! surface temperature imposed ; |
---|
| 82 | ! radiative cooling may be imposed (iflag_radia=0 in physiq.def) |
---|
| 83 | ! = 1 ==> forcing_radconv = .true. |
---|
| 84 | ! idem forcing_type = 0, but the imposed radiative cooling |
---|
| 85 | ! is set to 0 (hence, if iflag_radia=0 in physiq.def, |
---|
| 86 | ! then there is no radiative cooling at all) |
---|
| 87 | ! = 2 ==> forcing_toga = .true. |
---|
| 88 | ! initial profiles from TOGA-COARE IFA files |
---|
| 89 | ! LS convergence and SST imposed from TOGA-COARE IFA files |
---|
| 90 | ! = 3 ==> forcing_GCM2SCM = .true. |
---|
| 91 | ! initial profiles from the GCM output |
---|
| 92 | ! LS convergence imposed from the GCM output |
---|
| 93 | ! = 4 ==> forcing_twpi = .true. |
---|
| 94 | ! initial profiles from TWPICE nc files |
---|
| 95 | ! LS convergence and SST imposed from TWPICE nc files |
---|
| 96 | ! = 5 ==> forcing_rico = .true. |
---|
| 97 | ! initial profiles from RICO idealized |
---|
| 98 | ! LS convergence imposed from RICO (cst) |
---|
| 99 | ! = 6 ==> forcing_amma = .true. |
---|
[2191] | 100 | ! = 10 ==> forcing_case = .true. |
---|
| 101 | ! initial profiles from case.nc file |
---|
[2017] | 102 | ! = 40 ==> forcing_GCSSold = .true. |
---|
| 103 | ! initial profile from GCSS file |
---|
| 104 | ! LS convergence imposed from GCSS file |
---|
| 105 | ! = 50 ==> forcing_fire = .true. |
---|
| 106 | ! = 59 ==> forcing_sandu = .true. |
---|
| 107 | ! initial profiles from sanduref file: see prof.inp.001 |
---|
| 108 | ! SST varying with time and divergence constante: see ifa_sanduref.txt file |
---|
| 109 | ! Radiation has to be computed interactively |
---|
| 110 | ! = 60 ==> forcing_astex = .true. |
---|
| 111 | ! initial profiles from file: see prof.inp.001 |
---|
| 112 | ! SST,divergence,ug,vg,ufa,vfa varying with time : see ifa_astex.txt file |
---|
| 113 | ! Radiation has to be computed interactively |
---|
| 114 | ! = 61 ==> forcing_armcu = .true. |
---|
| 115 | ! initial profiles from file: see prof.inp.001 |
---|
| 116 | ! sensible and latent heat flux imposed: see ifa_arm_cu_1.txt |
---|
| 117 | ! large scale advective forcing & radiative tendencies applied below 1000m: see ifa_arm_cu_2.txt |
---|
| 118 | ! use geostrophic wind ug=10m/s vg=0m/s. Duration of the case 53100s |
---|
| 119 | ! Radiation to be switched off |
---|
[2716] | 120 | ! > 100 ==> forcing_case = .true. or forcing_case2 = .true. |
---|
| 121 | ! initial profiles from case.nc file |
---|
[2017] | 122 | ! |
---|
| 123 | forcing_type = 0 |
---|
| 124 | CALL getin('forcing_type',forcing_type) |
---|
| 125 | imp_fcg_gcssold = .false. |
---|
| 126 | ts_fcg_gcssold = .false. |
---|
| 127 | Tp_fcg_gcssold = .false. |
---|
| 128 | Tp_ini_gcssold = .false. |
---|
| 129 | xTurb_fcg_gcssold = .false. |
---|
| 130 | IF (forcing_type .eq.40) THEN |
---|
| 131 | CALL getin('imp_fcg',imp_fcg_gcssold) |
---|
| 132 | CALL getin('ts_fcg',ts_fcg_gcssold) |
---|
| 133 | CALL getin('tp_fcg',Tp_fcg_gcssold) |
---|
| 134 | CALL getin('tp_ini',Tp_ini_gcssold) |
---|
| 135 | CALL getin('turb_fcg',xTurb_fcg_gcssold) |
---|
| 136 | ENDIF |
---|
| 137 | |
---|
[2716] | 138 | !Parametres de forcage |
---|
[2191] | 139 | !Config Key = tend_t |
---|
| 140 | !Config Desc = forcage ou non par advection de T |
---|
| 141 | !Config Def = false |
---|
| 142 | !Config Help = forcage ou non par advection de T |
---|
| 143 | tend_t =0 |
---|
| 144 | CALL getin('tend_t',tend_t) |
---|
| 145 | |
---|
| 146 | !Config Key = tend_q |
---|
| 147 | !Config Desc = forcage ou non par advection de q |
---|
| 148 | !Config Def = false |
---|
| 149 | !Config Help = forcage ou non par advection de q |
---|
| 150 | tend_q =0 |
---|
| 151 | CALL getin('tend_q',tend_q) |
---|
| 152 | |
---|
| 153 | !Config Key = tend_u |
---|
| 154 | !Config Desc = forcage ou non par advection de u |
---|
| 155 | !Config Def = false |
---|
| 156 | !Config Help = forcage ou non par advection de u |
---|
| 157 | tend_u =0 |
---|
| 158 | CALL getin('tend_u',tend_u) |
---|
| 159 | |
---|
| 160 | !Config Key = tend_v |
---|
| 161 | !Config Desc = forcage ou non par advection de v |
---|
| 162 | !Config Def = false |
---|
| 163 | !Config Help = forcage ou non par advection de v |
---|
| 164 | tend_v =0 |
---|
| 165 | CALL getin('tend_v',tend_v) |
---|
| 166 | |
---|
| 167 | !Config Key = tend_w |
---|
| 168 | !Config Desc = forcage ou non par vitesse verticale |
---|
| 169 | !Config Def = false |
---|
| 170 | !Config Help = forcage ou non par vitesse verticale |
---|
| 171 | tend_w =0 |
---|
| 172 | CALL getin('tend_w',tend_w) |
---|
| 173 | |
---|
| 174 | !Config Key = tend_rayo |
---|
| 175 | !Config Desc = forcage ou non par dtrad |
---|
| 176 | !Config Def = false |
---|
| 177 | !Config Help = forcage ou non par dtrad |
---|
| 178 | tend_rayo =0 |
---|
| 179 | CALL getin('tend_rayo',tend_rayo) |
---|
| 180 | |
---|
| 181 | |
---|
| 182 | !Config Key = nudge_t |
---|
| 183 | !Config Desc = constante de nudging de T |
---|
| 184 | !Config Def = false |
---|
| 185 | !Config Help = constante de nudging de T |
---|
| 186 | nudge_t =0. |
---|
| 187 | CALL getin('nudge_t',nudge_t) |
---|
| 188 | |
---|
| 189 | !Config Key = nudge_q |
---|
| 190 | !Config Desc = constante de nudging de q |
---|
| 191 | !Config Def = false |
---|
| 192 | !Config Help = constante de nudging de q |
---|
| 193 | nudge_q =0. |
---|
| 194 | CALL getin('nudge_q',nudge_q) |
---|
| 195 | |
---|
| 196 | !Config Key = nudge_u |
---|
| 197 | !Config Desc = constante de nudging de u |
---|
| 198 | !Config Def = false |
---|
| 199 | !Config Help = constante de nudging de u |
---|
| 200 | nudge_u =0. |
---|
| 201 | CALL getin('nudge_u',nudge_u) |
---|
| 202 | |
---|
| 203 | !Config Key = nudge_v |
---|
| 204 | !Config Desc = constante de nudging de v |
---|
| 205 | !Config Def = false |
---|
| 206 | !Config Help = constante de nudging de v |
---|
| 207 | nudge_v =0. |
---|
| 208 | CALL getin('nudge_v',nudge_v) |
---|
| 209 | |
---|
| 210 | !Config Key = nudge_w |
---|
| 211 | !Config Desc = constante de nudging de w |
---|
| 212 | !Config Def = false |
---|
| 213 | !Config Help = constante de nudging de w |
---|
| 214 | nudge_w =0. |
---|
| 215 | CALL getin('nudge_w',nudge_w) |
---|
| 216 | |
---|
| 217 | |
---|
[2181] | 218 | !Config Key = iflag_nudge |
---|
| 219 | !Config Desc = atmospheric nudging ttype (decimal code) |
---|
| 220 | !Config Def = 0 |
---|
| 221 | !Config Help = 0 ==> no nudging |
---|
| 222 | ! If digit number n of iflag_nudge is set, then nudging of type n is on |
---|
| 223 | ! If digit number n of iflag_nudge is not set, then nudging of type n is off |
---|
| 224 | ! (digits are numbered from the right) |
---|
| 225 | iflag_nudge = 0 |
---|
| 226 | CALL getin('iflag_nudge',iflag_nudge) |
---|
| 227 | |
---|
[2017] | 228 | !Config Key = ok_flux_surf |
---|
| 229 | !Config Desc = forcage ou non par les flux de surface |
---|
| 230 | !Config Def = false |
---|
| 231 | !Config Help = forcage ou non par les flux de surface |
---|
[2019] | 232 | ok_flux_surf =.false. |
---|
[2017] | 233 | CALL getin('ok_flux_surf',ok_flux_surf) |
---|
| 234 | |
---|
[3780] | 235 | !Config Key = ok_forc_tsurf |
---|
| 236 | !Config Desc = forcage ou non par la Ts |
---|
| 237 | !Config Def = false |
---|
| 238 | !Config Help = forcage ou non par la Ts |
---|
| 239 | ok_forc_tsurf=.false. |
---|
| 240 | CALL getin('ok_forc_tsurf',ok_forc_tsurf) |
---|
| 241 | |
---|
[2126] | 242 | !Config Key = ok_prescr_ust |
---|
| 243 | !Config Desc = ustar impose ou non |
---|
| 244 | !Config Def = false |
---|
| 245 | !Config Help = ustar impose ou non |
---|
| 246 | ok_prescr_ust = .false. |
---|
| 247 | CALL getin('ok_prescr_ust',ok_prescr_ust) |
---|
| 248 | |
---|
[3780] | 249 | |
---|
| 250 | !Config Key = ok_prescr_beta |
---|
| 251 | !Config Desc = betaevap impose ou non |
---|
| 252 | !Config Def = false |
---|
| 253 | !Config Help = betaevap impose ou non |
---|
| 254 | ok_prescr_beta = .false. |
---|
| 255 | CALL getin('ok_prescr_beta',ok_prescr_beta) |
---|
| 256 | |
---|
[2017] | 257 | !Config Key = ok_old_disvert |
---|
| 258 | !Config Desc = utilisation de l ancien programme disvert0 (dans 1DUTILS.h) |
---|
| 259 | !Config Def = false |
---|
| 260 | !Config Help = utilisation de l ancien programme disvert0 (dans 1DUTILS.h) |
---|
| 261 | ok_old_disvert = .FALSE. |
---|
| 262 | CALL getin('ok_old_disvert',ok_old_disvert) |
---|
| 263 | |
---|
| 264 | !Config Key = time_ini |
---|
| 265 | !Config Desc = meaningless in this case |
---|
| 266 | !Config Def = 0. |
---|
| 267 | !Config Help = |
---|
[4650] | 268 | time_ini = 0. |
---|
[2017] | 269 | CALL getin('time_ini',time_ini) |
---|
| 270 | |
---|
| 271 | !Config Key = rlat et rlon |
---|
| 272 | !Config Desc = latitude et longitude |
---|
| 273 | !Config Def = 0.0 0.0 |
---|
| 274 | !Config Help = fixe la position de la colonne |
---|
| 275 | xlat = 0. |
---|
| 276 | xlon = 0. |
---|
| 277 | CALL getin('rlat',xlat) |
---|
| 278 | CALL getin('rlon',xlon) |
---|
| 279 | |
---|
| 280 | !Config Key = airephy |
---|
| 281 | !Config Desc = Grid cell area |
---|
| 282 | !Config Def = 1.e11 |
---|
| 283 | !Config Help = |
---|
| 284 | airefi = 1.e11 |
---|
| 285 | CALL getin('airephy',airefi) |
---|
| 286 | |
---|
| 287 | !Config Key = nat_surf |
---|
| 288 | !Config Desc = surface type |
---|
| 289 | !Config Def = 0 (ocean) |
---|
| 290 | !Config Help = 0=ocean,1=land,2=glacier,3=banquise |
---|
| 291 | nat_surf = 0. |
---|
| 292 | CALL getin('nat_surf',nat_surf) |
---|
| 293 | |
---|
| 294 | !Config Key = tsurf |
---|
| 295 | !Config Desc = surface temperature |
---|
| 296 | !Config Def = 290. |
---|
[3780] | 297 | !Config Help = surface temperature |
---|
[2017] | 298 | tsurf = 290. |
---|
| 299 | CALL getin('tsurf',tsurf) |
---|
| 300 | |
---|
| 301 | !Config Key = psurf |
---|
| 302 | !Config Desc = surface pressure |
---|
| 303 | !Config Def = 102400. |
---|
| 304 | !Config Help = |
---|
| 305 | psurf = 102400. |
---|
| 306 | CALL getin('psurf',psurf) |
---|
| 307 | |
---|
| 308 | !Config Key = zsurf |
---|
| 309 | !Config Desc = surface altitude |
---|
| 310 | !Config Def = 0. |
---|
| 311 | !Config Help = |
---|
| 312 | zsurf = 0. |
---|
| 313 | CALL getin('zsurf',zsurf) |
---|
[3780] | 314 | ! EV pour accord avec format standard |
---|
| 315 | CALL getin('zorog',zsurf) |
---|
[2017] | 316 | |
---|
[3780] | 317 | |
---|
[2017] | 318 | !Config Key = rugos |
---|
| 319 | !Config Desc = coefficient de frottement |
---|
| 320 | !Config Def = 0.0001 |
---|
| 321 | !Config Help = calcul du Cdrag |
---|
| 322 | rugos = 0.0001 |
---|
| 323 | CALL getin('rugos',rugos) |
---|
[3659] | 324 | ! FH/2020/04/08/confinement: Pour le nouveau format standard, la rugosite s'appelle z0 |
---|
| 325 | CALL getin('z0',rugos) |
---|
[2017] | 326 | |
---|
[2672] | 327 | !Config Key = rugosh |
---|
| 328 | !Config Desc = coefficient de frottement |
---|
| 329 | !Config Def = rugos |
---|
| 330 | !Config Help = calcul du Cdrag |
---|
| 331 | rugosh = rugos |
---|
| 332 | CALL getin('rugosh',rugosh) |
---|
| 333 | |
---|
| 334 | |
---|
| 335 | |
---|
| 336 | !Config Key = snowmass |
---|
| 337 | !Config Desc = mass de neige de la surface en kg/m2 |
---|
| 338 | !Config Def = 0.0000 |
---|
| 339 | !Config Help = snowmass |
---|
| 340 | snowmass = 0.0000 |
---|
| 341 | CALL getin('snowmass',snowmass) |
---|
| 342 | |
---|
[2017] | 343 | !Config Key = wtsurf et wqsurf |
---|
| 344 | !Config Desc = ??? |
---|
| 345 | !Config Def = 0.0 0.0 |
---|
| 346 | !Config Help = |
---|
| 347 | wtsurf = 0.0 |
---|
| 348 | wqsurf = 0.0 |
---|
| 349 | CALL getin('wtsurf',wtsurf) |
---|
| 350 | CALL getin('wqsurf',wqsurf) |
---|
| 351 | |
---|
| 352 | !Config Key = albedo |
---|
| 353 | !Config Desc = albedo |
---|
| 354 | !Config Def = 0.09 |
---|
| 355 | !Config Help = |
---|
| 356 | albedo = 0.09 |
---|
| 357 | CALL getin('albedo',albedo) |
---|
| 358 | |
---|
| 359 | !Config Key = agesno |
---|
| 360 | !Config Desc = age de la neige |
---|
| 361 | !Config Def = 30.0 |
---|
| 362 | !Config Help = |
---|
| 363 | xagesno = 30.0 |
---|
| 364 | CALL getin('agesno',xagesno) |
---|
| 365 | |
---|
| 366 | !Config Key = restart_runoff |
---|
| 367 | !Config Desc = age de la neige |
---|
| 368 | !Config Def = 30.0 |
---|
| 369 | !Config Help = |
---|
| 370 | restart_runoff = 0.0 |
---|
| 371 | CALL getin('restart_runoff',restart_runoff) |
---|
| 372 | |
---|
| 373 | !Config Key = qsolinp |
---|
| 374 | !Config Desc = initial bucket water content (kg/m2) when land (5std) |
---|
| 375 | !Config Def = 30.0 |
---|
| 376 | !Config Help = |
---|
| 377 | qsolinp = 1. |
---|
| 378 | CALL getin('qsolinp',qsolinp) |
---|
| 379 | |
---|
[3780] | 380 | |
---|
| 381 | |
---|
| 382 | !Config Key = betaevap |
---|
| 383 | !Config Desc = beta for actual evaporation when prescribed |
---|
| 384 | !Config Def = 1.0 |
---|
| 385 | !Config Help = |
---|
| 386 | betaevap = 1. |
---|
| 387 | CALL getin('betaevap',betaevap) |
---|
| 388 | |
---|
[2017] | 389 | !Config Key = zpicinp |
---|
| 390 | !Config Desc = denivellation orographie |
---|
[2668] | 391 | !Config Def = 0. |
---|
[2017] | 392 | !Config Help = input brise |
---|
[2668] | 393 | zpicinp = 0. |
---|
[2017] | 394 | CALL getin('zpicinp',zpicinp) |
---|
| 395 | !Config key = nudge_tsoil |
---|
| 396 | !Config Desc = activation of soil temperature nudging |
---|
| 397 | !Config Def = .FALSE. |
---|
| 398 | !Config Help = ... |
---|
| 399 | |
---|
| 400 | nudge_tsoil=.FALSE. |
---|
| 401 | CALL getin('nudge_tsoil',nudge_tsoil) |
---|
| 402 | |
---|
| 403 | !Config key = isoil_nudge |
---|
| 404 | !Config Desc = level number where soil temperature is nudged |
---|
| 405 | !Config Def = 3 |
---|
| 406 | !Config Help = ... |
---|
| 407 | |
---|
| 408 | isoil_nudge=3 |
---|
| 409 | CALL getin('isoil_nudge',isoil_nudge) |
---|
| 410 | |
---|
| 411 | !Config key = Tsoil_nudge |
---|
| 412 | !Config Desc = target temperature for tsoil(isoil_nudge) |
---|
| 413 | !Config Def = 300. |
---|
| 414 | !Config Help = ... |
---|
| 415 | |
---|
| 416 | Tsoil_nudge=300. |
---|
| 417 | CALL getin('Tsoil_nudge',Tsoil_nudge) |
---|
| 418 | |
---|
| 419 | !Config key = tau_soil_nudge |
---|
| 420 | !Config Desc = nudging relaxation time for tsoil |
---|
| 421 | !Config Def = 3600. |
---|
| 422 | !Config Help = ... |
---|
| 423 | |
---|
| 424 | tau_soil_nudge=3600. |
---|
| 425 | CALL getin('tau_soil_nudge',tau_soil_nudge) |
---|
| 426 | |
---|
[2716] | 427 | !---------------------------------------------------------- |
---|
| 428 | ! Param??tres de for??age pour les forcages communs: |
---|
| 429 | ! Pour les forcages communs: ces entiers valent 0 ou 1 |
---|
| 430 | ! tadv= advection tempe, tadvv= adv tempe verticale, tadvh= adv tempe horizontale |
---|
| 431 | ! qadv= advection q, qadvv= adv q verticale, qadvh= adv q horizontale |
---|
| 432 | ! trad= 0 (rayonnement actif) ou 1 (prescrit par tend_rad) ou adv (prescir et contenu dans les tadv) |
---|
| 433 | ! forcages en omega, w, vent geostrophique ou ustar |
---|
| 434 | ! Parametres de nudging en u,v,t,q valent 0 ou 1 ou le temps de nudging |
---|
| 435 | !---------------------------------------------------------- |
---|
[2017] | 436 | |
---|
[2716] | 437 | !Config Key = tadv |
---|
| 438 | !Config Desc = forcage ou non par advection totale de T |
---|
| 439 | !Config Def = false |
---|
| 440 | !Config Help = forcage ou non par advection totale de T |
---|
| 441 | tadv =0 |
---|
| 442 | CALL getin('tadv',tadv) |
---|
[2017] | 443 | |
---|
[2716] | 444 | !Config Key = tadvv |
---|
| 445 | !Config Desc = forcage ou non par advection verticale de T |
---|
| 446 | !Config Def = false |
---|
| 447 | !Config Help = forcage ou non par advection verticale de T |
---|
| 448 | tadvv =0 |
---|
| 449 | CALL getin('tadvv',tadvv) |
---|
[2017] | 450 | |
---|
[2716] | 451 | !Config Key = tadvh |
---|
| 452 | !Config Desc = forcage ou non par advection horizontale de T |
---|
| 453 | !Config Def = false |
---|
| 454 | !Config Help = forcage ou non par advection horizontale de T |
---|
| 455 | tadvh =0 |
---|
| 456 | CALL getin('tadvh',tadvh) |
---|
| 457 | |
---|
| 458 | !Config Key = thadv |
---|
| 459 | !Config Desc = forcage ou non par advection totale de Theta |
---|
| 460 | !Config Def = false |
---|
| 461 | !Config Help = forcage ou non par advection totale de Theta |
---|
| 462 | thadv =0 |
---|
| 463 | CALL getin('thadv',thadv) |
---|
| 464 | |
---|
| 465 | !Config Key = thadvv |
---|
| 466 | !Config Desc = forcage ou non par advection verticale de Theta |
---|
| 467 | !Config Def = false |
---|
| 468 | !Config Help = forcage ou non par advection verticale de Theta |
---|
| 469 | thadvv =0 |
---|
| 470 | CALL getin('thadvv',thadvv) |
---|
| 471 | |
---|
| 472 | !Config Key = thadvh |
---|
| 473 | !Config Desc = forcage ou non par advection horizontale de Theta |
---|
| 474 | !Config Def = false |
---|
| 475 | !Config Help = forcage ou non par advection horizontale de Theta |
---|
| 476 | thadvh =0 |
---|
| 477 | CALL getin('thadvh',thadvh) |
---|
| 478 | |
---|
| 479 | !Config Key = qadv |
---|
| 480 | !Config Desc = forcage ou non par advection totale de Q |
---|
| 481 | !Config Def = false |
---|
| 482 | !Config Help = forcage ou non par advection totale de Q |
---|
| 483 | qadv =0 |
---|
| 484 | CALL getin('qadv',qadv) |
---|
| 485 | |
---|
| 486 | !Config Key = qadvv |
---|
| 487 | !Config Desc = forcage ou non par advection verticale de Q |
---|
| 488 | !Config Def = false |
---|
| 489 | !Config Help = forcage ou non par advection verticale de Q |
---|
| 490 | qadvv =0 |
---|
| 491 | CALL getin('qadvv',qadvv) |
---|
| 492 | |
---|
| 493 | !Config Key = qadvh |
---|
| 494 | !Config Desc = forcage ou non par advection horizontale de Q |
---|
| 495 | !Config Def = false |
---|
| 496 | !Config Help = forcage ou non par advection horizontale de Q |
---|
| 497 | qadvh =0 |
---|
| 498 | CALL getin('qadvh',qadvh) |
---|
| 499 | |
---|
| 500 | !Config Key = trad |
---|
| 501 | !Config Desc = forcage ou non par tendance radiative |
---|
| 502 | !Config Def = false |
---|
| 503 | !Config Help = forcage ou non par tendance radiative |
---|
| 504 | trad =0 |
---|
| 505 | CALL getin('trad',trad) |
---|
| 506 | |
---|
| 507 | !Config Key = forc_omega |
---|
| 508 | !Config Desc = forcage ou non par omega |
---|
| 509 | !Config Def = false |
---|
| 510 | !Config Help = forcage ou non par omega |
---|
| 511 | forc_omega =0 |
---|
| 512 | CALL getin('forc_omega',forc_omega) |
---|
| 513 | |
---|
[2920] | 514 | !Config Key = forc_u |
---|
| 515 | !Config Desc = forcage ou non par u |
---|
| 516 | !Config Def = false |
---|
| 517 | !Config Help = forcage ou non par u |
---|
| 518 | forc_u =0 |
---|
| 519 | CALL getin('forc_u',forc_u) |
---|
| 520 | |
---|
| 521 | !Config Key = forc_v |
---|
| 522 | !Config Desc = forcage ou non par v |
---|
| 523 | !Config Def = false |
---|
| 524 | !Config Help = forcage ou non par v |
---|
| 525 | forc_v =0 |
---|
| 526 | CALL getin('forc_v',forc_v) |
---|
[2716] | 527 | !Config Key = forc_w |
---|
| 528 | !Config Desc = forcage ou non par w |
---|
| 529 | !Config Def = false |
---|
| 530 | !Config Help = forcage ou non par w |
---|
| 531 | forc_w =0 |
---|
| 532 | CALL getin('forc_w',forc_w) |
---|
| 533 | |
---|
| 534 | !Config Key = forc_geo |
---|
| 535 | !Config Desc = forcage ou non par geo |
---|
| 536 | !Config Def = false |
---|
| 537 | !Config Help = forcage ou non par geo |
---|
| 538 | forc_geo =0 |
---|
| 539 | CALL getin('forc_geo',forc_geo) |
---|
| 540 | |
---|
| 541 | ! Meme chose que ok_precr_ust |
---|
| 542 | !Config Key = forc_ustar |
---|
| 543 | !Config Desc = forcage ou non par ustar |
---|
| 544 | !Config Def = false |
---|
| 545 | !Config Help = forcage ou non par ustar |
---|
| 546 | forc_ustar =0 |
---|
| 547 | CALL getin('forc_ustar',forc_ustar) |
---|
| 548 | IF (forc_ustar .EQ. 1) ok_prescr_ust=.true. |
---|
| 549 | |
---|
[3780] | 550 | |
---|
[2716] | 551 | !Config Key = nudging_u |
---|
| 552 | !Config Desc = forcage ou non par nudging sur u |
---|
| 553 | !Config Def = false |
---|
| 554 | !Config Help = forcage ou non par nudging sur u |
---|
| 555 | nudging_u =0 |
---|
| 556 | CALL getin('nudging_u',nudging_u) |
---|
| 557 | |
---|
| 558 | !Config Key = nudging_v |
---|
| 559 | !Config Desc = forcage ou non par nudging sur v |
---|
| 560 | !Config Def = false |
---|
| 561 | !Config Help = forcage ou non par nudging sur v |
---|
| 562 | nudging_v =0 |
---|
| 563 | CALL getin('nudging_v',nudging_v) |
---|
| 564 | |
---|
| 565 | !Config Key = nudging_w |
---|
| 566 | !Config Desc = forcage ou non par nudging sur w |
---|
| 567 | !Config Def = false |
---|
| 568 | !Config Help = forcage ou non par nudging sur w |
---|
| 569 | nudging_w =0 |
---|
| 570 | CALL getin('nudging_w',nudging_w) |
---|
| 571 | |
---|
[3593] | 572 | ! RELIQUE ANCIENS FORMAT. ECRASE PAR LE SUIVANT |
---|
[2716] | 573 | !Config Key = nudging_q |
---|
| 574 | !Config Desc = forcage ou non par nudging sur q |
---|
| 575 | !Config Def = false |
---|
| 576 | !Config Help = forcage ou non par nudging sur q |
---|
[3593] | 577 | nudging_qv =0 |
---|
| 578 | CALL getin('nudging_q',nudging_qv) |
---|
| 579 | CALL getin('nudging_qv',nudging_qv) |
---|
[2716] | 580 | |
---|
[3594] | 581 | p_nudging_u=11000. |
---|
| 582 | p_nudging_v=11000. |
---|
| 583 | p_nudging_t=11000. |
---|
| 584 | p_nudging_qv=11000. |
---|
| 585 | CALL getin('p_nudging_u',p_nudging_u) |
---|
| 586 | CALL getin('p_nudging_v',p_nudging_v) |
---|
| 587 | CALL getin('p_nudging_t',p_nudging_t) |
---|
| 588 | CALL getin('p_nudging_qv',p_nudging_qv) |
---|
[3593] | 589 | |
---|
[2716] | 590 | !Config Key = nudging_t |
---|
| 591 | !Config Desc = forcage ou non par nudging sur t |
---|
| 592 | !Config Def = false |
---|
| 593 | !Config Help = forcage ou non par nudging sur t |
---|
| 594 | nudging_t =0 |
---|
| 595 | CALL getin('nudging_t',nudging_t) |
---|
| 596 | |
---|
| 597 | |
---|
| 598 | |
---|
[2017] | 599 | write(lunout,*)' +++++++++++++++++++++++++++++++++++++++' |
---|
| 600 | write(lunout,*)' Configuration des parametres du gcm1D: ' |
---|
| 601 | write(lunout,*)' +++++++++++++++++++++++++++++++++++++++' |
---|
| 602 | write(lunout,*)' restart = ', restart |
---|
| 603 | write(lunout,*)' forcing_type = ', forcing_type |
---|
| 604 | write(lunout,*)' time_ini = ', time_ini |
---|
| 605 | write(lunout,*)' rlat = ', xlat |
---|
| 606 | write(lunout,*)' rlon = ', xlon |
---|
| 607 | write(lunout,*)' airephy = ', airefi |
---|
| 608 | write(lunout,*)' nat_surf = ', nat_surf |
---|
| 609 | write(lunout,*)' tsurf = ', tsurf |
---|
| 610 | write(lunout,*)' psurf = ', psurf |
---|
| 611 | write(lunout,*)' zsurf = ', zsurf |
---|
| 612 | write(lunout,*)' rugos = ', rugos |
---|
[2672] | 613 | write(lunout,*)' snowmass=', snowmass |
---|
[2017] | 614 | write(lunout,*)' wtsurf = ', wtsurf |
---|
| 615 | write(lunout,*)' wqsurf = ', wqsurf |
---|
| 616 | write(lunout,*)' albedo = ', albedo |
---|
| 617 | write(lunout,*)' xagesno = ', xagesno |
---|
| 618 | write(lunout,*)' restart_runoff = ', restart_runoff |
---|
| 619 | write(lunout,*)' qsolinp = ', qsolinp |
---|
| 620 | write(lunout,*)' zpicinp = ', zpicinp |
---|
| 621 | write(lunout,*)' nudge_tsoil = ', nudge_tsoil |
---|
| 622 | write(lunout,*)' isoil_nudge = ', isoil_nudge |
---|
| 623 | write(lunout,*)' Tsoil_nudge = ', Tsoil_nudge |
---|
| 624 | write(lunout,*)' tau_soil_nudge = ', tau_soil_nudge |
---|
[2716] | 625 | write(lunout,*)' tadv = ', tadv |
---|
| 626 | write(lunout,*)' tadvv = ', tadvv |
---|
| 627 | write(lunout,*)' tadvh = ', tadvh |
---|
| 628 | write(lunout,*)' thadv = ', thadv |
---|
| 629 | write(lunout,*)' thadvv = ', thadvv |
---|
| 630 | write(lunout,*)' thadvh = ', thadvh |
---|
| 631 | write(lunout,*)' qadv = ', qadv |
---|
| 632 | write(lunout,*)' qadvv = ', qadvv |
---|
| 633 | write(lunout,*)' qadvh = ', qadvh |
---|
| 634 | write(lunout,*)' trad = ', trad |
---|
| 635 | write(lunout,*)' forc_omega = ', forc_omega |
---|
| 636 | write(lunout,*)' forc_w = ', forc_w |
---|
| 637 | write(lunout,*)' forc_geo = ', forc_geo |
---|
| 638 | write(lunout,*)' forc_ustar = ', forc_ustar |
---|
| 639 | write(lunout,*)' nudging_u = ', nudging_u |
---|
| 640 | write(lunout,*)' nudging_v = ', nudging_v |
---|
| 641 | write(lunout,*)' nudging_t = ', nudging_t |
---|
[3593] | 642 | write(lunout,*)' nudging_qv = ', nudging_qv |
---|
[2017] | 643 | IF (forcing_type .eq.40) THEN |
---|
| 644 | write(lunout,*) '--- Forcing type GCSS Old --- with:' |
---|
| 645 | write(lunout,*)'imp_fcg',imp_fcg_gcssold |
---|
| 646 | write(lunout,*)'ts_fcg',ts_fcg_gcssold |
---|
| 647 | write(lunout,*)'tp_fcg',Tp_fcg_gcssold |
---|
| 648 | write(lunout,*)'tp_ini',Tp_ini_gcssold |
---|
| 649 | write(lunout,*)'xturb_fcg',xTurb_fcg_gcssold |
---|
| 650 | ENDIF |
---|
| 651 | |
---|
| 652 | write(lunout,*)' +++++++++++++++++++++++++++++++++++++++' |
---|
| 653 | write(lunout,*) |
---|
[2019] | 654 | ! |
---|
[2017] | 655 | RETURN |
---|
| 656 | END |
---|
| 657 | ! |
---|
| 658 | ! $Id: dyn1deta0.F 1279 2010/07/30 A Lahellec$ |
---|
| 659 | ! |
---|
[2019] | 660 | ! |
---|
| 661 | SUBROUTINE dyn1deta0(fichnom,plev,play,phi,phis,presnivs, & |
---|
[2017] | 662 | & ucov,vcov,temp,q,omega2) |
---|
| 663 | USE dimphy |
---|
| 664 | USE mod_grid_phy_lmdz |
---|
| 665 | USE mod_phys_lmdz_para |
---|
| 666 | USE iophy |
---|
| 667 | USE phys_state_var_mod |
---|
| 668 | USE iostart |
---|
| 669 | USE write_field_phy |
---|
| 670 | USE infotrac |
---|
| 671 | use control_mod |
---|
[2597] | 672 | USE comconst_mod, ONLY: im, jm, lllm |
---|
[2603] | 673 | USE logic_mod, ONLY: fxyhypb, ysinus |
---|
[2601] | 674 | USE temps_mod, ONLY: annee_ref, day_ini, day_ref, itau_dyn |
---|
[2017] | 675 | |
---|
| 676 | IMPLICIT NONE |
---|
[2019] | 677 | !======================================================= |
---|
| 678 | ! Ecriture du fichier de redemarrage sous format NetCDF |
---|
| 679 | !======================================================= |
---|
| 680 | ! Declarations: |
---|
| 681 | ! ------------- |
---|
[2239] | 682 | include "dimensions.h" |
---|
[4603] | 683 | !!#include "control.h" |
---|
[5084] | 684 | include "netcdf.inc" |
---|
[2017] | 685 | |
---|
[2019] | 686 | ! Arguments: |
---|
| 687 | ! ---------- |
---|
[2017] | 688 | CHARACTER*(*) fichnom |
---|
[2019] | 689 | !Al1 plev tronque pour .nc mais plev(klev+1):=0 |
---|
[2139] | 690 | real :: plev(klon,klev+1),play (klon,klev),phi(klon,klev) |
---|
[2017] | 691 | real :: presnivs(klon,klev) |
---|
| 692 | real :: ucov(klon,klev),vcov(klon,klev),temp(klon,klev) |
---|
| 693 | real :: q(klon,klev,nqtot),omega2(klon,klev) |
---|
[2019] | 694 | ! real :: ug(klev),vg(klev),fcoriolis |
---|
[2017] | 695 | real :: phis(klon) |
---|
| 696 | |
---|
[2019] | 697 | ! Variables locales pour NetCDF: |
---|
| 698 | ! ------------------------------ |
---|
[2017] | 699 | INTEGER iq |
---|
| 700 | INTEGER length |
---|
| 701 | PARAMETER (length = 100) |
---|
| 702 | REAL tab_cntrl(length) ! tableau des parametres du run |
---|
| 703 | character*4 nmq(nqtot) |
---|
| 704 | character*12 modname |
---|
| 705 | character*80 abort_message |
---|
| 706 | LOGICAL found |
---|
| 707 | |
---|
| 708 | modname = 'dyn1deta0 : ' |
---|
[2933] | 709 | !! nmq(1)="vap" |
---|
| 710 | !! nmq(2)="cond" |
---|
| 711 | !! do iq=3,nqtot |
---|
| 712 | !! write(nmq(iq),'("tra",i1)') iq-2 |
---|
| 713 | !! enddo |
---|
| 714 | DO iq = 1,nqtot |
---|
[4046] | 715 | nmq(iq) = trim(tracers(iq)%name) |
---|
[2933] | 716 | ENDDO |
---|
[2017] | 717 | print*,'in dyn1deta0 ',fichnom,klon,klev,nqtot |
---|
| 718 | CALL open_startphy(fichnom) |
---|
| 719 | print*,'after open startphy ',fichnom,nmq |
---|
| 720 | |
---|
[2019] | 721 | ! |
---|
| 722 | ! Lecture des parametres de controle: |
---|
| 723 | ! |
---|
[2017] | 724 | CALL get_var("controle",tab_cntrl) |
---|
| 725 | |
---|
| 726 | |
---|
| 727 | im = tab_cntrl(1) |
---|
| 728 | jm = tab_cntrl(2) |
---|
| 729 | lllm = tab_cntrl(3) |
---|
| 730 | day_ref = tab_cntrl(4) |
---|
| 731 | annee_ref = tab_cntrl(5) |
---|
[2019] | 732 | ! rad = tab_cntrl(6) |
---|
| 733 | ! omeg = tab_cntrl(7) |
---|
| 734 | ! g = tab_cntrl(8) |
---|
| 735 | ! cpp = tab_cntrl(9) |
---|
| 736 | ! kappa = tab_cntrl(10) |
---|
| 737 | ! daysec = tab_cntrl(11) |
---|
| 738 | ! dtvr = tab_cntrl(12) |
---|
| 739 | ! etot0 = tab_cntrl(13) |
---|
| 740 | ! ptot0 = tab_cntrl(14) |
---|
| 741 | ! ztot0 = tab_cntrl(15) |
---|
| 742 | ! stot0 = tab_cntrl(16) |
---|
| 743 | ! ang0 = tab_cntrl(17) |
---|
| 744 | ! pa = tab_cntrl(18) |
---|
| 745 | ! preff = tab_cntrl(19) |
---|
| 746 | ! |
---|
| 747 | ! clon = tab_cntrl(20) |
---|
| 748 | ! clat = tab_cntrl(21) |
---|
| 749 | ! grossismx = tab_cntrl(22) |
---|
| 750 | ! grossismy = tab_cntrl(23) |
---|
| 751 | ! |
---|
[2017] | 752 | IF ( tab_cntrl(24).EQ.1. ) THEN |
---|
[2019] | 753 | fxyhypb =.true. |
---|
| 754 | ! dzoomx = tab_cntrl(25) |
---|
| 755 | ! dzoomy = tab_cntrl(26) |
---|
| 756 | ! taux = tab_cntrl(28) |
---|
| 757 | ! tauy = tab_cntrl(29) |
---|
[2017] | 758 | ELSE |
---|
[2019] | 759 | fxyhypb = .false. |
---|
| 760 | ysinus = .false. |
---|
| 761 | IF( tab_cntrl(27).EQ.1. ) ysinus =.true. |
---|
[2017] | 762 | ENDIF |
---|
| 763 | |
---|
| 764 | day_ini = tab_cntrl(30) |
---|
| 765 | itau_dyn = tab_cntrl(31) |
---|
[2019] | 766 | ! ................................................................. |
---|
| 767 | ! |
---|
| 768 | ! |
---|
| 769 | ! PRINT*,'rad,omeg,g,cpp,kappa',rad,omeg,g,cpp,kappa |
---|
| 770 | !Al1 |
---|
| 771 | Print*,'day_ref,annee_ref,day_ini,itau_dyn', & |
---|
[2017] | 772 | & day_ref,annee_ref,day_ini,itau_dyn |
---|
| 773 | |
---|
[2019] | 774 | ! Lecture des champs |
---|
| 775 | ! |
---|
[2017] | 776 | CALL get_field("play",play,found) |
---|
| 777 | IF (.NOT. found) PRINT*, modname//'Le champ <Play> est absent' |
---|
| 778 | CALL get_field("phi",phi,found) |
---|
| 779 | IF (.NOT. found) PRINT*, modname//'Le champ <Phi> est absent' |
---|
| 780 | CALL get_field("phis",phis,found) |
---|
| 781 | IF (.NOT. found) PRINT*, modname//'Le champ <Phis> est absent' |
---|
| 782 | CALL get_field("presnivs",presnivs,found) |
---|
| 783 | IF (.NOT. found) PRINT*, modname//'Le champ <Presnivs> est absent' |
---|
| 784 | CALL get_field("ucov",ucov,found) |
---|
| 785 | IF (.NOT. found) PRINT*, modname//'Le champ <ucov> est absent' |
---|
| 786 | CALL get_field("vcov",vcov,found) |
---|
| 787 | IF (.NOT. found) PRINT*, modname//'Le champ <vcov> est absent' |
---|
| 788 | CALL get_field("temp",temp,found) |
---|
| 789 | IF (.NOT. found) PRINT*, modname//'Le champ <temp> est absent' |
---|
| 790 | CALL get_field("omega2",omega2,found) |
---|
| 791 | IF (.NOT. found) PRINT*, modname//'Le champ <omega2> est absent' |
---|
[2139] | 792 | plev(1,klev+1)=0. |
---|
| 793 | CALL get_field("plev",plev(:,1:klev),found) |
---|
| 794 | IF (.NOT. found) PRINT*, modname//'Le champ <Plev> est absent' |
---|
[2017] | 795 | |
---|
| 796 | Do iq=1,nqtot |
---|
| 797 | CALL get_field("q"//nmq(iq),q(:,:,iq),found) |
---|
[2019] | 798 | IF (.NOT.found)PRINT*, modname//'Le champ <q'//nmq//'> est absent' |
---|
[2017] | 799 | EndDo |
---|
| 800 | |
---|
| 801 | CALL close_startphy |
---|
[2019] | 802 | print*,' close startphy',fichnom,play(1,1),play(1,klev),temp(1,klev) |
---|
| 803 | ! |
---|
[2017] | 804 | RETURN |
---|
| 805 | END |
---|
| 806 | ! |
---|
| 807 | ! $Id: dyn1dredem.F 1279 2010/07/29 A Lahellec$ |
---|
| 808 | ! |
---|
[2019] | 809 | ! |
---|
| 810 | SUBROUTINE dyn1dredem(fichnom,plev,play,phi,phis,presnivs, & |
---|
[2017] | 811 | & ucov,vcov,temp,q,omega2) |
---|
| 812 | USE dimphy |
---|
| 813 | USE mod_grid_phy_lmdz |
---|
| 814 | USE mod_phys_lmdz_para |
---|
| 815 | USE phys_state_var_mod |
---|
| 816 | USE iostart |
---|
| 817 | USE infotrac |
---|
| 818 | use control_mod |
---|
[2597] | 819 | USE comconst_mod, ONLY: cpp, daysec, dtvr, g, kappa, omeg, rad |
---|
[2603] | 820 | USE logic_mod, ONLY: fxyhypb, ysinus |
---|
[2601] | 821 | USE temps_mod, ONLY: annee_ref,day_end,day_ref,itau_dyn,itaufin |
---|
[2017] | 822 | |
---|
| 823 | IMPLICIT NONE |
---|
[2019] | 824 | !======================================================= |
---|
| 825 | ! Ecriture du fichier de redemarrage sous format NetCDF |
---|
| 826 | !======================================================= |
---|
| 827 | ! Declarations: |
---|
| 828 | ! ------------- |
---|
[2239] | 829 | include "dimensions.h" |
---|
[4603] | 830 | !!#include "control.h" |
---|
[5084] | 831 | include "netcdf.inc" |
---|
[2017] | 832 | |
---|
[2019] | 833 | ! Arguments: |
---|
| 834 | ! ---------- |
---|
[2017] | 835 | CHARACTER*(*) fichnom |
---|
[2019] | 836 | !Al1 plev tronque pour .nc mais plev(klev+1):=0 |
---|
[2017] | 837 | real :: plev(klon,klev),play (klon,klev),phi(klon,klev) |
---|
| 838 | real :: presnivs(klon,klev) |
---|
| 839 | real :: ucov(klon,klev),vcov(klon,klev),temp(klon,klev) |
---|
| 840 | real :: q(klon,klev,nqtot) |
---|
| 841 | real :: omega2(klon,klev),rho(klon,klev+1) |
---|
[2019] | 842 | ! real :: ug(klev),vg(klev),fcoriolis |
---|
[2017] | 843 | real :: phis(klon) |
---|
| 844 | |
---|
[2019] | 845 | ! Variables locales pour NetCDF: |
---|
| 846 | ! ------------------------------ |
---|
[2017] | 847 | INTEGER nid |
---|
| 848 | INTEGER ierr |
---|
| 849 | INTEGER iq,l |
---|
| 850 | INTEGER length |
---|
| 851 | PARAMETER (length = 100) |
---|
| 852 | REAL tab_cntrl(length) ! tableau des parametres du run |
---|
| 853 | character*4 nmq(nqtot) |
---|
| 854 | character*20 modname |
---|
| 855 | character*80 abort_message |
---|
[2019] | 856 | ! |
---|
[3513] | 857 | INTEGER pass |
---|
[2017] | 858 | |
---|
| 859 | CALL open_restartphy(fichnom) |
---|
| 860 | print*,'redm1 ',fichnom,klon,klev,nqtot |
---|
[2933] | 861 | !! nmq(1)="vap" |
---|
| 862 | !! nmq(2)="cond" |
---|
| 863 | !! nmq(3)="tra1" |
---|
| 864 | !! nmq(4)="tra2" |
---|
| 865 | DO iq = 1,nqtot |
---|
[4046] | 866 | nmq(iq) = trim(tracers(iq)%name) |
---|
[2933] | 867 | ENDDO |
---|
[2017] | 868 | |
---|
[3540] | 869 | ! modname = 'dyn1dredem' |
---|
| 870 | ! ierr = NF_OPEN(fichnom, NF_WRITE, nid) |
---|
| 871 | ! IF (ierr .NE. NF_NOERR) THEN |
---|
| 872 | ! abort_message="Pb. d ouverture "//fichnom |
---|
| 873 | ! CALL abort_gcm('Modele 1D',abort_message,1) |
---|
| 874 | ! ENDIF |
---|
[2017] | 875 | |
---|
| 876 | DO l=1,length |
---|
| 877 | tab_cntrl(l) = 0. |
---|
| 878 | ENDDO |
---|
| 879 | tab_cntrl(1) = FLOAT(iim) |
---|
| 880 | tab_cntrl(2) = FLOAT(jjm) |
---|
| 881 | tab_cntrl(3) = FLOAT(llm) |
---|
| 882 | tab_cntrl(4) = FLOAT(day_ref) |
---|
| 883 | tab_cntrl(5) = FLOAT(annee_ref) |
---|
| 884 | tab_cntrl(6) = rad |
---|
| 885 | tab_cntrl(7) = omeg |
---|
| 886 | tab_cntrl(8) = g |
---|
| 887 | tab_cntrl(9) = cpp |
---|
| 888 | tab_cntrl(10) = kappa |
---|
| 889 | tab_cntrl(11) = daysec |
---|
| 890 | tab_cntrl(12) = dtvr |
---|
[2019] | 891 | ! tab_cntrl(13) = etot0 |
---|
| 892 | ! tab_cntrl(14) = ptot0 |
---|
| 893 | ! tab_cntrl(15) = ztot0 |
---|
| 894 | ! tab_cntrl(16) = stot0 |
---|
| 895 | ! tab_cntrl(17) = ang0 |
---|
| 896 | ! tab_cntrl(18) = pa |
---|
| 897 | ! tab_cntrl(19) = preff |
---|
| 898 | ! |
---|
| 899 | ! ..... parametres pour le zoom ...... |
---|
[2017] | 900 | |
---|
[2019] | 901 | ! tab_cntrl(20) = clon |
---|
| 902 | ! tab_cntrl(21) = clat |
---|
| 903 | ! tab_cntrl(22) = grossismx |
---|
| 904 | ! tab_cntrl(23) = grossismy |
---|
| 905 | ! |
---|
[2017] | 906 | IF ( fxyhypb ) THEN |
---|
| 907 | tab_cntrl(24) = 1. |
---|
[2019] | 908 | ! tab_cntrl(25) = dzoomx |
---|
| 909 | ! tab_cntrl(26) = dzoomy |
---|
[2017] | 910 | tab_cntrl(27) = 0. |
---|
[2019] | 911 | ! tab_cntrl(28) = taux |
---|
| 912 | ! tab_cntrl(29) = tauy |
---|
[2017] | 913 | ELSE |
---|
| 914 | tab_cntrl(24) = 0. |
---|
[2019] | 915 | ! tab_cntrl(25) = dzoomx |
---|
| 916 | ! tab_cntrl(26) = dzoomy |
---|
[2017] | 917 | tab_cntrl(27) = 0. |
---|
| 918 | tab_cntrl(28) = 0. |
---|
| 919 | tab_cntrl(29) = 0. |
---|
| 920 | IF( ysinus ) tab_cntrl(27) = 1. |
---|
| 921 | ENDIF |
---|
[2019] | 922 | !Al1 iday_end -> day_end |
---|
[2017] | 923 | tab_cntrl(30) = FLOAT(day_end) |
---|
| 924 | tab_cntrl(31) = FLOAT(itau_dyn + itaufin) |
---|
[2019] | 925 | ! |
---|
[3513] | 926 | DO pass=1,2 |
---|
| 927 | CALL put_var(pass,"controle","Param. de controle Dyn1D",tab_cntrl) |
---|
[2019] | 928 | ! |
---|
[2017] | 929 | |
---|
[2019] | 930 | ! Ecriture/extension de la coordonnee temps |
---|
[2017] | 931 | |
---|
| 932 | |
---|
[2019] | 933 | ! Ecriture des champs |
---|
| 934 | ! |
---|
[3513] | 935 | CALL put_field(pass,"plev","p interfaces sauf la nulle",plev) |
---|
| 936 | CALL put_field(pass,"play","",play) |
---|
| 937 | CALL put_field(pass,"phi","geopotentielle",phi) |
---|
| 938 | CALL put_field(pass,"phis","geopotentiell de surface",phis) |
---|
| 939 | CALL put_field(pass,"presnivs","",presnivs) |
---|
| 940 | CALL put_field(pass,"ucov","",ucov) |
---|
| 941 | CALL put_field(pass,"vcov","",vcov) |
---|
| 942 | CALL put_field(pass,"temp","",temp) |
---|
| 943 | CALL put_field(pass,"omega2","",omega2) |
---|
[2017] | 944 | |
---|
| 945 | Do iq=1,nqtot |
---|
[3513] | 946 | CALL put_field(pass,"q"//nmq(iq),"eau vap ou condens et traceurs", & |
---|
[2019] | 947 | & q(:,:,iq)) |
---|
[2017] | 948 | EndDo |
---|
[3513] | 949 | IF (pass==1) CALL enddef_restartphy |
---|
| 950 | IF (pass==2) CALL close_restartphy |
---|
[2017] | 951 | |
---|
[3513] | 952 | |
---|
| 953 | ENDDO |
---|
| 954 | |
---|
[2019] | 955 | ! |
---|
[2017] | 956 | RETURN |
---|
| 957 | END |
---|
| 958 | SUBROUTINE gr_fi_dyn(nfield,ngrid,im,jm,pfi,pdyn) |
---|
| 959 | IMPLICIT NONE |
---|
| 960 | !======================================================================= |
---|
| 961 | ! passage d'un champ de la grille scalaire a la grille physique |
---|
| 962 | !======================================================================= |
---|
| 963 | |
---|
| 964 | !----------------------------------------------------------------------- |
---|
| 965 | ! declarations: |
---|
| 966 | ! ------------- |
---|
| 967 | |
---|
| 968 | INTEGER im,jm,ngrid,nfield |
---|
| 969 | REAL pdyn(im,jm,nfield) |
---|
| 970 | REAL pfi(ngrid,nfield) |
---|
| 971 | |
---|
| 972 | INTEGER i,j,ifield,ig |
---|
| 973 | |
---|
| 974 | !----------------------------------------------------------------------- |
---|
| 975 | ! calcul: |
---|
| 976 | ! ------- |
---|
| 977 | |
---|
| 978 | DO ifield=1,nfield |
---|
| 979 | ! traitement des poles |
---|
| 980 | DO i=1,im |
---|
| 981 | pdyn(i,1,ifield)=pfi(1,ifield) |
---|
| 982 | pdyn(i,jm,ifield)=pfi(ngrid,ifield) |
---|
| 983 | ENDDO |
---|
| 984 | |
---|
| 985 | ! traitement des point normaux |
---|
| 986 | DO j=2,jm-1 |
---|
| 987 | ig=2+(j-2)*(im-1) |
---|
| 988 | CALL SCOPY(im-1,pfi(ig,ifield),1,pdyn(1,j,ifield),1) |
---|
| 989 | pdyn(im,j,ifield)=pdyn(1,j,ifield) |
---|
| 990 | ENDDO |
---|
| 991 | ENDDO |
---|
| 992 | |
---|
| 993 | RETURN |
---|
| 994 | END |
---|
| 995 | |
---|
| 996 | |
---|
| 997 | |
---|
| 998 | SUBROUTINE abort_gcm(modname, message, ierr) |
---|
| 999 | |
---|
| 1000 | USE IOIPSL |
---|
| 1001 | ! |
---|
| 1002 | ! Stops the simulation cleanly, closing files and printing various |
---|
| 1003 | ! comments |
---|
| 1004 | ! |
---|
| 1005 | ! Input: modname = name of calling program |
---|
| 1006 | ! message = stuff to print |
---|
| 1007 | ! ierr = severity of situation ( = 0 normal ) |
---|
| 1008 | |
---|
| 1009 | character(len=*) modname |
---|
| 1010 | integer ierr |
---|
| 1011 | character(len=*) message |
---|
| 1012 | |
---|
| 1013 | write(*,*) 'in abort_gcm' |
---|
| 1014 | call histclo |
---|
| 1015 | ! call histclo(2) |
---|
| 1016 | ! call histclo(3) |
---|
| 1017 | ! call histclo(4) |
---|
| 1018 | ! call histclo(5) |
---|
| 1019 | write(*,*) 'out of histclo' |
---|
| 1020 | write(*,*) 'Stopping in ', modname |
---|
| 1021 | write(*,*) 'Reason = ',message |
---|
| 1022 | call getin_dump |
---|
| 1023 | ! |
---|
| 1024 | if (ierr .eq. 0) then |
---|
| 1025 | write(*,*) 'Everything is cool' |
---|
| 1026 | else |
---|
| 1027 | write(*,*) 'Houston, we have a problem ', ierr |
---|
| 1028 | endif |
---|
| 1029 | STOP |
---|
| 1030 | END |
---|
| 1031 | REAL FUNCTION fq_sat(kelvin, millibar) |
---|
| 1032 | ! |
---|
| 1033 | IMPLICIT none |
---|
| 1034 | !====================================================================== |
---|
| 1035 | ! Autheur(s): Z.X. Li (LMD/CNRS) |
---|
| 1036 | ! Objet: calculer la vapeur d'eau saturante (formule Centre Euro.) |
---|
| 1037 | !====================================================================== |
---|
| 1038 | ! Arguments: |
---|
| 1039 | ! kelvin---input-R: temperature en Kelvin |
---|
| 1040 | ! millibar--input-R: pression en mb |
---|
| 1041 | ! |
---|
| 1042 | ! fq_sat----output-R: vapeur d'eau saturante en kg/kg |
---|
| 1043 | !====================================================================== |
---|
| 1044 | ! |
---|
| 1045 | REAL kelvin, millibar |
---|
| 1046 | ! |
---|
| 1047 | REAL r2es |
---|
| 1048 | PARAMETER (r2es=611.14 *18.0153/28.9644) |
---|
| 1049 | ! |
---|
| 1050 | REAL r3les, r3ies, r3es |
---|
| 1051 | PARAMETER (R3LES=17.269) |
---|
| 1052 | PARAMETER (R3IES=21.875) |
---|
| 1053 | ! |
---|
| 1054 | REAL r4les, r4ies, r4es |
---|
| 1055 | PARAMETER (R4LES=35.86) |
---|
| 1056 | PARAMETER (R4IES=7.66) |
---|
| 1057 | ! |
---|
| 1058 | REAL rtt |
---|
| 1059 | PARAMETER (rtt=273.16) |
---|
| 1060 | ! |
---|
| 1061 | REAL retv |
---|
| 1062 | PARAMETER (retv=28.9644/18.0153 - 1.0) |
---|
| 1063 | ! |
---|
| 1064 | REAL zqsat |
---|
| 1065 | REAL temp, pres |
---|
| 1066 | ! ------------------------------------------------------------------ |
---|
| 1067 | ! |
---|
| 1068 | ! |
---|
| 1069 | temp = kelvin |
---|
| 1070 | pres = millibar * 100.0 |
---|
| 1071 | ! write(*,*)'kelvin,millibar=',kelvin,millibar |
---|
| 1072 | ! write(*,*)'temp,pres=',temp,pres |
---|
| 1073 | ! |
---|
| 1074 | IF (temp .LE. rtt) THEN |
---|
| 1075 | r3es = r3ies |
---|
| 1076 | r4es = r4ies |
---|
| 1077 | ELSE |
---|
| 1078 | r3es = r3les |
---|
| 1079 | r4es = r4les |
---|
| 1080 | ENDIF |
---|
| 1081 | ! |
---|
| 1082 | zqsat=r2es/pres * EXP ( r3es*(temp-rtt) / (temp-r4es) ) |
---|
| 1083 | zqsat=MIN(0.5,ZQSAT) |
---|
| 1084 | zqsat=zqsat/(1.-retv *zqsat) |
---|
| 1085 | ! |
---|
| 1086 | fq_sat = zqsat |
---|
| 1087 | ! |
---|
| 1088 | RETURN |
---|
| 1089 | END |
---|
| 1090 | |
---|
| 1091 | SUBROUTINE gr_dyn_fi(nfield,im,jm,ngrid,pdyn,pfi) |
---|
| 1092 | IMPLICIT NONE |
---|
| 1093 | !======================================================================= |
---|
| 1094 | ! passage d'un champ de la grille scalaire a la grille physique |
---|
| 1095 | !======================================================================= |
---|
| 1096 | |
---|
| 1097 | !----------------------------------------------------------------------- |
---|
| 1098 | ! declarations: |
---|
| 1099 | ! ------------- |
---|
| 1100 | |
---|
| 1101 | INTEGER im,jm,ngrid,nfield |
---|
| 1102 | REAL pdyn(im,jm,nfield) |
---|
| 1103 | REAL pfi(ngrid,nfield) |
---|
| 1104 | |
---|
| 1105 | INTEGER j,ifield,ig |
---|
| 1106 | |
---|
| 1107 | !----------------------------------------------------------------------- |
---|
| 1108 | ! calcul: |
---|
| 1109 | ! ------- |
---|
| 1110 | |
---|
[2019] | 1111 | IF(ngrid.NE.2+(jm-2)*(im-1).AND.ngrid.NE.1) & |
---|
| 1112 | & STOP 'probleme de dim' |
---|
[2017] | 1113 | ! traitement des poles |
---|
| 1114 | CALL SCOPY(nfield,pdyn,im*jm,pfi,ngrid) |
---|
| 1115 | CALL SCOPY(nfield,pdyn(1,jm,1),im*jm,pfi(ngrid,1),ngrid) |
---|
| 1116 | |
---|
| 1117 | ! traitement des point normaux |
---|
| 1118 | DO ifield=1,nfield |
---|
| 1119 | DO j=2,jm-1 |
---|
| 1120 | ig=2+(j-2)*(im-1) |
---|
| 1121 | CALL SCOPY(im-1,pdyn(1,j,ifield),1,pfi(ig,ifield),1) |
---|
| 1122 | ENDDO |
---|
| 1123 | ENDDO |
---|
| 1124 | |
---|
| 1125 | RETURN |
---|
| 1126 | END |
---|
| 1127 | |
---|
| 1128 | SUBROUTINE disvert0(pa,preff,ap,bp,dpres,presnivs,nivsigs,nivsig) |
---|
| 1129 | |
---|
| 1130 | ! Ancienne version disvert dont on a modifie nom pour utiliser |
---|
| 1131 | ! le disvert de dyn3d (qui permet d'utiliser grille avec ab,bp imposes) |
---|
| 1132 | ! (MPL 18092012) |
---|
| 1133 | ! |
---|
| 1134 | ! Auteur : P. Le Van . |
---|
| 1135 | ! |
---|
| 1136 | IMPLICIT NONE |
---|
| 1137 | |
---|
[2239] | 1138 | include "dimensions.h" |
---|
| 1139 | include "paramet.h" |
---|
[2017] | 1140 | ! |
---|
| 1141 | !======================================================================= |
---|
| 1142 | ! |
---|
| 1143 | ! |
---|
| 1144 | ! s = sigma ** kappa : coordonnee verticale |
---|
| 1145 | ! dsig(l) : epaisseur de la couche l ds la coord. s |
---|
| 1146 | ! sig(l) : sigma a l'interface des couches l et l-1 |
---|
| 1147 | ! ds(l) : distance entre les couches l et l-1 en coord.s |
---|
| 1148 | ! |
---|
| 1149 | !======================================================================= |
---|
| 1150 | ! |
---|
| 1151 | REAL pa,preff |
---|
| 1152 | REAL ap(llmp1),bp(llmp1),dpres(llm),nivsigs(llm),nivsig(llmp1) |
---|
| 1153 | REAL presnivs(llm) |
---|
| 1154 | ! |
---|
| 1155 | ! declarations: |
---|
| 1156 | ! ------------- |
---|
| 1157 | ! |
---|
| 1158 | REAL sig(llm+1),dsig(llm) |
---|
| 1159 | ! |
---|
| 1160 | INTEGER l |
---|
| 1161 | REAL snorm |
---|
| 1162 | REAL alpha,beta,gama,delta,deltaz,h |
---|
| 1163 | INTEGER np,ierr |
---|
| 1164 | REAL pi,x |
---|
| 1165 | |
---|
| 1166 | !----------------------------------------------------------------------- |
---|
| 1167 | ! |
---|
| 1168 | pi=2.*ASIN(1.) |
---|
| 1169 | |
---|
[2019] | 1170 | OPEN(99,file='sigma.def',status='old',form='formatted', & |
---|
| 1171 | & iostat=ierr) |
---|
[2017] | 1172 | |
---|
| 1173 | !----------------------------------------------------------------------- |
---|
| 1174 | ! cas 1 on lit les options dans sigma.def: |
---|
| 1175 | ! ---------------------------------------- |
---|
| 1176 | |
---|
| 1177 | IF (ierr.eq.0) THEN |
---|
| 1178 | |
---|
| 1179 | print*,'WARNING!!! on lit les options dans sigma.def' |
---|
| 1180 | READ(99,*) deltaz |
---|
| 1181 | READ(99,*) h |
---|
| 1182 | READ(99,*) beta |
---|
| 1183 | READ(99,*) gama |
---|
| 1184 | READ(99,*) delta |
---|
| 1185 | READ(99,*) np |
---|
| 1186 | CLOSE(99) |
---|
| 1187 | alpha=deltaz/(llm*h) |
---|
| 1188 | ! |
---|
| 1189 | |
---|
| 1190 | DO 1 l = 1, llm |
---|
[2019] | 1191 | dsig(l) = (alpha+(1.-alpha)*exp(-beta*(llm-l)))* & |
---|
| 1192 | & ( (tanh(gama*l)/tanh(gama*llm))**np + & |
---|
| 1193 | & (1.-l/FLOAT(llm))*delta ) |
---|
[2017] | 1194 | 1 CONTINUE |
---|
| 1195 | |
---|
| 1196 | sig(1)=1. |
---|
| 1197 | DO 101 l=1,llm-1 |
---|
| 1198 | sig(l+1)=sig(l)*(1.-dsig(l))/(1.+dsig(l)) |
---|
| 1199 | 101 CONTINUE |
---|
| 1200 | sig(llm+1)=0. |
---|
| 1201 | |
---|
| 1202 | DO 2 l = 1, llm |
---|
| 1203 | dsig(l) = sig(l)-sig(l+1) |
---|
| 1204 | 2 CONTINUE |
---|
| 1205 | ! |
---|
| 1206 | |
---|
| 1207 | ELSE |
---|
| 1208 | !----------------------------------------------------------------------- |
---|
| 1209 | ! cas 2 ancienne discretisation (LMD5...): |
---|
| 1210 | ! ---------------------------------------- |
---|
| 1211 | |
---|
| 1212 | PRINT*,'WARNING!!! Ancienne discretisation verticale' |
---|
| 1213 | |
---|
| 1214 | h=7. |
---|
| 1215 | snorm = 0. |
---|
| 1216 | DO l = 1, llm |
---|
| 1217 | x = 2.*asin(1.) * (FLOAT(l)-0.5) / float(llm+1) |
---|
| 1218 | dsig(l) = 1.0 + 7.0 * SIN(x)**2 |
---|
| 1219 | snorm = snorm + dsig(l) |
---|
| 1220 | ENDDO |
---|
| 1221 | snorm = 1./snorm |
---|
| 1222 | DO l = 1, llm |
---|
| 1223 | dsig(l) = dsig(l)*snorm |
---|
| 1224 | ENDDO |
---|
| 1225 | sig(llm+1) = 0. |
---|
| 1226 | DO l = llm, 1, -1 |
---|
| 1227 | sig(l) = sig(l+1) + dsig(l) |
---|
| 1228 | ENDDO |
---|
| 1229 | |
---|
| 1230 | ENDIF |
---|
| 1231 | |
---|
| 1232 | |
---|
| 1233 | DO l=1,llm |
---|
| 1234 | nivsigs(l) = FLOAT(l) |
---|
| 1235 | ENDDO |
---|
| 1236 | |
---|
| 1237 | DO l=1,llmp1 |
---|
| 1238 | nivsig(l)= FLOAT(l) |
---|
| 1239 | ENDDO |
---|
| 1240 | |
---|
| 1241 | ! |
---|
| 1242 | ! .... Calculs de ap(l) et de bp(l) .... |
---|
| 1243 | ! ......................................... |
---|
| 1244 | ! |
---|
| 1245 | ! |
---|
| 1246 | ! ..... pa et preff sont lus sur les fichiers start par lectba ..... |
---|
| 1247 | ! |
---|
| 1248 | |
---|
| 1249 | bp(llmp1) = 0. |
---|
| 1250 | |
---|
| 1251 | DO l = 1, llm |
---|
| 1252 | !c |
---|
| 1253 | !cc ap(l) = 0. |
---|
| 1254 | !cc bp(l) = sig(l) |
---|
| 1255 | |
---|
| 1256 | bp(l) = EXP( 1. -1./( sig(l)*sig(l)) ) |
---|
| 1257 | ap(l) = pa * ( sig(l) - bp(l) ) |
---|
| 1258 | ! |
---|
| 1259 | ENDDO |
---|
| 1260 | ap(llmp1) = pa * ( sig(llmp1) - bp(llmp1) ) |
---|
| 1261 | |
---|
| 1262 | PRINT *,' BP ' |
---|
| 1263 | PRINT *, bp |
---|
| 1264 | PRINT *,' AP ' |
---|
| 1265 | PRINT *, ap |
---|
| 1266 | |
---|
[2096] | 1267 | DO l = 1, llm |
---|
| 1268 | dpres(l) = bp(l) - bp(l+1) |
---|
| 1269 | presnivs(l) = 0.5 *( ap(l)+bp(l)*preff + ap(l+1)+bp(l+1)*preff ) |
---|
| 1270 | ENDDO |
---|
[2017] | 1271 | |
---|
| 1272 | PRINT *,' PRESNIVS ' |
---|
| 1273 | PRINT *,presnivs |
---|
| 1274 | |
---|
| 1275 | RETURN |
---|
| 1276 | END |
---|
| 1277 | |
---|
[3780] | 1278 | !!====================================================================== |
---|
| 1279 | ! SUBROUTINE read_tsurf1d(knon,sst_out) |
---|
| 1280 | ! |
---|
| 1281 | !! This subroutine specifies the surface temperature to be used in 1D simulations |
---|
| 1282 | ! |
---|
| 1283 | ! USE dimphy, ONLY : klon |
---|
| 1284 | ! |
---|
| 1285 | ! INTEGER, INTENT(IN) :: knon ! nomber of points on compressed grid |
---|
| 1286 | ! REAL, DIMENSION(klon), INTENT(OUT) :: sst_out ! tsurf used to force the single-column model |
---|
| 1287 | ! |
---|
| 1288 | ! INTEGER :: i |
---|
| 1289 | !! COMMON defined in lmdz1d.F: |
---|
| 1290 | ! real ts_cur |
---|
| 1291 | ! common /sst_forcing/ts_cur |
---|
[2017] | 1292 | |
---|
[3780] | 1293 | ! DO i = 1, knon |
---|
| 1294 | ! sst_out(i) = ts_cur |
---|
| 1295 | ! ENDDO |
---|
| 1296 | ! |
---|
| 1297 | ! END SUBROUTINE read_tsurf1d |
---|
| 1298 | ! |
---|
[2017] | 1299 | !=============================================================== |
---|
| 1300 | subroutine advect_vert(llm,w,dt,q,plev) |
---|
| 1301 | !=============================================================== |
---|
| 1302 | ! Schema amont pour l'advection verticale en 1D |
---|
| 1303 | ! w est la vitesse verticale dp/dt en Pa/s |
---|
| 1304 | ! Traitement en volumes finis |
---|
| 1305 | ! d / dt ( zm q ) = delta_z ( omega q ) |
---|
| 1306 | ! d / dt ( zm ) = delta_z ( omega ) |
---|
| 1307 | ! avec zm = delta_z ( p ) |
---|
| 1308 | ! si * designe la valeur au pas de temps t+dt |
---|
| 1309 | ! zm*(l) q*(l) - zm(l) q(l) = w(l+1) q(l+1) - w(l) q(l) |
---|
| 1310 | ! zm*(l) -zm(l) = w(l+1) - w(l) |
---|
| 1311 | ! avec w=omega * dt |
---|
| 1312 | !--------------------------------------------------------------- |
---|
| 1313 | implicit none |
---|
| 1314 | ! arguments |
---|
| 1315 | integer llm |
---|
| 1316 | real w(llm+1),q(llm),plev(llm+1),dt |
---|
| 1317 | |
---|
| 1318 | ! local |
---|
| 1319 | integer l |
---|
| 1320 | real zwq(llm+1),zm(llm+1),zw(llm+1) |
---|
| 1321 | real qold |
---|
| 1322 | |
---|
| 1323 | !--------------------------------------------------------------- |
---|
| 1324 | |
---|
| 1325 | do l=1,llm |
---|
| 1326 | zw(l)=dt*w(l) |
---|
| 1327 | zm(l)=plev(l)-plev(l+1) |
---|
| 1328 | zwq(l)=q(l)*zw(l) |
---|
| 1329 | enddo |
---|
| 1330 | zwq(llm+1)=0. |
---|
| 1331 | zw(llm+1)=0. |
---|
| 1332 | |
---|
| 1333 | do l=1,llm |
---|
| 1334 | qold=q(l) |
---|
| 1335 | q(l)=(q(l)*zm(l)+zwq(l+1)-zwq(l))/(zm(l)+zw(l+1)-zw(l)) |
---|
| 1336 | print*,'ADV Q ',zm(l),zw(l),zwq(l),qold,q(l) |
---|
| 1337 | enddo |
---|
| 1338 | |
---|
| 1339 | |
---|
| 1340 | return |
---|
| 1341 | end |
---|
| 1342 | |
---|
| 1343 | !=============================================================== |
---|
| 1344 | |
---|
| 1345 | |
---|
[2019] | 1346 | SUBROUTINE advect_va(llm,omega,d_t_va,d_q_va,d_u_va,d_v_va, & |
---|
| 1347 | & q,temp,u,v,play) |
---|
[2017] | 1348 | !itlmd |
---|
| 1349 | !---------------------------------------------------------------------- |
---|
| 1350 | ! Calcul de l'advection verticale (ascendance et subsidence) de |
---|
[2716] | 1351 | ! temperature et d'humidite. Hypothese : ce qui rentre de l'exterieur |
---|
| 1352 | ! a les memes caracteristiques que l'air de la colonne 1D (WTG) ou |
---|
[2017] | 1353 | ! sans WTG rajouter une advection horizontale |
---|
| 1354 | !---------------------------------------------------------------------- |
---|
| 1355 | implicit none |
---|
[4603] | 1356 | #include "YOMCST.h" |
---|
[2017] | 1357 | ! argument |
---|
| 1358 | integer llm |
---|
| 1359 | real omega(llm+1),d_t_va(llm), d_q_va(llm,3) |
---|
| 1360 | real d_u_va(llm), d_v_va(llm) |
---|
| 1361 | real q(llm,3),temp(llm) |
---|
| 1362 | real u(llm),v(llm) |
---|
| 1363 | real play(llm) |
---|
| 1364 | ! interne |
---|
| 1365 | integer l |
---|
| 1366 | real alpha,omgdown,omgup |
---|
| 1367 | |
---|
| 1368 | do l= 1,llm |
---|
| 1369 | if(l.eq.1) then |
---|
| 1370 | !si omgup pour la couche 1, alors tendance nulle |
---|
| 1371 | omgdown=max(omega(2),0.0) |
---|
[2019] | 1372 | alpha = rkappa*temp(l)*(1.+q(l,1)*rv/rd)/(play(l)*(1.+q(l,1))) |
---|
| 1373 | d_t_va(l)= alpha*(omgdown)-omgdown*(temp(l)-temp(l+1)) & |
---|
| 1374 | & /(play(l)-play(l+1)) |
---|
[2017] | 1375 | |
---|
[2019] | 1376 | d_q_va(l,:)= -omgdown*(q(l,:)-q(l+1,:))/(play(l)-play(l+1)) |
---|
[2017] | 1377 | |
---|
[2019] | 1378 | d_u_va(l)= -omgdown*(u(l)-u(l+1))/(play(l)-play(l+1)) |
---|
| 1379 | d_v_va(l)= -omgdown*(v(l)-v(l+1))/(play(l)-play(l+1)) |
---|
[2017] | 1380 | |
---|
| 1381 | |
---|
| 1382 | elseif(l.eq.llm) then |
---|
| 1383 | omgup=min(omega(l),0.0) |
---|
[2019] | 1384 | alpha = rkappa*temp(l)*(1.+q(l,1)*rv/rd)/(play(l)*(1.+q(l,1))) |
---|
| 1385 | d_t_va(l)= alpha*(omgup)- & |
---|
| 1386 | |
---|
[2017] | 1387 | !bug? & omgup*(temp(l-1)-temp(l))/(play(l-1)-plev(l)) |
---|
| 1388 | & omgup*(temp(l-1)-temp(l))/(play(l-1)-play(l)) |
---|
| 1389 | d_q_va(l,:)= -omgup*(q(l-1,:)-q(l,:))/(play(l-1)-play(l)) |
---|
| 1390 | d_u_va(l)= -omgup*(u(l-1)-u(l))/(play(l-1)-play(l)) |
---|
| 1391 | d_v_va(l)= -omgup*(v(l-1)-v(l))/(play(l-1)-play(l)) |
---|
| 1392 | |
---|
| 1393 | else |
---|
| 1394 | omgup=min(omega(l),0.0) |
---|
| 1395 | omgdown=max(omega(l+1),0.0) |
---|
[2019] | 1396 | alpha = rkappa*temp(l)*(1.+q(l,1)*rv/rd)/(play(l)*(1.+q(l,1))) |
---|
| 1397 | d_t_va(l)= alpha*(omgup+omgdown)-omgdown*(temp(l)-temp(l+1)) & |
---|
| 1398 | & /(play(l)-play(l+1))- & |
---|
[2017] | 1399 | !bug? & omgup*(temp(l-1)-temp(l))/(play(l-1)-plev(l)) |
---|
| 1400 | & omgup*(temp(l-1)-temp(l))/(play(l-1)-play(l)) |
---|
| 1401 | ! print*, ' ??? ' |
---|
| 1402 | |
---|
[2019] | 1403 | d_q_va(l,:)= -omgdown*(q(l,:)-q(l+1,:)) & |
---|
| 1404 | & /(play(l)-play(l+1))- & |
---|
| 1405 | & omgup*(q(l-1,:)-q(l,:))/(play(l-1)-play(l)) |
---|
| 1406 | d_u_va(l)= -omgdown*(u(l)-u(l+1)) & |
---|
| 1407 | & /(play(l)-play(l+1))- & |
---|
| 1408 | & omgup*(u(l-1)-u(l))/(play(l-1)-play(l)) |
---|
| 1409 | d_v_va(l)= -omgdown*(v(l)-v(l+1)) & |
---|
| 1410 | & /(play(l)-play(l+1))- & |
---|
[2017] | 1411 | & omgup*(v(l-1)-v(l))/(play(l-1)-play(l)) |
---|
| 1412 | |
---|
| 1413 | endif |
---|
| 1414 | |
---|
| 1415 | enddo |
---|
| 1416 | !fin itlmd |
---|
| 1417 | return |
---|
| 1418 | end |
---|
| 1419 | ! SUBROUTINE lstendH(llm,omega,d_t_va,d_q_va,d_u_va,d_v_va, |
---|
[2019] | 1420 | SUBROUTINE lstendH(llm,nqtot,omega,d_t_va,d_q_va, & |
---|
| 1421 | & q,temp,u,v,play) |
---|
[2017] | 1422 | !itlmd |
---|
| 1423 | !---------------------------------------------------------------------- |
---|
| 1424 | ! Calcul de l'advection verticale (ascendance et subsidence) de |
---|
[2716] | 1425 | ! temperature et d'humidite. Hypothese : ce qui rentre de l'exterieur |
---|
| 1426 | ! a les memes caracteristiques que l'air de la colonne 1D (WTG) ou |
---|
[2017] | 1427 | ! sans WTG rajouter une advection horizontale |
---|
| 1428 | !---------------------------------------------------------------------- |
---|
| 1429 | implicit none |
---|
[4603] | 1430 | #include "YOMCST.h" |
---|
[2017] | 1431 | ! argument |
---|
| 1432 | integer llm,nqtot |
---|
| 1433 | real omega(llm+1),d_t_va(llm), d_q_va(llm,nqtot) |
---|
| 1434 | ! real d_u_va(llm), d_v_va(llm) |
---|
| 1435 | real q(llm,nqtot),temp(llm) |
---|
| 1436 | real u(llm),v(llm) |
---|
| 1437 | real play(llm) |
---|
| 1438 | real cor(llm) |
---|
| 1439 | ! real dph(llm),dudp(llm),dvdp(llm),dqdp(llm),dtdp(llm) |
---|
| 1440 | real dph(llm),dqdp(llm),dtdp(llm) |
---|
| 1441 | ! interne |
---|
| 1442 | integer k |
---|
| 1443 | real omdn,omup |
---|
| 1444 | |
---|
| 1445 | ! dudp=0. |
---|
| 1446 | ! dvdp=0. |
---|
| 1447 | dqdp=0. |
---|
| 1448 | dtdp=0. |
---|
| 1449 | ! d_u_va=0. |
---|
| 1450 | ! d_v_va=0. |
---|
| 1451 | |
---|
| 1452 | cor(:) = rkappa*temp*(1.+q(:,1)*rv/rd)/(play*(1.+q(:,1))) |
---|
| 1453 | |
---|
[2933] | 1454 | |
---|
[2017] | 1455 | do k=2,llm-1 |
---|
| 1456 | |
---|
| 1457 | dph (k-1) = (play(k )- play(k-1 )) |
---|
| 1458 | ! dudp (k-1) = (u (k )- u (k-1 ))/dph(k-1) |
---|
| 1459 | ! dvdp (k-1) = (v (k )- v (k-1 ))/dph(k-1) |
---|
| 1460 | dqdp (k-1) = (q (k,1)- q (k-1,1))/dph(k-1) |
---|
| 1461 | dtdp (k-1) = (temp(k )- temp(k-1 ))/dph(k-1) |
---|
| 1462 | |
---|
| 1463 | enddo |
---|
| 1464 | |
---|
| 1465 | ! dudp ( llm ) = dudp ( llm-1 ) |
---|
| 1466 | ! dvdp ( llm ) = dvdp ( llm-1 ) |
---|
| 1467 | dqdp ( llm ) = dqdp ( llm-1 ) |
---|
| 1468 | dtdp ( llm ) = dtdp ( llm-1 ) |
---|
| 1469 | |
---|
| 1470 | do k=2,llm-1 |
---|
| 1471 | omdn=max(0.0,omega(k+1)) |
---|
| 1472 | omup=min(0.0,omega( k )) |
---|
| 1473 | |
---|
| 1474 | ! d_u_va(k) = -omdn*dudp(k)-omup*dudp(k-1) |
---|
| 1475 | ! d_v_va(k) = -omdn*dvdp(k)-omup*dvdp(k-1) |
---|
| 1476 | d_q_va(k,1)= -omdn*dqdp(k)-omup*dqdp(k-1) |
---|
[2019] | 1477 | d_t_va(k) = -omdn*dtdp(k)-omup*dtdp(k-1)+(omup+omdn)*cor(k) |
---|
[2017] | 1478 | enddo |
---|
| 1479 | |
---|
| 1480 | omdn=max(0.0,omega( 2 )) |
---|
| 1481 | omup=min(0.0,omega(llm)) |
---|
| 1482 | ! d_u_va( 1 ) = -omdn*dudp( 1 ) |
---|
| 1483 | ! d_u_va(llm) = -omup*dudp(llm) |
---|
| 1484 | ! d_v_va( 1 ) = -omdn*dvdp( 1 ) |
---|
| 1485 | ! d_v_va(llm) = -omup*dvdp(llm) |
---|
| 1486 | d_q_va( 1 ,1) = -omdn*dqdp( 1 ) |
---|
| 1487 | d_q_va(llm,1) = -omup*dqdp(llm) |
---|
| 1488 | d_t_va( 1 ) = -omdn*dtdp( 1 )+omdn*cor( 1 ) |
---|
| 1489 | d_t_va(llm) = -omup*dtdp(llm)!+omup*cor(llm) |
---|
| 1490 | |
---|
| 1491 | ! if(abs(rlat(1))>10.) then |
---|
| 1492 | ! Calculate the tendency due agestrophic motions |
---|
| 1493 | ! du_age = fcoriolis*(v-vg) |
---|
| 1494 | ! dv_age = fcoriolis*(ug-u) |
---|
| 1495 | ! endif |
---|
| 1496 | |
---|
| 1497 | ! call writefield_phy('d_t_va',d_t_va,llm) |
---|
| 1498 | |
---|
| 1499 | return |
---|
| 1500 | end |
---|
| 1501 | |
---|
| 1502 | !====================================================================== |
---|
| 1503 | |
---|
[2181] | 1504 | ! Subroutines for nudging |
---|
| 1505 | |
---|
| 1506 | Subroutine Nudge_RHT_init (paprs,pplay,t,q,t_targ,rh_targ) |
---|
| 1507 | ! ======================================================== |
---|
| 1508 | USE dimphy |
---|
| 1509 | |
---|
| 1510 | implicit none |
---|
| 1511 | |
---|
| 1512 | ! ======================================================== |
---|
| 1513 | REAL paprs(klon,klevp1) |
---|
| 1514 | REAL pplay(klon,klev) |
---|
| 1515 | ! |
---|
| 1516 | ! Variables d'etat |
---|
| 1517 | REAL t(klon,klev) |
---|
| 1518 | REAL q(klon,klev) |
---|
| 1519 | ! |
---|
| 1520 | ! Profiles cible |
---|
| 1521 | REAL t_targ(klon,klev) |
---|
| 1522 | REAL rh_targ(klon,klev) |
---|
| 1523 | ! |
---|
| 1524 | INTEGER k,i |
---|
| 1525 | REAL zx_qs |
---|
| 1526 | |
---|
| 1527 | ! Declaration des constantes et des fonctions thermodynamiques |
---|
| 1528 | ! |
---|
| 1529 | include "YOMCST.h" |
---|
| 1530 | include "YOETHF.h" |
---|
| 1531 | ! |
---|
| 1532 | ! ---------------------------------------- |
---|
| 1533 | ! Statement functions |
---|
| 1534 | include "FCTTRE.h" |
---|
| 1535 | ! ---------------------------------------- |
---|
| 1536 | ! |
---|
| 1537 | DO k = 1,klev |
---|
| 1538 | DO i = 1,klon |
---|
| 1539 | t_targ(i,k) = t(i,k) |
---|
| 1540 | IF (t(i,k).LT.RTT) THEN |
---|
| 1541 | zx_qs = qsats(t(i,k))/(pplay(i,k)) |
---|
| 1542 | ELSE |
---|
| 1543 | zx_qs = qsatl(t(i,k))/(pplay(i,k)) |
---|
| 1544 | ENDIF |
---|
| 1545 | rh_targ(i,k) = q(i,k)/zx_qs |
---|
| 1546 | ENDDO |
---|
| 1547 | ENDDO |
---|
| 1548 | print *, 't_targ',t_targ |
---|
| 1549 | print *, 'rh_targ',rh_targ |
---|
| 1550 | ! |
---|
| 1551 | ! |
---|
| 1552 | RETURN |
---|
| 1553 | END |
---|
| 1554 | |
---|
| 1555 | Subroutine Nudge_UV_init (paprs,pplay,u,v,u_targ,v_targ) |
---|
| 1556 | ! ======================================================== |
---|
| 1557 | USE dimphy |
---|
| 1558 | |
---|
| 1559 | implicit none |
---|
| 1560 | |
---|
| 1561 | ! ======================================================== |
---|
| 1562 | REAL paprs(klon,klevp1) |
---|
| 1563 | REAL pplay(klon,klev) |
---|
| 1564 | ! |
---|
| 1565 | ! Variables d'etat |
---|
| 1566 | REAL u(klon,klev) |
---|
| 1567 | REAL v(klon,klev) |
---|
| 1568 | ! |
---|
| 1569 | ! Profiles cible |
---|
| 1570 | REAL u_targ(klon,klev) |
---|
| 1571 | REAL v_targ(klon,klev) |
---|
| 1572 | ! |
---|
| 1573 | INTEGER k,i |
---|
| 1574 | ! |
---|
| 1575 | DO k = 1,klev |
---|
| 1576 | DO i = 1,klon |
---|
| 1577 | u_targ(i,k) = u(i,k) |
---|
| 1578 | v_targ(i,k) = v(i,k) |
---|
| 1579 | ENDDO |
---|
| 1580 | ENDDO |
---|
| 1581 | print *, 'u_targ',u_targ |
---|
| 1582 | print *, 'v_targ',v_targ |
---|
| 1583 | ! |
---|
| 1584 | ! |
---|
| 1585 | RETURN |
---|
| 1586 | END |
---|
| 1587 | |
---|
| 1588 | Subroutine Nudge_RHT (dtime,paprs,pplay,t_targ,rh_targ,t,q, & |
---|
| 1589 | & d_t,d_q) |
---|
| 1590 | ! ======================================================== |
---|
| 1591 | USE dimphy |
---|
| 1592 | |
---|
| 1593 | implicit none |
---|
| 1594 | |
---|
| 1595 | ! ======================================================== |
---|
| 1596 | REAL dtime |
---|
| 1597 | REAL paprs(klon,klevp1) |
---|
| 1598 | REAL pplay(klon,klev) |
---|
| 1599 | ! |
---|
| 1600 | ! Variables d'etat |
---|
| 1601 | REAL t(klon,klev) |
---|
| 1602 | REAL q(klon,klev) |
---|
| 1603 | ! |
---|
| 1604 | ! Tendances |
---|
| 1605 | REAL d_t(klon,klev) |
---|
| 1606 | REAL d_q(klon,klev) |
---|
| 1607 | ! |
---|
| 1608 | ! Profiles cible |
---|
| 1609 | REAL t_targ(klon,klev) |
---|
| 1610 | REAL rh_targ(klon,klev) |
---|
| 1611 | ! |
---|
| 1612 | ! Temps de relaxation |
---|
| 1613 | REAL tau |
---|
| 1614 | !c DATA tau /3600./ |
---|
| 1615 | !! DATA tau /5400./ |
---|
| 1616 | DATA tau /1800./ |
---|
| 1617 | ! |
---|
| 1618 | INTEGER k,i |
---|
[2455] | 1619 | REAL zx_qs, rh, tnew, d_rh, rhnew |
---|
[2181] | 1620 | |
---|
| 1621 | ! Declaration des constantes et des fonctions thermodynamiques |
---|
| 1622 | ! |
---|
| 1623 | include "YOMCST.h" |
---|
| 1624 | include "YOETHF.h" |
---|
| 1625 | ! |
---|
| 1626 | ! ---------------------------------------- |
---|
| 1627 | ! Statement functions |
---|
| 1628 | include "FCTTRE.h" |
---|
| 1629 | ! ---------------------------------------- |
---|
| 1630 | ! |
---|
| 1631 | print *,'dtime, tau ',dtime,tau |
---|
| 1632 | print *, 't_targ',t_targ |
---|
| 1633 | print *, 'rh_targ',rh_targ |
---|
| 1634 | print *,'temp ',t |
---|
| 1635 | print *,'hum ',q |
---|
[2455] | 1636 | ! |
---|
[2181] | 1637 | DO k = 1,klev |
---|
| 1638 | DO i = 1,klon |
---|
[2455] | 1639 | IF (paprs(i,1)-pplay(i,k) .GT. 10000.) THEN |
---|
[2181] | 1640 | IF (t(i,k).LT.RTT) THEN |
---|
| 1641 | zx_qs = qsats(t(i,k))/(pplay(i,k)) |
---|
| 1642 | ELSE |
---|
| 1643 | zx_qs = qsatl(t(i,k))/(pplay(i,k)) |
---|
| 1644 | ENDIF |
---|
| 1645 | rh = q(i,k)/zx_qs |
---|
| 1646 | ! |
---|
| 1647 | d_t(i,k) = d_t(i,k) + 1./tau*(t_targ(i,k)-t(i,k)) |
---|
| 1648 | d_rh = 1./tau*(rh_targ(i,k)-rh) |
---|
| 1649 | ! |
---|
[2455] | 1650 | tnew = t(i,k)+d_t(i,k)*dtime |
---|
| 1651 | !jyg< |
---|
| 1652 | ! Formule pour q : |
---|
| 1653 | ! d_q = (1/tau) [rh_targ*qsat(T_new) - q] |
---|
| 1654 | ! |
---|
| 1655 | ! Cette formule remplace d_q = (1/tau) [rh_targ - rh] qsat(T_new) |
---|
[2716] | 1656 | ! qui n'etait pas correcte. |
---|
[2455] | 1657 | ! |
---|
[2181] | 1658 | IF (tnew.LT.RTT) THEN |
---|
| 1659 | zx_qs = qsats(tnew)/(pplay(i,k)) |
---|
| 1660 | ELSE |
---|
| 1661 | zx_qs = qsatl(tnew)/(pplay(i,k)) |
---|
| 1662 | ENDIF |
---|
[2455] | 1663 | !! d_q(i,k) = d_q(i,k) + d_rh*zx_qs |
---|
| 1664 | d_q(i,k) = d_q(i,k) + (1./tau)*(rh_targ(i,k)*zx_qs - q(i,k)) |
---|
| 1665 | rhnew = (q(i,k)+d_q(i,k)*dtime)/zx_qs |
---|
[2181] | 1666 | ! |
---|
[2455] | 1667 | print *,' k,d_t,rh,d_rh,rhnew,d_q ', & |
---|
| 1668 | k,d_t(i,k),rh,d_rh,rhnew,d_q(i,k) |
---|
| 1669 | ENDIF |
---|
[2181] | 1670 | ! |
---|
| 1671 | ENDDO |
---|
| 1672 | ENDDO |
---|
| 1673 | ! |
---|
| 1674 | RETURN |
---|
| 1675 | END |
---|
| 1676 | |
---|
| 1677 | Subroutine Nudge_UV (dtime,paprs,pplay,u_targ,v_targ,u,v, & |
---|
| 1678 | & d_u,d_v) |
---|
| 1679 | ! ======================================================== |
---|
| 1680 | USE dimphy |
---|
| 1681 | |
---|
| 1682 | implicit none |
---|
| 1683 | |
---|
| 1684 | ! ======================================================== |
---|
| 1685 | REAL dtime |
---|
| 1686 | REAL paprs(klon,klevp1) |
---|
| 1687 | REAL pplay(klon,klev) |
---|
| 1688 | ! |
---|
| 1689 | ! Variables d'etat |
---|
| 1690 | REAL u(klon,klev) |
---|
| 1691 | REAL v(klon,klev) |
---|
| 1692 | ! |
---|
| 1693 | ! Tendances |
---|
| 1694 | REAL d_u(klon,klev) |
---|
| 1695 | REAL d_v(klon,klev) |
---|
| 1696 | ! |
---|
| 1697 | ! Profiles cible |
---|
| 1698 | REAL u_targ(klon,klev) |
---|
| 1699 | REAL v_targ(klon,klev) |
---|
| 1700 | ! |
---|
| 1701 | ! Temps de relaxation |
---|
| 1702 | REAL tau |
---|
| 1703 | !c DATA tau /3600./ |
---|
[2933] | 1704 | ! DATA tau /5400./ |
---|
| 1705 | DATA tau /43200./ |
---|
[2181] | 1706 | ! |
---|
| 1707 | INTEGER k,i |
---|
| 1708 | |
---|
| 1709 | ! |
---|
[4292] | 1710 | !print *,'dtime, tau ',dtime,tau |
---|
| 1711 | !print *, 'u_targ',u_targ |
---|
| 1712 | !print *, 'v_targ',v_targ |
---|
| 1713 | !print *,'zonal velocity ',u |
---|
| 1714 | !print *,'meridional velocity ',v |
---|
[2181] | 1715 | DO k = 1,klev |
---|
| 1716 | DO i = 1,klon |
---|
[2933] | 1717 | !CR: nudging everywhere |
---|
| 1718 | ! IF (paprs(i,1)-pplay(i,k) .GT. 10000.) THEN |
---|
[2181] | 1719 | ! |
---|
| 1720 | d_u(i,k) = d_u(i,k) + 1./tau*(u_targ(i,k)-u(i,k)) |
---|
| 1721 | d_v(i,k) = d_v(i,k) + 1./tau*(v_targ(i,k)-v(i,k)) |
---|
| 1722 | ! |
---|
[4292] | 1723 | ! print *,' k,u,d_u,v,d_v ', & |
---|
| 1724 | ! k,u(i,k),d_u(i,k),v(i,k),d_v(i,k) |
---|
[2933] | 1725 | ! ENDIF |
---|
[2181] | 1726 | ! |
---|
| 1727 | ENDDO |
---|
| 1728 | ENDDO |
---|
| 1729 | ! |
---|
| 1730 | RETURN |
---|
| 1731 | END |
---|
| 1732 | |
---|
[2716] | 1733 | !===================================================================== |
---|
| 1734 | SUBROUTINE interp2_case_vertical(play,nlev_cas,plev_prof_cas & |
---|
| 1735 | & ,t_prof_cas,th_prof_cas,thv_prof_cas,thl_prof_cas & |
---|
| 1736 | & ,qv_prof_cas,ql_prof_cas,qi_prof_cas,u_prof_cas,v_prof_cas & |
---|
| 1737 | & ,ug_prof_cas,vg_prof_cas,vitw_prof_cas,omega_prof_cas & |
---|
| 1738 | & ,du_prof_cas,hu_prof_cas,vu_prof_cas,dv_prof_cas,hv_prof_cas,vv_prof_cas & |
---|
| 1739 | & ,dt_prof_cas,ht_prof_cas,vt_prof_cas,dtrad_prof_cas,dq_prof_cas,hq_prof_cas,vq_prof_cas & |
---|
| 1740 | & ,dth_prof_cas,hth_prof_cas,vth_prof_cas & |
---|
| 1741 | ! |
---|
| 1742 | & ,t_mod_cas,theta_mod_cas,thv_mod_cas,thl_mod_cas & |
---|
| 1743 | & ,qv_mod_cas,ql_mod_cas,qi_mod_cas,u_mod_cas,v_mod_cas & |
---|
| 1744 | & ,ug_mod_cas,vg_mod_cas,w_mod_cas,omega_mod_cas & |
---|
| 1745 | & ,du_mod_cas,hu_mod_cas,vu_mod_cas,dv_mod_cas,hv_mod_cas,vv_mod_cas & |
---|
| 1746 | & ,dt_mod_cas,ht_mod_cas,vt_mod_cas,dtrad_mod_cas,dq_mod_cas,hq_mod_cas,vq_mod_cas & |
---|
| 1747 | & ,dth_mod_cas,hth_mod_cas,vth_mod_cas,mxcalc) |
---|
| 1748 | |
---|
| 1749 | implicit none |
---|
| 1750 | |
---|
[4603] | 1751 | #include "YOMCST.h" |
---|
| 1752 | #include "dimensions.h" |
---|
[2683] | 1753 | |
---|
[2716] | 1754 | !------------------------------------------------------------------------- |
---|
| 1755 | ! Vertical interpolation of generic case forcing data onto mod_casel levels |
---|
| 1756 | !------------------------------------------------------------------------- |
---|
| 1757 | |
---|
| 1758 | integer nlevmax |
---|
| 1759 | parameter (nlevmax=41) |
---|
| 1760 | integer nlev_cas,mxcalc |
---|
| 1761 | ! real play(llm), plev_prof(nlevmax) |
---|
| 1762 | ! real t_prof(nlevmax),q_prof(nlevmax) |
---|
| 1763 | ! real u_prof(nlevmax),v_prof(nlevmax), w_prof(nlevmax) |
---|
| 1764 | ! real ht_prof(nlevmax),vt_prof(nlevmax) |
---|
| 1765 | ! real hq_prof(nlevmax),vq_prof(nlevmax) |
---|
| 1766 | |
---|
| 1767 | real play(llm), plev_prof_cas(nlev_cas) |
---|
| 1768 | real t_prof_cas(nlev_cas),th_prof_cas(nlev_cas),thv_prof_cas(nlev_cas),thl_prof_cas(nlev_cas) |
---|
| 1769 | real qv_prof_cas(nlev_cas),ql_prof_cas(nlev_cas),qi_prof_cas(nlev_cas) |
---|
| 1770 | real u_prof_cas(nlev_cas),v_prof_cas(nlev_cas) |
---|
| 1771 | real ug_prof_cas(nlev_cas),vg_prof_cas(nlev_cas), vitw_prof_cas(nlev_cas),omega_prof_cas(nlev_cas) |
---|
| 1772 | real du_prof_cas(nlev_cas),hu_prof_cas(nlev_cas),vu_prof_cas(nlev_cas) |
---|
| 1773 | real dv_prof_cas(nlev_cas),hv_prof_cas(nlev_cas),vv_prof_cas(nlev_cas) |
---|
| 1774 | real dt_prof_cas(nlev_cas),ht_prof_cas(nlev_cas),vt_prof_cas(nlev_cas),dtrad_prof_cas(nlev_cas) |
---|
| 1775 | real dth_prof_cas(nlev_cas),hth_prof_cas(nlev_cas),vth_prof_cas(nlev_cas) |
---|
| 1776 | real dq_prof_cas(nlev_cas),hq_prof_cas(nlev_cas),vq_prof_cas(nlev_cas) |
---|
| 1777 | |
---|
| 1778 | real t_mod_cas(llm),theta_mod_cas(llm),thv_mod_cas(llm),thl_mod_cas(llm) |
---|
| 1779 | real qv_mod_cas(llm),ql_mod_cas(llm),qi_mod_cas(llm) |
---|
| 1780 | real u_mod_cas(llm),v_mod_cas(llm) |
---|
| 1781 | real ug_mod_cas(llm),vg_mod_cas(llm), w_mod_cas(llm),omega_mod_cas(llm) |
---|
| 1782 | real du_mod_cas(llm),hu_mod_cas(llm),vu_mod_cas(llm) |
---|
| 1783 | real dv_mod_cas(llm),hv_mod_cas(llm),vv_mod_cas(llm) |
---|
| 1784 | real dt_mod_cas(llm),ht_mod_cas(llm),vt_mod_cas(llm),dtrad_mod_cas(llm) |
---|
| 1785 | real dth_mod_cas(llm),hth_mod_cas(llm),vth_mod_cas(llm) |
---|
| 1786 | real dq_mod_cas(llm),hq_mod_cas(llm),vq_mod_cas(llm) |
---|
| 1787 | |
---|
| 1788 | integer l,k,k1,k2 |
---|
| 1789 | real frac,frac1,frac2,fact |
---|
| 1790 | |
---|
[3541] | 1791 | ! do l = 1, llm |
---|
| 1792 | ! print *,'debut interp2, play=',l,play(l) |
---|
| 1793 | ! enddo |
---|
[2716] | 1794 | ! do l = 1, nlev_cas |
---|
| 1795 | ! print *,'debut interp2, plev_prof_cas=',l,play(l),plev_prof_cas(l) |
---|
| 1796 | ! enddo |
---|
| 1797 | |
---|
| 1798 | do l = 1, llm |
---|
| 1799 | |
---|
| 1800 | if (play(l).ge.plev_prof_cas(nlev_cas)) then |
---|
| 1801 | |
---|
| 1802 | mxcalc=l |
---|
[3541] | 1803 | ! print *,'debut interp2, mxcalc=',mxcalc |
---|
[2716] | 1804 | k1=0 |
---|
| 1805 | k2=0 |
---|
| 1806 | |
---|
| 1807 | if (play(l).le.plev_prof_cas(1)) then |
---|
| 1808 | |
---|
| 1809 | do k = 1, nlev_cas-1 |
---|
| 1810 | if (play(l).le.plev_prof_cas(k).and. play(l).gt.plev_prof_cas(k+1)) then |
---|
| 1811 | k1=k |
---|
| 1812 | k2=k+1 |
---|
| 1813 | endif |
---|
| 1814 | enddo |
---|
| 1815 | |
---|
| 1816 | if (k1.eq.0 .or. k2.eq.0) then |
---|
| 1817 | write(*,*) 'PB! k1, k2 = ',k1,k2 |
---|
| 1818 | write(*,*) 'l,play(l) = ',l,play(l)/100 |
---|
| 1819 | do k = 1, nlev_cas-1 |
---|
| 1820 | write(*,*) 'k,plev_prof_cas(k) = ',k,plev_prof_cas(k)/100 |
---|
| 1821 | enddo |
---|
| 1822 | endif |
---|
| 1823 | |
---|
| 1824 | frac = (plev_prof_cas(k2)-play(l))/(plev_prof_cas(k2)-plev_prof_cas(k1)) |
---|
| 1825 | t_mod_cas(l)= t_prof_cas(k2) - frac*(t_prof_cas(k2)-t_prof_cas(k1)) |
---|
| 1826 | theta_mod_cas(l)= th_prof_cas(k2) - frac*(th_prof_cas(k2)-th_prof_cas(k1)) |
---|
[3008] | 1827 | if(theta_mod_cas(l).NE.0) t_mod_cas(l)= theta_mod_cas(l)*(play(l)/100000.)**(RD/RCPD) |
---|
[2716] | 1828 | thv_mod_cas(l)= thv_prof_cas(k2) - frac*(thv_prof_cas(k2)-thv_prof_cas(k1)) |
---|
| 1829 | thl_mod_cas(l)= thl_prof_cas(k2) - frac*(thl_prof_cas(k2)-thl_prof_cas(k1)) |
---|
| 1830 | qv_mod_cas(l)= qv_prof_cas(k2) - frac*(qv_prof_cas(k2)-qv_prof_cas(k1)) |
---|
| 1831 | ql_mod_cas(l)= ql_prof_cas(k2) - frac*(ql_prof_cas(k2)-ql_prof_cas(k1)) |
---|
| 1832 | qi_mod_cas(l)= qi_prof_cas(k2) - frac*(qi_prof_cas(k2)-qi_prof_cas(k1)) |
---|
| 1833 | u_mod_cas(l)= u_prof_cas(k2) - frac*(u_prof_cas(k2)-u_prof_cas(k1)) |
---|
| 1834 | v_mod_cas(l)= v_prof_cas(k2) - frac*(v_prof_cas(k2)-v_prof_cas(k1)) |
---|
| 1835 | ug_mod_cas(l)= ug_prof_cas(k2) - frac*(ug_prof_cas(k2)-ug_prof_cas(k1)) |
---|
| 1836 | vg_mod_cas(l)= vg_prof_cas(k2) - frac*(vg_prof_cas(k2)-vg_prof_cas(k1)) |
---|
| 1837 | w_mod_cas(l)= vitw_prof_cas(k2) - frac*(vitw_prof_cas(k2)-vitw_prof_cas(k1)) |
---|
| 1838 | omega_mod_cas(l)= omega_prof_cas(k2) - frac*(omega_prof_cas(k2)-omega_prof_cas(k1)) |
---|
| 1839 | du_mod_cas(l)= du_prof_cas(k2) - frac*(du_prof_cas(k2)-du_prof_cas(k1)) |
---|
| 1840 | hu_mod_cas(l)= hu_prof_cas(k2) - frac*(hu_prof_cas(k2)-hu_prof_cas(k1)) |
---|
| 1841 | vu_mod_cas(l)= vu_prof_cas(k2) - frac*(vu_prof_cas(k2)-vu_prof_cas(k1)) |
---|
| 1842 | dv_mod_cas(l)= dv_prof_cas(k2) - frac*(dv_prof_cas(k2)-dv_prof_cas(k1)) |
---|
| 1843 | hv_mod_cas(l)= hv_prof_cas(k2) - frac*(hv_prof_cas(k2)-hv_prof_cas(k1)) |
---|
| 1844 | vv_mod_cas(l)= vv_prof_cas(k2) - frac*(vv_prof_cas(k2)-vv_prof_cas(k1)) |
---|
| 1845 | dt_mod_cas(l)= dt_prof_cas(k2) - frac*(dt_prof_cas(k2)-dt_prof_cas(k1)) |
---|
| 1846 | ht_mod_cas(l)= ht_prof_cas(k2) - frac*(ht_prof_cas(k2)-ht_prof_cas(k1)) |
---|
| 1847 | vt_mod_cas(l)= vt_prof_cas(k2) - frac*(vt_prof_cas(k2)-vt_prof_cas(k1)) |
---|
| 1848 | dth_mod_cas(l)= dth_prof_cas(k2) - frac*(dth_prof_cas(k2)-dth_prof_cas(k1)) |
---|
| 1849 | hth_mod_cas(l)= hth_prof_cas(k2) - frac*(hth_prof_cas(k2)-hth_prof_cas(k1)) |
---|
| 1850 | vth_mod_cas(l)= vth_prof_cas(k2) - frac*(vth_prof_cas(k2)-vth_prof_cas(k1)) |
---|
| 1851 | dq_mod_cas(l)= dq_prof_cas(k2) - frac*(dq_prof_cas(k2)-dq_prof_cas(k1)) |
---|
| 1852 | hq_mod_cas(l)= hq_prof_cas(k2) - frac*(hq_prof_cas(k2)-hq_prof_cas(k1)) |
---|
| 1853 | vq_mod_cas(l)= vq_prof_cas(k2) - frac*(vq_prof_cas(k2)-vq_prof_cas(k1)) |
---|
[2933] | 1854 | dtrad_mod_cas(l)= dtrad_prof_cas(k2) - frac*(dtrad_prof_cas(k2)-dtrad_prof_cas(k1)) |
---|
[2716] | 1855 | |
---|
| 1856 | else !play>plev_prof_cas(1) |
---|
| 1857 | |
---|
| 1858 | k1=1 |
---|
| 1859 | k2=2 |
---|
| 1860 | print *,'interp2_vert, k1,k2=',plev_prof_cas(k1),plev_prof_cas(k2) |
---|
| 1861 | frac1 = (play(l)-plev_prof_cas(k2))/(plev_prof_cas(k1)-plev_prof_cas(k2)) |
---|
| 1862 | frac2 = (play(l)-plev_prof_cas(k1))/(plev_prof_cas(k1)-plev_prof_cas(k2)) |
---|
| 1863 | t_mod_cas(l)= frac1*t_prof_cas(k1) - frac2*t_prof_cas(k2) |
---|
| 1864 | theta_mod_cas(l)= frac1*th_prof_cas(k1) - frac2*th_prof_cas(k2) |
---|
[3008] | 1865 | if(theta_mod_cas(l).NE.0) t_mod_cas(l)= theta_mod_cas(l)*(play(l)/100000.)**(RD/RCPD) |
---|
[2716] | 1866 | thv_mod_cas(l)= frac1*thv_prof_cas(k1) - frac2*thv_prof_cas(k2) |
---|
| 1867 | thl_mod_cas(l)= frac1*thl_prof_cas(k1) - frac2*thl_prof_cas(k2) |
---|
| 1868 | qv_mod_cas(l)= frac1*qv_prof_cas(k1) - frac2*qv_prof_cas(k2) |
---|
| 1869 | ql_mod_cas(l)= frac1*ql_prof_cas(k1) - frac2*ql_prof_cas(k2) |
---|
| 1870 | qi_mod_cas(l)= frac1*qi_prof_cas(k1) - frac2*qi_prof_cas(k2) |
---|
| 1871 | u_mod_cas(l)= frac1*u_prof_cas(k1) - frac2*u_prof_cas(k2) |
---|
| 1872 | v_mod_cas(l)= frac1*v_prof_cas(k1) - frac2*v_prof_cas(k2) |
---|
| 1873 | ug_mod_cas(l)= frac1*ug_prof_cas(k1) - frac2*ug_prof_cas(k2) |
---|
| 1874 | vg_mod_cas(l)= frac1*vg_prof_cas(k1) - frac2*vg_prof_cas(k2) |
---|
| 1875 | w_mod_cas(l)= frac1*vitw_prof_cas(k1) - frac2*vitw_prof_cas(k2) |
---|
| 1876 | omega_mod_cas(l)= frac1*omega_prof_cas(k1) - frac2*omega_prof_cas(k2) |
---|
| 1877 | du_mod_cas(l)= frac1*du_prof_cas(k1) - frac2*du_prof_cas(k2) |
---|
| 1878 | hu_mod_cas(l)= frac1*hu_prof_cas(k1) - frac2*hu_prof_cas(k2) |
---|
| 1879 | vu_mod_cas(l)= frac1*vu_prof_cas(k1) - frac2*vu_prof_cas(k2) |
---|
| 1880 | dv_mod_cas(l)= frac1*dv_prof_cas(k1) - frac2*dv_prof_cas(k2) |
---|
| 1881 | hv_mod_cas(l)= frac1*hv_prof_cas(k1) - frac2*hv_prof_cas(k2) |
---|
| 1882 | vv_mod_cas(l)= frac1*vv_prof_cas(k1) - frac2*vv_prof_cas(k2) |
---|
| 1883 | dt_mod_cas(l)= frac1*dt_prof_cas(k1) - frac2*dt_prof_cas(k2) |
---|
| 1884 | ht_mod_cas(l)= frac1*ht_prof_cas(k1) - frac2*ht_prof_cas(k2) |
---|
| 1885 | vt_mod_cas(l)= frac1*vt_prof_cas(k1) - frac2*vt_prof_cas(k2) |
---|
| 1886 | dth_mod_cas(l)= frac1*dth_prof_cas(k1) - frac2*dth_prof_cas(k2) |
---|
| 1887 | hth_mod_cas(l)= frac1*hth_prof_cas(k1) - frac2*hth_prof_cas(k2) |
---|
| 1888 | vth_mod_cas(l)= frac1*vth_prof_cas(k1) - frac2*vth_prof_cas(k2) |
---|
| 1889 | dq_mod_cas(l)= frac1*dq_prof_cas(k1) - frac2*dq_prof_cas(k2) |
---|
| 1890 | hq_mod_cas(l)= frac1*hq_prof_cas(k1) - frac2*hq_prof_cas(k2) |
---|
| 1891 | vq_mod_cas(l)= frac1*vq_prof_cas(k1) - frac2*vq_prof_cas(k2) |
---|
[2933] | 1892 | dtrad_mod_cas(l)= frac1*dtrad_prof_cas(k1) - frac2*dtrad_prof_cas(k2) |
---|
[2716] | 1893 | |
---|
| 1894 | endif ! play.le.plev_prof_cas(1) |
---|
| 1895 | |
---|
| 1896 | else ! above max altitude of forcing file |
---|
| 1897 | |
---|
| 1898 | !jyg |
---|
| 1899 | fact=20.*(plev_prof_cas(nlev_cas)-play(l))/plev_prof_cas(nlev_cas) !jyg |
---|
| 1900 | fact = max(fact,0.) !jyg |
---|
| 1901 | fact = exp(-fact) !jyg |
---|
| 1902 | t_mod_cas(l)= t_prof_cas(nlev_cas) !jyg |
---|
| 1903 | theta_mod_cas(l)= th_prof_cas(nlev_cas) !jyg |
---|
| 1904 | thv_mod_cas(l)= thv_prof_cas(nlev_cas) !jyg |
---|
| 1905 | thl_mod_cas(l)= thl_prof_cas(nlev_cas) !jyg |
---|
| 1906 | qv_mod_cas(l)= qv_prof_cas(nlev_cas)*fact !jyg |
---|
| 1907 | ql_mod_cas(l)= ql_prof_cas(nlev_cas)*fact !jyg |
---|
| 1908 | qi_mod_cas(l)= qi_prof_cas(nlev_cas)*fact !jyg |
---|
| 1909 | u_mod_cas(l)= u_prof_cas(nlev_cas)*fact !jyg |
---|
| 1910 | v_mod_cas(l)= v_prof_cas(nlev_cas)*fact !jyg |
---|
| 1911 | ug_mod_cas(l)= ug_prof_cas(nlev_cas)*fact !jyg |
---|
| 1912 | vg_mod_cas(l)= vg_prof_cas(nlev_cas)*fact !jyg |
---|
| 1913 | w_mod_cas(l)= 0.0 !jyg |
---|
[2990] | 1914 | omega_mod_cas(l)= 0.0 !jyg |
---|
[2716] | 1915 | du_mod_cas(l)= du_prof_cas(nlev_cas)*fact |
---|
| 1916 | hu_mod_cas(l)= hu_prof_cas(nlev_cas)*fact !jyg |
---|
| 1917 | vu_mod_cas(l)= vu_prof_cas(nlev_cas)*fact !jyg |
---|
| 1918 | dv_mod_cas(l)= dv_prof_cas(nlev_cas)*fact |
---|
| 1919 | hv_mod_cas(l)= hv_prof_cas(nlev_cas)*fact !jyg |
---|
| 1920 | vv_mod_cas(l)= vv_prof_cas(nlev_cas)*fact !jyg |
---|
| 1921 | dt_mod_cas(l)= dt_prof_cas(nlev_cas) |
---|
| 1922 | ht_mod_cas(l)= ht_prof_cas(nlev_cas) !jyg |
---|
| 1923 | vt_mod_cas(l)= vt_prof_cas(nlev_cas) !jyg |
---|
| 1924 | dth_mod_cas(l)= dth_prof_cas(nlev_cas) |
---|
| 1925 | hth_mod_cas(l)= hth_prof_cas(nlev_cas) !jyg |
---|
| 1926 | vth_mod_cas(l)= vth_prof_cas(nlev_cas) !jyg |
---|
| 1927 | dq_mod_cas(l)= dq_prof_cas(nlev_cas)*fact |
---|
| 1928 | hq_mod_cas(l)= hq_prof_cas(nlev_cas)*fact !jyg |
---|
| 1929 | vq_mod_cas(l)= vq_prof_cas(nlev_cas)*fact !jyg |
---|
[2933] | 1930 | dtrad_mod_cas(l)= dtrad_prof_cas(nlev_cas)*fact !jyg |
---|
[2716] | 1931 | |
---|
| 1932 | endif ! play |
---|
| 1933 | |
---|
| 1934 | enddo ! l |
---|
| 1935 | |
---|
| 1936 | return |
---|
| 1937 | end |
---|
| 1938 | !***************************************************************************** |
---|
| 1939 | |
---|
| 1940 | |
---|
[3008] | 1941 | |
---|
| 1942 | |
---|