Changeset 5202 for LMDZ6/branches/cirrus/libf/phylmdiso/isotopes_mod.F90
- Timestamp:
- Sep 20, 2024, 12:32:04 PM (7 weeks ago)
- Location:
- LMDZ6/branches/cirrus
- Files:
-
- 2 edited
-
. (modified) (1 prop)
-
libf/phylmdiso/isotopes_mod.F90 (modified) (11 diffs)
Legend:
- Unmodified
- Added
- Removed
-
LMDZ6/branches/cirrus
- Property svn:mergeinfo changed
-
LMDZ6/branches/cirrus/libf/phylmdiso/isotopes_mod.F90
r4491 r5202 20 20 21 21 !--- Variables not depending on isotopes 22 REAL, SAVE :: pxtmelt, pxtice, pxtmin, pxtmax 23 !$OMP THREADPRIVATE(pxtmelt, pxtice, pxtmin, pxtmax) 24 REAL, SAVE :: tdifexp, tv0cin, thumxt1 25 !$OMP THREADPRIVATE(tdifexp, tv0cin, thumxt1) 22 REAL, SAVE :: thumxt1 23 !$OMP THREADPRIVATE(thumxt1) 26 24 INTEGER, SAVE :: ntot 27 25 !$OMP THREADPRIVATE(ntot) … … 30 28 REAL, SAVE :: P_veg 31 29 !$OMP THREADPRIVATE(P_veg) 32 REAL, SAVE :: musi, lambda_sursat 33 !$OMP THREADPRIVATE(musi, lambda_sursat) 34 REAL, SAVE :: Kd 35 !$OMP THREADPRIVATE(Kd) 36 REAL, SAVE :: rh_cste_surf_cond, T_cste_surf_cond 37 !$OMP THREADPRIVATE(rh_cste_surf_cond, T_cste_surf_cond) 30 REAL, SAVE :: lambda_sursat 31 !$OMP THREADPRIVATE(lambda_sursat) 38 32 LOGICAL, SAVE :: bidouille_anti_divergence ! T: regularly, xteau <- q to avoid slow drifts 39 33 !$OMP THREADPRIVATE(bidouille_anti_divergence) … … 54 48 REAL, SAVE :: sstlatcrit, dsstlatcrit 55 49 !$OMP THREADPRIVATE(sstlatcrit, dsstlatcrit) 56 REAL, SAVE :: deltaO18_oce57 !$OMP THREADPRIVATE(deltaO18_oce)58 50 INTEGER, SAVE :: albedo_prescrit ! 0: default ; 1: constant albedo 59 51 !$OMP THREADPRIVATE(albedo_prescrit) … … 88 80 REAL, SAVE :: fac_modif_evaoce 89 81 !$OMP THREADPRIVATE(fac_modif_evaoce) 82 REAL, SAVE :: deltaO18_oce 83 !$OMP THREADPRIVATE(deltaO18_oce) 90 84 INTEGER, SAVE :: ok_bidouille_wake 91 85 !$OMP THREADPRIVATE(ok_bidouille_wake) … … 106 100 alpha_liq_sol, Rdefault, Rmethox 107 101 !$OMP THREADPRIVATE(alpha_liq_sol, Rdefault, Rmethox) 108 REAL, SAVE :: fac_coeff_eq17_liq, fac_coeff_eq17_ice 109 !$OMP THREADPRIVATE(fac_coeff_eq17_liq, fac_coeff_eq17_ice) 102 ! REAL, SAVE :: fac_coeff_eq17_liq, fac_coeff_eq17_ice 103 !!$OMP THREADPRIVATE(fac_coeff_eq17_liq, fac_coeff_eq17_ice) 104 105 !--- H2[18]O reference 106 REAL, PARAMETER :: fac_enrichoce18=0.0005 107 REAL, PARAMETER :: alpha_liq_sol_O18=1.00291 108 REAL, PARAMETER :: talph1_O18=1137. 109 REAL, PARAMETER :: talph2_O18=-0.4156 110 REAL, PARAMETER :: talph3_O18=-2.0667E-3 111 REAL, PARAMETER :: talps1_O18=11.839 112 REAL, PARAMETER :: talps2_O18=-0.028244 113 REAL, PARAMETER :: tdifrel_O18=1./0.9723 114 REAL, PARAMETER :: tkcin0_O18=0.006 115 REAL, PARAMETER :: tkcin1_O18=0.000285 116 REAL, PARAMETER :: tkcin2_O18=0.00082 117 REAL, PARAMETER :: fac_coeff_eq17_liq=0.529 118 REAL, PARAMETER :: fac_coeff_eq17_ice=0.529 119 120 !---- Parameters that do not depend on the nature of water isotopes: 121 REAL, PARAMETER :: pxtmelt = 273.15 ! temperature at which ice formation starts 122 REAL, PARAMETER :: pxtice = 273.15-10.0 ! -- temperature at which all condensate is ice: 123 REAL, PARAMETER :: pxtmin = 273.15 - 120.0 ! On ne calcule qu'au dessus de -120°C 124 REAL, PARAMETER :: pxtmax = 273.15 + 60.0 ! On ne calcule qu'au dessus de +60°C 125 REAL, PARAMETER :: tdifexp = 0.58 ! -- a constant for alpha_eff for equilibrium below cloud base: 126 REAL, PARAMETER :: tv0cin = 7.0 ! wind threshold (m/s) for smooth/rough regime (Merlivat and Jouzel 1979) 127 REAL, PARAMETER :: musi=1.0 ! facteurs lambda et mu dans Si=musi-lambda*T 128 REAL, PARAMETER :: Kd=2.5e-9 ! m2/s ! diffusion dans le sol 129 REAL, PARAMETER :: rh_cste_surf_cond = 0.6 ! cas où cste_surf_cond: on met rhs ou/et Ts cste pour voir 130 REAL, PARAMETER :: T_cste_surf_cond = 288.0 131 110 132 111 133 !--- Negligible lower thresholds: no need to check for absurd values under these lower limits … … 140 162 INTEGER :: ixt 141 163 142 !--- H2[18]O reference 143 REAL :: fac_enrichoce18, alpha_liq_sol_O18, & 144 talph1_O18, talph2_O18, talph3_O18, talps1_O18, talps2_O18, & 145 tkcin0_O18, tkcin1_O18, tkcin2_O18, tdifrel_O18 146 164 147 165 !--- For H2[17]O 148 166 REAL :: fac_kcin, pente_MWL … … 152 170 LOGICAL, PARAMETER :: ok_nocinsat = .FALSE. ! if T: no sursaturation effect for ice 153 171 LOGICAL, PARAMETER :: Rdefault_smow = .FALSE. ! if T: Rdefault=smow; if F: nul 172 LOGICAL, PARAMETER :: tnat1 = .TRUE. ! If T: all tnats are 1. 154 173 155 174 !--- For [3]H … … 157 176 158 177 CHARACTER(LEN=maxlen) :: modname, sxt 159 REAL, ALLOCATABLE :: tmp(:)160 178 161 179 modname = 'iso_init' … … 214 232 CALL get_in('lat_max_albedo', lat_max_albedo, 100.) 215 233 END IF 216 deltaO18_oce=0.0217 234 CALL get_in('deltaP_BL', deltaP_BL, 10.0) 218 235 CALL get_in('ruissellement_pluie', ruissellement_pluie, 0) … … 249 266 CALL get_in('ok_prod_nucl_tritium', ok_prod_nucl_tritium, .FALSE., .FALSE.) 250 267 251 !-------------------------------------------------------------- 252 ! Parameters that do not depend on the nature of water isotopes: 253 !-------------------------------------------------------------- 254 ! -- temperature at which ice condensate starts to form (valeur ECHAM?): 255 pxtmelt = 273.15 256 257 ! -- temperature at which all condensate is ice: 258 pxtice = 273.15-10.0 259 260 !- -- test PHASE 261 ! pxtmelt = 273.15 - 10.0 262 ! pxtice = 273.15 - 30.0 263 264 ! -- minimum temperature to calculate fractionation coeff 265 pxtmin = 273.15 - 120.0 ! On ne calcule qu'au dessus de -120°C 266 pxtmax = 273.15 + 60.0 ! On ne calcule qu'au dessus de +60°C 267 ! Remarque: les coeffs ont ete mesures seulement jusq'à -40! 268 269 ! -- a constant for alpha_eff for equilibrium below cloud base: 270 tdifexp = 0.58 271 tv0cin = 7.0 272 273 ! facteurs lambda et mu dans Si=musi-lambda*T 274 musi=1.0 275 if (ok_nocinsat) lambda_sursat = 0.0 ! no sursaturation effect 276 277 ! diffusion dans le sol 278 Kd=2.5e-9 ! m2/s 279 280 ! cas où cste_surf_cond: on met rhs ou/et Ts cste pour voir 281 rh_cste_surf_cond = 0.6 282 T_cste_surf_cond = 288.0 268 ! Ocean composition 269 CALL get_in('deltaO18_oce', deltaO18_oce, 0.0) 283 270 284 271 CALL msg('iso_O18, iso_HDO, iso_eau = '//TRIM(strStack(int2str([iso_O18, iso_HDO, iso_eau]))), modname) 285 272 286 273 !-------------------------------------------------------------- 287 ! Parameters that depend on the nature of water isotopes:274 ! Isotope fractionation factors and a few isotopic constants 288 275 !-------------------------------------------------------------- 289 IF(getKey('tnat', tnat, isoName)) CALL abort_physic(modname, 'can''t get tnat', 1) 290 IF(getKey('toce', toce, isoName)) CALL abort_physic(modname, 'can''t get toce', 1) 291 IF(getKey('tcorr', tcorr, isoName)) CALL abort_physic(modname, 'can''t get tcorr', 1) 292 IF(getKey('talph1', talph1, isoName)) CALL abort_physic(modname, 'can''t get talph1', 1) 293 IF(getKey('talph2', talph2, isoName)) CALL abort_physic(modname, 'can''t get talph2', 1) 294 IF(getKey('talph3', talph3, isoName)) CALL abort_physic(modname, 'can''t get talph3', 1) 295 IF(getKey('talps1', talps1, isoName)) CALL abort_physic(modname, 'can''t get talps1', 1) 296 IF(getKey('talps2', talps2, isoName)) CALL abort_physic(modname, 'can''t get talps2', 1) 297 IF(getKey('tkcin0', tkcin0, isoName)) CALL abort_physic(modname, 'can''t get tkcin0', 1) 298 IF(getKey('tkcin1', tkcin1, isoName)) CALL abort_physic(modname, 'can''t get tkcin1', 1) 299 IF(getKey('tkcin2', tkcin2, isoName)) CALL abort_physic(modname, 'can''t get tkcin2', 1) 300 IF(getKey('tdifrel', tdifrel, isoName)) CALL abort_physic(modname, 'can''t get tdifrel', 1) 301 IF(getKey('alpha_liq_sol', alpha_liq_sol, isoName)) CALL abort_physic(modname, 'can''t get alpha_liq_sol', 1) 302 IF(getKey('Rdefault',Rdefault,isoName)) CALL abort_physic(modname, 'can''t get Rdefault',1) 303 IF(getKey('Rmethox', Rmethox, isoName)) CALL abort_physic(modname, 'can''t get Rmethox', 1) 276 ALLOCATE(tkcin0(niso)) 277 ALLOCATE(tkcin1(niso)) 278 ALLOCATE(tkcin2(niso)) 279 ALLOCATE(tnat(niso)) 280 ALLOCATE(tdifrel(niso)) 281 ALLOCATE(toce(niso)) 282 ALLOCATE(tcorr(niso)) 283 ALLOCATE(talph1(niso)) 284 ALLOCATE(talph2(niso)) 285 ALLOCATE(talph3(niso)) 286 ALLOCATE(talps1(niso)) 287 ALLOCATE(talps2(niso)) 288 ALLOCATE(alpha_liq_sol(niso)) 289 ALLOCATE(Rdefault(niso)) 290 ALLOCATE(Rmethox(niso)) 291 292 do ixt=1,niso 293 if (ixt.eq.iso_HTO) then ! Tritium 294 tkcin0(ixt) = 0.01056 295 tkcin1(ixt) = 0.0005016 296 tkcin2(ixt) = 0.0014432 297 if (tnat1) then 298 tnat(ixt)=1 299 else 300 tnat(ixt)=0. 301 endif 302 toce(ixt)=4.0E-19 ! rapport T/H = 0.2 TU Dreisigacker and Roether 1978 303 tcorr(ixt)=1. 304 tdifrel(ixt)=1./0.968 305 talph1(ixt)=46480. 306 talph2(ixt)=-103.87 307 talph3(ixt)=0. 308 talps1(ixt)=46480. 309 talps2(ixt)=-103.87 310 alpha_liq_sol(ixt)=1. 311 Rmethox(ixt)=0.0 312 endif 313 if (ixt.eq.iso_O17) then ! O17 314 pente_MWL=0.528 315 tdifrel(ixt)=1./0.98555 ! valeur utilisée en 1D et dans modèle de LdG ! tdifrel(ixt)=1./0.985452 ! donné par Amaelle 316 fac_kcin= (tdifrel(ixt)-1.0)/(tdifrel_O18-1.0) ! fac_kcin=0.5145 ! donné par Amaelle 317 tkcin0(ixt) = tkcin0_O18*fac_kcin 318 tkcin1(ixt) = tkcin1_O18*fac_kcin 319 tkcin2(ixt) = tkcin2_O18*fac_kcin 320 if (tnat1) then 321 tnat(ixt)=1 322 else 323 tnat(ixt)=0.004/100. ! O17 représente 0.004% de l'oxygène 324 endif 325 toce(ixt)=tnat(ixt)*(1.0+deltaO18_oce/1000.0)**pente_MWL 326 tcorr(ixt)=1.0+fac_enrichoce18*pente_MWL ! donné par Amaelle 327 talph1(ixt)=talph1_O18 328 talph2(ixt)=talph2_O18 329 talph3(ixt)=talph3_O18 330 talps1(ixt)=talps1_O18 331 talps2(ixt)=talps2_O18 332 alpha_liq_sol(ixt)=(alpha_liq_sol_O18)**fac_coeff_eq17_liq 333 Rdefault(ixt)=tnat(ixt)*(-3.15/1000.0+1.0) 334 Rmethox(ixt)=(230./1000.+1.)*tnat(ixt) !Zahn et al 2006 335 endif 336 if (ixt.eq.iso_O18) then ! Oxygene18 337 tkcin0(ixt) = tkcin0_O18 338 tkcin1(ixt) = tkcin1_O18 339 tkcin2(ixt) = tkcin2_O18 340 if (tnat1) then 341 tnat(ixt)=1 342 else 343 tnat(ixt)=2005.2E-6 344 endif 345 toce(ixt)=tnat(ixt)*(1.0+deltaO18_oce/1000.0) 346 tcorr(ixt)=1.0+fac_enrichoce18 347 tdifrel(ixt)=tdifrel_O18 348 talph1(ixt)=talph1_O18 349 talph2(ixt)=talph2_O18 350 talph3(ixt)=talph3_O18 351 talps1(ixt)=talps1_O18 352 talps2(ixt)=talps2_O18 353 alpha_liq_sol(ixt)=alpha_liq_sol_O18 354 Rdefault(ixt)=tnat(ixt)*(-6.0/1000.0+1.0) 355 Rmethox(ixt)=(130./1000.+1.)*tnat(ixt) !Zahn et al 2006 356 endif 357 if (ixt.eq.iso_HDO) then ! Deuterium 358 pente_MWL=8.0 359 tdifrel(ixt)=1./0.9755 ! fac_kcin=0.88 360 fac_kcin= (tdifrel(ixt)-1)/(tdifrel_O18-1) 361 tkcin0(ixt) = tkcin0_O18*fac_kcin 362 tkcin1(ixt) = tkcin1_O18*fac_kcin 363 tkcin2(ixt) = tkcin2_O18*fac_kcin 364 if (tnat1) then 365 tnat(ixt)=1 366 else 367 tnat(ixt)=155.76E-6 368 endif 369 toce(ixt)=tnat(ixt)*(1.0+pente_MWL*deltaO18_oce/1000.0) 370 tcorr(ixt)=1.0+fac_enrichoce18*pente_MWL 371 talph1(ixt)=24844. 372 talph2(ixt)=-76.248 373 talph3(ixt)=52.612E-3 374 talps1(ixt)=16288. 375 talps2(ixt)=-0.0934 376 !ZXalpha_liq_sol=1.0192 ! Weston, Ralph, 1955 377 alpha_liq_sol(ixt)=1.0212 378 ! valeur de Lehmann & Siegenthaler, 1991, Journal of 379 ! Glaciology, vol 37, p 23 380 Rdefault(ixt)=tnat(ixt)*((-6.0*pente_MWL+10.0)/1000.0+1.0) 381 Rmethox(ixt)=tnat(ixt)*(-25.0/1000.+1.) ! Zahn et al 2006 382 endif 383 if (ixt.eq.iso_eau) then ! Oxygene16 384 tkcin0(ixt) = 0.0 385 tkcin1(ixt) = 0.0 386 tkcin2(ixt) = 0.0 387 tnat(ixt)=1. 388 toce(ixt)=tnat(ixt) 389 tcorr(ixt)=1.0 390 tdifrel(ixt)=1. 391 talph1(ixt)=0. 392 talph2(ixt)=0. 393 talph3(ixt)=0. 394 talps1(ixt)=0. 395 talph3(ixt)=0. 396 alpha_liq_sol(ixt)=1. 397 Rdefault(ixt)=tnat(ixt)*1.0 398 Rmethox(ixt)=1.0 399 endif 400 enddo ! ixt=1,niso 304 401 305 402 IF(.NOT.Rdefault_smow) then … … 308 405 ENDIF 309 406 write(*,*) 'Rdefault=',Rdefault 407 write(*,*) 'toce=',toce 310 408 311 409 !--- Sensitivity test: no kinetic effect in sfc evaporation
Note: See TracChangeset
for help on using the changeset viewer.
