- Timestamp:
- Mar 1, 2023, 6:22:39 PM (16 months ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
LMDZ6/branches/Portage_acc/libf/phylmdiso/isotopes_mod.F90
r4124 r4446 3 3 4 4 MODULE isotopes_mod 5 USE infotrac_phy, ONLY: ntraciso,niso,indnum_fn_num,use_iso, & 6 & niso_possibles 7 IMPLICIT NONE 8 SAVE 9 10 ! contient toutes les variables isotopiques et leur initialisation 11 ! les routines specifiquement isotopiques sont dans 12 ! isotopes_routines_mod pour éviter dépendance circulaire avec 13 ! isotopes_verif_mod. 14 15 16 ! indices des isotopes 17 integer, save :: iso_eau,iso_HDO,iso_O18,iso_O17,iso_HTO ! indices de 1 à niso: les isos n'existant pas sont mis à 0 18 !$OMP THREADPRIVATE(iso_eau,iso_HDO,iso_O18,iso_O17,iso_HTO) 19 20 integer :: iso_eau_possible,iso_HDO_possible,iso_O18_possible,iso_O17_possible,iso_HTO_possible ! indices de 1 à niso_possibles: ils correspondent aux tableaux définis dans infotrac: 21 ! tnom_iso=(/'eau','HDO','O18','O17','HTO'/) 22 ! ce sont ces indices qui doivent être utilisés avec use_iso, puisque use_iso est défini comme DIMENSION(niso_possibles) 23 parameter (iso_eau_possible=1) 24 parameter (iso_HDO_possible=2) 25 parameter (iso_O18_possible=3) 26 parameter (iso_O17_possible=4) 27 parameter (iso_HTO_possible=5) 28 29 integer, save :: ntracisoOR 5 USE strings_mod, ONLY: msg, real2str, int2str, bool2str, maxlen, strIdx, strStack 6 USE infotrac_phy, ONLY: isoName 7 IMPLICIT NONE 8 INTERFACE get_in; MODULE PROCEDURE getinp_s, getinp_i, getinp_r, getinp_l; END INTERFACE get_in 9 SAVE 10 11 !--- Contains all isotopic variables + their initialization 12 !--- Isotopes-specific routines are in isotopes_routines_mod to avoid circular dependencies with isotopes_verif_mod. 13 14 !--- Isotopes indices (in [1,niso] ; non-existing => 0 index) 15 INTEGER, SAVE :: iso_eau, iso_HDO, iso_O18, iso_O17, iso_HTO 16 !$OMP THREADPRIVATE(iso_eau, iso_HDO, iso_O18, iso_O17, iso_HTO) 17 18 INTEGER, SAVE :: ntracisoOR 30 19 !$OMP THREADPRIVATE(ntracisoOR) 31 20 32 ! variables indépendantes des isotopes 33 34 real, save :: pxtmelt,pxtice,pxtmin,pxtmax 35 !$OMP THREADPRIVATE(pxtmelt,pxtice,pxtmin,pxtmax) 36 real, save :: tdifexp, tv0cin, thumxt1 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 37 25 !$OMP THREADPRIVATE(tdifexp, tv0cin, thumxt1) 38 integer, save ::ntot26 INTEGER, SAVE :: ntot 39 27 !$OMP THREADPRIVATE(ntot) 40 real, save ::h_land_ice28 REAL, SAVE :: h_land_ice 41 29 !$OMP THREADPRIVATE(h_land_ice) 42 real, save ::P_veg30 REAL, SAVE :: P_veg 43 31 !$OMP THREADPRIVATE(P_veg) 44 real, save :: musi,lambda_sursat45 !$OMP THREADPRIVATE( lambda_sursat)46 real, save ::Kd32 REAL, SAVE :: musi, lambda_sursat 33 !$OMP THREADPRIVATE(musi, lambda_sursat) 34 REAL, SAVE :: Kd 47 35 !$OMP THREADPRIVATE(Kd) 48 real, save :: rh_cste_surf_cond,T_cste_surf_cond 49 !$OMP THREADPRIVATE(rh_cste_surf_cond,T_cste_surf_cond) 50 51 logical, save :: bidouille_anti_divergence 52 ! si true, rappel régulier de xteau vers q, pour éviter dérives lentes 36 REAL, SAVE :: rh_cste_surf_cond, T_cste_surf_cond 37 !$OMP THREADPRIVATE(rh_cste_surf_cond, T_cste_surf_cond) 38 LOGICAL, SAVE :: bidouille_anti_divergence ! T: regularly, xteau <- q to avoid slow drifts 53 39 !$OMP THREADPRIVATE(bidouille_anti_divergence) 54 logical, save :: essai_convergence 55 ! si false, on fait rigoureusement comme dans LMDZ sans isotopes, 56 ! meme si c'est génant pour les isotopes 40 LOGICAL, SAVE :: essai_convergence ! F: as in LMDZ without isotopes (bad for isotopes) 57 41 !$OMP THREADPRIVATE(essai_convergence) 58 integer, save :: initialisation_iso 59 ! 0: dans fichier 60 ! 1: R=0 61 ! 2: R selon distill rayleigh 62 ! 3: R=Rsmow 42 INTEGER, SAVE :: initialisation_iso ! 0: file ; 1: R=0 ; 2: R=distill. Rayleigh ; 3: R=Rsmow 63 43 !$OMP THREADPRIVATE(initialisation_iso) 64 integer, save :: modif_SST ! 0 par defaut, 1 si on veut modifier la sst 65 ! 2 et 3: profils de SST 44 INTEGER, SAVE :: modif_SST ! 0: default ; 1: modified SST ; 2, 3: SST profiles 66 45 !$OMP THREADPRIVATE(modif_SST) 67 real, save :: deltaTtest ! modif de la SST, uniforme. 46 REAL, SAVE :: deltaTtest ! Uniform modification of the SST 68 47 !$OMP THREADPRIVATE(deltaTtest) 69 integer, save :: modif_sic ! on met des trous dans glace de mer 48 INTEGER, SAVE :: modif_sic ! Holes in the Sea Ice 70 49 !$OMP THREADPRIVATE(modif_sic) 71 real, save :: deltasic ! fraction de trous minimale 50 REAL, SAVE :: deltasic ! Minimal holes fraction 72 51 !$OMP THREADPRIVATE(deltasic) 73 real, save ::deltaTtestpoles52 REAL, SAVE :: deltaTtestpoles 74 53 !$OMP THREADPRIVATE(deltaTtestpoles) 75 real, save :: sstlatcrit 76 !$OMP THREADPRIVATE(sstlatcrit) 77 real, save :: dsstlatcrit 78 !$OMP THREADPRIVATE(dsstlatcrit) 79 real, save :: deltaO18_oce 54 REAL, SAVE :: sstlatcrit, dsstlatcrit 55 !$OMP THREADPRIVATE(sstlatcrit, dsstlatcrit) 56 REAL, SAVE :: deltaO18_oce 80 57 !$OMP THREADPRIVATE(deltaO18_oce) 81 integer, save :: albedo_prescrit ! 0 par defaut 82 ! 1 si on veut garder albedo constant 58 INTEGER, SAVE :: albedo_prescrit ! 0: default ; 1: constant albedo 83 59 !$OMP THREADPRIVATE(albedo_prescrit) 84 real, save :: lon_min_albedo,lon_max_albedo 85 !$OMP THREADPRIVATE(lon_min_albedo,lon_max_albedo) 86 real, save :: lat_min_albedo,lat_max_albedo 87 !$OMP THREADPRIVATE(lat_min_albedo,lat_max_albedo) 88 real, save :: deltaP_BL,tdifexp_sol 60 REAL, SAVE :: lon_min_albedo, lon_max_albedo, lat_min_albedo, lat_max_albedo 61 !$OMP THREADPRIVATE(lon_min_albedo, lon_max_albedo, lat_min_albedo, lat_max_albedo) 62 REAL, SAVE :: deltaP_BL,tdifexp_sol 89 63 !$OMP THREADPRIVATE(deltaP_BL,tdifexp_sol) 90 integer, save :: ruissellement_pluie,alphak_stewart91 !$OMP THREADPRIVATE(ruissellement_pluie, alphak_stewart)92 integer, save ::calendrier_guide64 INTEGER, SAVE :: ruissellement_pluie, alphak_stewart 65 !$OMP THREADPRIVATE(ruissellement_pluie, alphak_stewart) 66 INTEGER, SAVE :: calendrier_guide 93 67 !$OMP THREADPRIVATE(calendrier_guide) 94 integer, save ::cste_surf_cond68 INTEGER, SAVE :: cste_surf_cond 95 69 !$OMP THREADPRIVATE(cste_surf_cond) 96 real, save ::mixlen70 REAL, SAVE :: mixlen 97 71 !$OMP THREADPRIVATE(mixlen) 98 integer, save ::evap_cont_cste72 INTEGER, SAVE :: evap_cont_cste 99 73 !$OMP THREADPRIVATE(evap_cont_cste) 100 real, save :: deltaO18_evap_cont,d_evap_cont101 !$OMP THREADPRIVATE(deltaO18_evap_cont, d_evap_cont)102 integer, save :: nudge_qsol,region_nudge_qsol103 !$OMP THREADPRIVATE(nudge_qsol, region_nudge_qsol)104 integer, save:: nlevmaxO1774 REAL, SAVE :: deltaO18_evap_cont, d_evap_cont 75 !$OMP THREADPRIVATE(deltaO18_evap_cont, d_evap_cont) 76 INTEGER, SAVE :: nudge_qsol, region_nudge_qsol 77 !$OMP THREADPRIVATE(nudge_qsol, region_nudge_qsol) 78 INTEGER, SAVE :: nlevmaxO17 105 79 !$OMP THREADPRIVATE(nlevmaxO17) 106 integer, save :: no_pce 107 ! real, save :: slope_limiterxy,slope_limiterz 80 INTEGER, SAVE :: no_pce 108 81 !$OMP THREADPRIVATE(no_pce) 109 real, save ::A_satlim82 REAL, SAVE :: A_satlim 110 83 !$OMP THREADPRIVATE(A_satlim) 111 integer, save :: ok_restrict_A_satlim,modif_ratqs112 !$OMP THREADPRIVATE(ok_restrict_A_satlim, modif_ratqs)113 real, save :: Pcrit_ratqs,ratqsbasnew114 !$OMP THREADPRIVATE(Pcrit_ratqs, ratqsbasnew)115 real, save ::fac_modif_evaoce84 INTEGER, SAVE :: ok_restrict_A_satlim, modif_ratqs 85 !$OMP THREADPRIVATE(ok_restrict_A_satlim, modif_ratqs) 86 REAL, SAVE :: Pcrit_ratqs, ratqsbasnew 87 !$OMP THREADPRIVATE(Pcrit_ratqs, ratqsbasnew) 88 REAL, SAVE :: fac_modif_evaoce 116 89 !$OMP THREADPRIVATE(fac_modif_evaoce) 117 integer, save ::ok_bidouille_wake90 INTEGER, SAVE :: ok_bidouille_wake 118 91 !$OMP THREADPRIVATE(ok_bidouille_wake) 119 logical ::cond_temp_env92 LOGICAL, SAVE :: cond_temp_env 120 93 !$OMP THREADPRIVATE(cond_temp_env) 121 94 122 123 ! variables tableaux fn de niso 124 real, ALLOCATABLE, DIMENSION(:), save :: tnat, toce, tcorr 125 !$OMP THREADPRIVATE(tnat, toce, tcorr) 126 real, ALLOCATABLE, DIMENSION(:), save :: tdifrel 127 !$OMP THREADPRIVATE(tdifrel) 128 real, ALLOCATABLE, DIMENSION(:), save :: talph1, talph2, talph3 129 !$OMP THREADPRIVATE(talph1, talph2, talph3) 130 real, ALLOCATABLE, DIMENSION(:), save :: talps1, talps2 131 !$OMP THREADPRIVATE(talps1, talps2) 132 real, ALLOCATABLE, DIMENSION(:), save :: tkcin0, tkcin1, tkcin2 95 !--- Vectors of length "niso" 96 REAL, ALLOCATABLE, DIMENSION(:), SAVE :: & 97 tnat, toce, tcorr, tdifrel 98 !$OMP THREADPRIVATE(tnat, toce, tcorr, tdifrel) 99 REAL, ALLOCATABLE, DIMENSION(:), SAVE :: & 100 talph1, talph2, talph3, talps1, talps2 101 !$OMP THREADPRIVATE(talph1, talph2, talph3, talps1, talps2) 102 REAL, ALLOCATABLE, DIMENSION(:), SAVE :: & 103 tkcin0, tkcin1, tkcin2 133 104 !$OMP THREADPRIVATE(tkcin0, tkcin1, tkcin2) 134 real, ALLOCATABLE, DIMENSION(:), save :: alpha_liq_sol 135 !$OMP THREADPRIVATE(alpha_liq_sol) 136 real, ALLOCATABLE, DIMENSION(:), save :: Rdefault, Rmethox 137 !$OMP THREADPRIVATE(Rdefault, Rmethox) 138 character*3, ALLOCATABLE, DIMENSION(:), save :: striso 139 !$OMP THREADPRIVATE(striso) 140 real, save :: fac_coeff_eq17_liq, fac_coeff_eq17_ice 105 REAL, ALLOCATABLE, DIMENSION(:), SAVE :: & 106 alpha_liq_sol, Rdefault, Rmethox 107 !$OMP THREADPRIVATE(alpha_liq_sol, Rdefault, Rmethox) 108 REAL, SAVE :: fac_coeff_eq17_liq, fac_coeff_eq17_ice 141 109 !$OMP THREADPRIVATE(fac_coeff_eq17_liq, fac_coeff_eq17_ice) 142 110 143 real ridicule ! valeur maximale pour qu'une variable de type 144 ! rapoport de mélange puisse être considérée comme négligeable. Si 145 ! négligeable, alors on ne verifie pas si sa compo iso esta bérrante. 146 parameter (ridicule=1e-12) 147 ! parameter (ridicule=1) 148 ! 149 real ridicule_rain ! valeur limite de ridicule pour les flux de pluies (rain, zrfl...) 150 parameter (ridicule_rain=1e-8) ! en kg/s <-> 1e-3mm/day 151 152 real ridicule_evap ! valeur limite de ridicule pour les evap 153 parameter (ridicule_evap=ridicule_rain*1e-2) ! en kg/s <-> 1e-3mm/day 154 155 real ridicule_qsol ! valeur limite de ridicule pour les qsol 156 parameter (ridicule_qsol=ridicule_rain) ! en kg <-> 1e-8kg 157 158 real ridicule_snow ! valeur limite de ridicule pour les snow 159 parameter (ridicule_snow=ridicule_qsol) ! en kg/s <-> 1e-8kg 160 161 real expb_max 162 parameter (expb_max=30.0) 163 164 ! spécifique au tritium: 165 166 167 logical, save :: ok_prod_nucl_tritium ! si oui, production de tritium par essais nucleaires 111 !--- Negligible lower thresholds: no need to check for absurd values under these lower limits 112 REAL, PARAMETER :: & 113 ridicule = 1e-12, & ! For mixing ratios 114 ridicule_rain = 1e-8, & ! For rain fluxes (rain, zrfl...) in kg/s <-> 1e-3 mm/day 115 ridicule_evap = ridicule_rain*1e-2, & ! For evaporations in kg/s <-> 1e-3 mm/day 116 ridicule_qsol = ridicule_rain, & ! For qsol in kg <-> 1e-8 kg 117 ridicule_snow = ridicule_qsol ! For snow in kg <-> 1e-8 kg 118 REAL, PARAMETER :: expb_max = 30.0 119 120 !--- Specific to HTO: 121 LOGICAL, SAVE :: ok_prod_nucl_tritium !--- TRUE => HTO production by nuclear tests 168 122 !$OMP THREADPRIVATE(ok_prod_nucl_tritium) 169 integer nessai 170 parameter (nessai=486) 171 integer, save :: day_nucl(nessai) 172 !$OMP THREADPRIVATE(day_nucl) 173 integer, save :: month_nucl(nessai) 174 !$OMP THREADPRIVATE(month_nucl) 175 integer, save :: year_nucl(nessai) 176 !$OMP THREADPRIVATE(year_nucl) 177 real, save :: lat_nucl(nessai) 178 !$OMP THREADPRIVATE(lat_nucl) 179 real, save :: lon_nucl(nessai) 180 !$OMP THREADPRIVATE(lon_nucl) 181 real, save :: zmin_nucl(nessai) 182 !$OMP THREADPRIVATE(zmin_nucl) 183 real, save :: zmax_nucl(nessai) 184 !$OMP THREADPRIVATE(zmax_nucl) 185 real, save :: HTO_nucl(nessai) 186 !$OMP THREADPRIVATE(HTO_nucl) 187 123 INTEGER, PARAMETER :: nessai = 486 124 INTEGER, DIMENSION(nessai), SAVE :: & 125 day_nucl, month_nucl, year_nucl 126 !$OMP THREADPRIVATE(day_nucl, month_nucl, year_nucl) 127 REAL, DIMENSION(nessai), SAVE :: & 128 lat_nucl, lon_nucl, zmin_nucl, zmax_nucl, HTO_nucl 129 !$OMP THREADPRIVATE(lat_nucl, lon_nucl, zmin_nucl, zmax_nucl, HTO_nucl) 130 188 131 189 132 CONTAINS 190 133 191 SUBROUTINE iso_init() 192 use ioipsl_getin_p_mod, ONLY : getin_p 193 implicit none 194 195 ! -- local variables: 196 197 integer ixt 198 ! référence O18 199 real fac_enrichoce18 200 real alpha_liq_sol_O18, & 201 & talph1_O18,talph2_O18,talph3_O18, & 202 & talps1_O18,talps2_O18, & 203 & tkcin0_O18,tkcin1_O18,tkcin2_O18, & 204 & tdifrel_O18 134 SUBROUTINE iso_init() 135 USE infotrac_phy, ONLY: ntiso, niso, getKey 136 USE strings_mod, ONLY: maxlen 137 IMPLICIT NONE 138 139 !=== Local variables: 140 INTEGER :: ixt 141 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 147 !--- For H2[17]O 148 REAL :: fac_kcin, pente_MWL 205 149 206 ! cas de l'O17 207 real fac_kcin 208 real pente_MWL 209 integer ierr 210 211 logical ok_nocinsat, ok_nocinsfc !sensi test 212 parameter (ok_nocinsfc=.FALSE.) ! if T: no kinetic effect in sfc evap 213 parameter (ok_nocinsat=.FALSE.) ! if T: no sursaturation effect for ice 214 logical Rdefault_smow 215 parameter (Rdefault_smow=.FALSE.) ! si T: Rdefault=smow; si F: nul 216 ! pour le tritium 217 integer iessai 218 219 write(*,*) 'iso_init 219: entree' 220 221 ! allocations mémoire 222 allocate (tnat(niso)) 223 allocate (toce(niso)) 224 allocate (tcorr(niso)) 225 allocate (tdifrel(niso)) 226 allocate (talph1(niso)) 227 allocate (talph2(niso)) 228 allocate (talph3(niso)) 229 allocate (talps1(niso)) 230 allocate (talps2(niso)) 231 allocate (tkcin0(niso)) 232 allocate (tkcin1(niso)) 233 allocate (tkcin2(niso)) 234 allocate (alpha_liq_sol(niso)) 235 allocate (Rdefault(niso)) 236 allocate (Rmethox(niso)) 237 allocate (striso(niso)) 238 239 240 !-------------------------------------------------------------- 241 ! General: 242 !-------------------------------------------------------------- 243 244 ! -- verif du nombre d'isotopes: 245 write(*,*) 'iso_init 64: niso=',niso 246 247 ! init de ntracisoOR: on écrasera en cas de nzone>0 si complications avec 248 ! ORCHIDEE 249 ntracisoOR=ntraciso 250 251 ! -- Type of water isotopes: 252 253 iso_eau=indnum_fn_num(1) 254 iso_HDO=indnum_fn_num(2) 255 iso_O18=indnum_fn_num(3) 256 iso_O17=indnum_fn_num(4) 257 iso_HTO=indnum_fn_num(5) 258 write(*,*) 'iso_init 59: iso_eau=',iso_eau 259 write(*,*) 'iso_HDO=',iso_HDO 260 write(*,*) 'iso_O18=',iso_O18 261 write(*,*) 'iso_O17=',iso_O17 262 write(*,*) 'iso_HTO=',iso_HTO 263 write(*,*) 'iso_init 251: use_iso=',use_iso 264 265 ! initialisation 266 lambda_sursat=0.004 267 thumxt1=0.75*1.2 268 ntot=20 269 h_land_ice=20. ! à comparer aux 3000mm de snow_max 270 P_veg=1.0 271 bidouille_anti_divergence=.false. 272 essai_convergence=.false. 273 initialisation_iso=0 274 modif_sst=0 275 modif_sic=0 276 deltaTtest=0.0 277 deltasic=0.1 278 deltaTtestpoles=0.0 279 sstlatcrit=30.0 280 deltaO18_oce=0.0 281 albedo_prescrit=0 282 lon_min_albedo=-200 283 lon_max_albedo=200 284 lat_min_albedo=-100 285 lat_max_albedo=100 286 deltaP_BL=10.0 287 ruissellement_pluie=0 288 alphak_stewart=1 289 tdifexp_sol=0.67 290 calendrier_guide=0 291 cste_surf_cond=0 292 mixlen=35.0 293 evap_cont_cste=0.0 294 deltaO18_evap_cont=0.0 295 d_evap_cont=0.0 296 nudge_qsol=0 297 region_nudge_qsol=1 298 nlevmaxO17=50 299 no_pce=0 300 A_satlim=1.0 301 ok_restrict_A_satlim=0 302 ! slope_limiterxy=2.0 303 ! slope_limiterz=2.0 304 modif_ratqs=0 305 Pcrit_ratqs=500.0 306 ratqsbasnew=0.05 307 308 fac_modif_evaoce=1.0 309 ok_bidouille_wake=0 310 cond_temp_env=.false. 311 ! si oui, la temperature de cond est celle de l'environnement, 312 ! pour eviter bugs quand temperature dans ascendances convs est 313 ! mal calculee 314 ok_prod_nucl_tritium=.false. 315 316 ! lecture des paramètres isotopiques: 317 ! pour que ça marche en openMP, il faut utiliser getin_p. Car le getin ne peut 318 ! être appelé que par un thread à la fois, et ça pose tout un tas de problème, 319 ! d'où tout un tas de magouilles comme dans conf_phys_m. A terme, tout le monde 320 ! lira par getin_p. 321 call getin_p('lambda',lambda_sursat) 322 call getin_p('thumxt1',thumxt1) 323 call getin_p('ntot',ntot) 324 call getin_p('h_land_ice',h_land_ice) 325 call getin_p('P_veg',P_veg) 326 call getin_p('bidouille_anti_divergence',bidouille_anti_divergence) 327 call getin_p('essai_convergence',essai_convergence) 328 call getin_p('initialisation_iso',initialisation_iso) 329 !if (nzone>0) then 330 !if (initialisation_iso.eq.0) then 331 ! call getin_p('initialisation_isotrac',initialisation_isotrac) 332 !endif !if (initialisation_iso.eq.0) then 333 !endif !if (nzone>0) 334 call getin_p('modif_sst',modif_sst) 335 if (modif_sst.ge.1) then 336 call getin_p('deltaTtest',deltaTtest) 337 if (modif_sst.ge.2) then 338 call getin_p('deltaTtestpoles',deltaTtestpoles) 339 call getin_p('sstlatcrit',sstlatcrit) 150 !--- Sensitivity tests 151 LOGICAL, PARAMETER :: ok_nocinsfc = .FALSE. ! if T: no kinetic effect in sfc evap 152 LOGICAL, PARAMETER :: ok_nocinsat = .FALSE. ! if T: no sursaturation effect for ice 153 LOGICAL, PARAMETER :: Rdefault_smow = .FALSE. ! if T: Rdefault=smow; if F: nul 154 155 !--- For [3]H 156 INTEGER :: iessai 157 158 CHARACTER(LEN=maxlen) :: modname, sxt 159 REAL, ALLOCATABLE :: tmp(:) 160 161 modname = 'iso_init' 162 CALL msg('219: entree', modname) 163 164 !-------------------------------------------------------------- 165 ! General: 166 !-------------------------------------------------------------- 167 168 !--- Check number of isotopes 169 CALL msg('64: niso = '//TRIM(int2str(niso)), modname) 170 171 !--- Init de ntracisoOR: on ecrasera en cas de traceurs de tagging isotopiques 172 ! (nzone>0) si complications avec ORCHIDEE 173 ntracisoOR = ntiso 174 175 !--- Type of water isotopes: 176 iso_eau = strIdx(isoName, 'H216O'); CALL msg('iso_eau='//int2str(iso_eau), modname) 177 iso_HDO = strIdx(isoName, 'HDO'); CALL msg('iso_HDO='//int2str(iso_HDO), modname) 178 iso_O18 = strIdx(isoName, 'H218O'); CALL msg('iso_O18='//int2str(iso_O18), modname) 179 iso_O17 = strIdx(isoName, 'H217O'); CALL msg('iso_O17='//int2str(iso_O17), modname) 180 iso_HTO = strIdx(isoName, 'HTO'); CALL msg('iso_HTO='//int2str(iso_HTO), modname) 181 182 !--- Initialiaation: reading the isotopic parameters. 183 CALL get_in('lambda', lambda_sursat, 0.004) 184 CALL get_in('thumxt1', thumxt1, 0.75*1.2) 185 CALL get_in('ntot', ntot, 20, .FALSE.) 186 CALL get_in('h_land_ice', h_land_ice, 20., .FALSE.) 187 CALL get_in('P_veg', P_veg, 1.0, .FALSE.) 188 CALL get_in('bidouille_anti_divergence', bidouille_anti_divergence, .FALSE.) 189 CALL get_in('essai_convergence', essai_convergence, .FALSE.) 190 CALL get_in('initialisation_iso', initialisation_iso, 0) 191 192 ! IF(nzone>0 .AND. initialisation_iso==0) & 193 ! CALL get_in('initialisation_isotrac',initialisation_isotrac) 194 CALL get_in('modif_sst', modif_sst, 0) 195 CALL get_in('deltaTtest', deltaTtest, 0.0) !--- For modif_sst>=1 196 CALL get_in('deltaTtestpoles',deltaTtestpoles, 0.0) !--- For modif_sst>=2 197 CALL get_in( 'sstlatcrit', sstlatcrit, 30.0) !--- For modif_sst>=3 198 CALL get_in('dsstlatcrit', dsstlatcrit, 0.0) !--- For modif_sst>=3 340 199 #ifdef ISOVERIF 341 !call iso_verif_positif(sstlatcrit,'iso_init 107') 342 if (sstlatcrit.lt.0.0) then 343 write(*,*) 'iso_init 270: sstlatcrit=',sstlatcrit 344 stop 345 endif 200 CALL msg('iso_init 270: sstlatcrit='//real2str( sstlatcrit), modname, sstlatcrit < 0.0) !--- For modif_sst>=2 201 CALL msg('iso_init 279: dsstlatcrit='//real2str(dsstlatcrit), modname, sstlatcrit < 0.0) !--- For modif_sst>=3 202 IF(modif_sst >= 2 .AND. sstlatcrit < 0.0) STOP 346 203 #endif 347 if (modif_sst.ge.3) then 348 call getin_p('dsstlatcrit',dsstlatcrit) 204 205 CALL get_in('modif_sic', modif_sic, 0) 206 IF(modif_sic >= 1) & 207 CALL get_in('deltasic', deltasic, 0.1) 208 209 CALL get_in('albedo_prescrit', albedo_prescrit, 0) 210 IF(albedo_prescrit == 1) THEN 211 CALL get_in('lon_min_albedo', lon_min_albedo, -200.) 212 CALL get_in('lon_max_albedo', lon_max_albedo, 200.) 213 CALL get_in('lat_min_albedo', lat_min_albedo, -100.) 214 CALL get_in('lat_max_albedo', lat_max_albedo, 100.) 215 END IF 216 deltaO18_oce=0.0 217 CALL get_in('deltaP_BL', deltaP_BL, 10.0) 218 CALL get_in('ruissellement_pluie', ruissellement_pluie, 0) 219 CALL get_in('alphak_stewart', alphak_stewart, 1) 220 CALL get_in('tdifexp_sol', tdifexp_sol, 0.67) 221 CALL get_in('calendrier_guide', calendrier_guide, 0) 222 CALL get_in('cste_surf_cond', cste_surf_cond, 0) 223 CALL get_in('mixlen', mixlen, 35.0) 224 CALL get_in('evap_cont_cste', evap_cont_cste, 0) 225 CALL get_in('deltaO18_evap_cont', deltaO18_evap_cont,0.0) 226 CALL get_in('d_evap_cont', d_evap_cont, 0.0) 227 CALL get_in('nudge_qsol', nudge_qsol, 0) 228 CALL get_in('region_nudge_qsol', region_nudge_qsol, 1) 229 nlevmaxO17 = 50 230 CALL msg('nlevmaxO17='//TRIM(int2str(nlevmaxO17))) 231 CALL get_in('no_pce', no_pce, 0) 232 CALL get_in('A_satlim', A_satlim, 1.0) 233 CALL get_in('ok_restrict_A_satlim', ok_restrict_A_satlim, 0) 349 234 #ifdef ISOVERIF 350 !call iso_verif_positif(dsstlatcrit,'iso_init 110') 351 if (sstlatcrit.lt.0.0) then 352 write(*,*) 'iso_init 279: dsstlatcrit=',dsstlatcrit 353 stop 354 endif 355 #endif 356 endif !if (modif_sst.ge.3) then 357 endif !if (modif_sst.ge.2) then 358 endif ! if (modif_sst.ge.1) then 359 call getin_p('modif_sic',modif_sic) 360 if (modif_sic.ge.1) then 361 call getin_p('deltasic',deltasic) 362 endif !if (modif_sic.ge.1) then 363 364 call getin_p('albedo_prescrit',albedo_prescrit) 365 call getin_p('lon_min_albedo',lon_min_albedo) 366 call getin_p('lon_max_albedo',lon_max_albedo) 367 call getin_p('lat_min_albedo',lat_min_albedo) 368 call getin_p('lat_max_albedo',lat_max_albedo) 369 call getin_p('deltaO18_oce',deltaO18_oce) 370 call getin_p('deltaP_BL',deltaP_BL) 371 call getin_p('ruissellement_pluie',ruissellement_pluie) 372 call getin_p('alphak_stewart',alphak_stewart) 373 call getin_p('tdifexp_sol',tdifexp_sol) 374 call getin_p('calendrier_guide',calendrier_guide) 375 call getin_p('cste_surf_cond',cste_surf_cond) 376 call getin_p('mixlen',mixlen) 377 call getin_p('evap_cont_cste',evap_cont_cste) 378 call getin_p('deltaO18_evap_cont',deltaO18_evap_cont) 379 call getin_p('d_evap_cont',d_evap_cont) 380 call getin_p('nudge_qsol',nudge_qsol) 381 call getin_p('region_nudge_qsol',region_nudge_qsol) 382 call getin_p('no_pce',no_pce) 383 call getin_p('A_satlim',A_satlim) 384 call getin_p('ok_restrict_A_satlim',ok_restrict_A_satlim) 385 #ifdef ISOVERIF 386 !call iso_verif_positif(1.0-A_satlim,'iso_init 158') 387 if (A_satlim.gt.1.0) then 388 write(*,*) 'iso_init 315: A_satlim=',A_satlim 389 stop 390 endif 391 #endif 392 ! call getin_p('slope_limiterxy',slope_limiterxy) 393 ! call getin_p('slope_limiterz',slope_limiterz) 394 call getin_p('modif_ratqs',modif_ratqs) 395 call getin_p('Pcrit_ratqs',Pcrit_ratqs) 396 call getin_p('ratqsbasnew',ratqsbasnew) 397 call getin_p('fac_modif_evaoce',fac_modif_evaoce) 398 call getin_p('ok_bidouille_wake',ok_bidouille_wake) 399 call getin_p('cond_temp_env',cond_temp_env) 400 if (use_iso(iso_HTO_possible)) then 401 ok_prod_nucl_tritium=.true. 402 call getin_p('ok_prod_nucl_tritium',ok_prod_nucl_tritium) 403 endif 404 405 write(*,*) 'lambda,thumxt1=',lambda_sursat,thumxt1 406 write(*,*) 'bidouille_anti_divergence=',bidouille_anti_divergence 407 write(*,*) 'essai_convergence=',essai_convergence 408 write(*,*) 'initialisation_iso=',initialisation_iso 409 write(*,*) 'modif_sst=',modif_sst 410 if (modif_sst.ge.1) then 411 write(*,*) 'deltaTtest=',deltaTtest 412 if (modif_sst.ge.2) then 413 write(*,*) 'deltaTtestpoles,sstlatcrit=', & 414 & deltaTtestpoles,sstlatcrit 415 if (modif_sst.ge.3) then 416 write(*,*) 'dsstlatcrit=',dsstlatcrit 417 endif !if (modif_sst.ge.3) then 418 endif !if (modif_sst.ge.2) then 419 endif !if (modif_sst.ge.1) then 420 write(*,*) 'modif_sic=',modif_sic 421 if (modif_sic.ge.1) then 422 write(*,*) 'deltasic=',deltasic 423 endif !if (modif_sic.ge.1) then 424 write(*,*) 'deltaO18_oce=',deltaO18_oce 425 write(*,*) 'albedo_prescrit=',albedo_prescrit 426 if (albedo_prescrit.eq.1) then 427 write(*,*) 'lon_min_albedo,lon_max_albedo=', & 428 & lon_min_albedo,lon_max_albedo 429 write(*,*) 'lat_min_albedo,lat_max_albedo=', & 430 & lat_min_albedo,lat_max_albedo 431 endif !if (albedo_prescrit.eq.1) then 432 write(*,*) 'deltaP_BL,ruissellement_pluie,alphak_stewart=', & 433 & deltaP_BL,ruissellement_pluie,alphak_stewart 434 write(*,*) 'cste_surf_cond=',cste_surf_cond 435 write(*,*) 'mixlen=',mixlen 436 write(*,*) 'tdifexp_sol=',tdifexp_sol 437 write(*,*) 'calendrier_guide=',calendrier_guide 438 write(*,*) 'evap_cont_cste=',evap_cont_cste 439 write(*,*) 'deltaO18_evap_cont,d_evap_cont=', & 440 & deltaO18_evap_cont,d_evap_cont 441 write(*,*) 'nudge_qsol,region_nudge_qsol=', & 442 & nudge_qsol,region_nudge_qsol 443 write(*,*) 'nlevmaxO17=',nlevmaxO17 444 write(*,*) 'no_pce=',no_pce 445 write(*,*) 'A_satlim=',A_satlim 446 write(*,*) 'ok_restrict_A_satlim=',ok_restrict_A_satlim 447 ! write(*,*) 'slope_limiterxy=',slope_limiterxy 448 ! write(*,*) 'slope_limiterz=',slope_limiterz 449 write(*,*) 'modif_ratqs=',modif_ratqs 450 write(*,*) 'Pcrit_ratqs=',Pcrit_ratqs 451 write(*,*) 'ratqsbasnew=',ratqsbasnew 452 write(*,*) 'fac_modif_evaoce=',fac_modif_evaoce 453 write(*,*) 'ok_bidouille_wake=',ok_bidouille_wake 454 write(*,*) 'cond_temp_env=',cond_temp_env 455 write(*,*) 'ok_prod_nucl_tritium=',ok_prod_nucl_tritium 456 457 458 !-------------------------------------------------------------- 459 ! Parameters that do not depend on the nature of water isotopes: 460 !-------------------------------------------------------------- 461 462 ! -- temperature at which ice condensate starts to form (valeur ECHAM?): 463 pxtmelt=273.15 464 ! pxtmelt=273.15-10.0 ! test PHASE 465 466 ! -- temperature at which all condensate is ice: 467 pxtice=273.15-10.0 468 ! pxtice=273.15-30.0 ! test PHASE 469 470 ! -- minimum temperature to calculate fractionation coeff 471 pxtmin=273.15-120.0 ! On ne calcule qu'au dessus de -120°C 472 pxtmax=273.15+60.0 ! On ne calcule qu'au dessus de +60°C 473 ! remarque: les coeffs ont été mesurés seulement jusq'à -40! 474 475 ! -- a constant for alpha_eff for equilibrium below cloud base: 476 tdifexp=0.58 477 tv0cin=7.0 478 479 ! facteurs lambda et mu dans Si=musi-lambda*T 480 musi=1.0 481 if (ok_nocinsat) then 482 lambda_sursat = 0.0 ! no sursaturation effect 483 endif 484 485 486 ! diffusion dans le sol 487 Kd=2.5e-9 ! m2/s 488 489 ! cas où cste_surf_cond: on met rhs ou/et Ts cste pour voir 490 rh_cste_surf_cond=0.6 491 T_cste_surf_cond=288.0 492 493 !-------------------------------------------------------------- 494 ! Parameters that depend on the nature of water isotopes: 495 !-------------------------------------------------------------- 496 ! ** constantes locales 497 fac_enrichoce18=0.0005 498 ! on a alors tcor018=1+fac_enrichoce18 499 ! tcorD=1+fac_enrichoce18*8 500 ! tcorO17=1+fac_enrichoce18*0.528 501 alpha_liq_sol_O18=1.00291 ! valeur de Lehmann & Siegenthaler, 1991, 502 ! Journal of Glaciology, vol 37, p 23 503 talph1_O18=1137. 504 talph2_O18=-0.4156 505 talph3_O18=-2.0667E-3 506 talps1_O18=11.839 507 talps2_O18=-0.028244 508 tkcin0_O18 = 0.006 509 tkcin1_O18 = 0.000285 510 tkcin2_O18 = 0.00082 511 tdifrel_O18= 1./0.9723 512 513 ! rapport des ln(alphaeq) entre O18 et O17 514 fac_coeff_eq17_liq=0.529 ! donné par Amaelle 515 ! fac_coeff_eq17_ice=0.528 ! slope MWL 516 fac_coeff_eq17_ice=0.529 517 518 519 write(*,*) 'iso_O18,iso_HDO,iso_eau=',iso_O18,iso_HDO,iso_eau 520 do 999 ixt = 1, niso 521 write(*,*) 'iso_init 80: ixt=',ixt 522 523 524 ! -- kinetic factor for surface evaporation: 525 ! (cf: kcin = tkcin0 if |V|<tv0cin 526 ! kcin = tkcin1*|Vsurf| + tkcin2 if |V|>tv0cin ) 527 ! (Rq: formula discontinuous for |V|=tv0cin... ) 528 529 ! -- main: 530 if (ixt.eq.iso_HTO) then ! Tritium 531 tkcin0(ixt) = 0.01056 532 tkcin1(ixt) = 0.0005016 533 tkcin2(ixt) = 0.0014432 534 tnat(ixt)=0. 535 !toce(ixt)=2.2222E-8 ! corrigé par Alex Cauquoin 536 !toce(ixt)=1.0E-18 ! rapport 3H/1H ocean 537 toce(ixt)=4.0E-19 ! rapport T/H = 0.2 TU Dreisigacker and Roether 1978 538 tcorr(ixt)=1. 539 tdifrel(ixt)=1./0.968 540 talph1(ixt)=46480. 541 talph2(ixt)=-103.87 542 talph3(ixt)=0. 543 talps1(ixt)=46480. 544 talps2(ixt)=-103.87 545 alpha_liq_sol(ixt)=1. 546 Rdefault(ixt)=0.0 547 Rmethox(ixt)=0.0 548 striso(ixt)='HTO' 549 endif 550 if (ixt.eq.iso_O17) then ! Deuterium 551 pente_MWL=0.528 552 ! tdifrel(ixt)=1./0.985452 ! donné par Amaelle 553 tdifrel(ixt)=1./0.98555 ! valeur utilisée en 1D et dans modèle de LdG 554 ! fac_kcin=0.5145 ! donné par Amaelle 555 fac_kcin= (tdifrel(ixt)-1.0)/(tdifrel_O18-1.0) 556 tkcin0(ixt) = tkcin0_O18*fac_kcin 557 tkcin1(ixt) = tkcin1_O18*fac_kcin 558 tkcin2(ixt) = tkcin2_O18*fac_kcin 559 tnat(ixt)=0.004/100. ! O17 représente 0.004% de l'oxygène 560 toce(ixt)=tnat(ixt)*(1.0+deltaO18_oce/1000.0)**pente_MWL 561 tcorr(ixt)=1.0+fac_enrichoce18*pente_MWL ! donné par Amaelle 562 talph1(ixt)=talph1_O18 563 talph2(ixt)=talph2_O18 564 talph3(ixt)=talph3_O18 565 talps1(ixt)=talps1_O18 566 talps2(ixt)=talps2_O18 567 alpha_liq_sol(ixt)=(alpha_liq_sol_O18)**fac_coeff_eq17_liq 568 if (Rdefault_smow) then 569 Rdefault(ixt)=tnat(ixt)*(-3.15/1000.0+1.0) 570 else 571 Rdefault(ixt)=0.0 572 endif 573 Rmethox(ixt)=(230./1000.+1.)*tnat(ixt) !Zahn et al 2006 574 striso(ixt)='O17' 575 endif 576 577 if (ixt.eq.iso_O18) then ! Oxygene18 578 tkcin0(ixt) = tkcin0_O18 579 tkcin1(ixt) = tkcin1_O18 580 tkcin2(ixt) = tkcin2_O18 581 tnat(ixt)=2005.2E-6 582 toce(ixt)=tnat(ixt)*(1.0+deltaO18_oce/1000.0) 583 tcorr(ixt)=1.0+fac_enrichoce18 584 tdifrel(ixt)=tdifrel_O18 585 talph1(ixt)=talph1_O18 586 talph2(ixt)=talph2_O18 587 talph3(ixt)=talph3_O18 588 talps1(ixt)=talps1_O18 589 talps2(ixt)=talps2_O18 590 alpha_liq_sol(ixt)=alpha_liq_sol_O18 591 if (Rdefault_smow) then 592 Rdefault(ixt)=tnat(ixt)*(-6.0/1000.0+1.0) 593 else 594 Rdefault(ixt)=0.0 595 endif 596 Rmethox(ixt)=(130./1000.+1.)*tnat(ixt) !Zahn et al 2006 597 ! write(*,*) 'iso_init 163: ZXalpha_liq_sol=',ZXalpha_liq_sol 598 striso(ixt)='O18' 599 write(*,*) 'isotopes_mod 519: ixt,striso(ixt)=',ixt,striso(ixt) 600 endif 601 602 if (ixt.eq.iso_HDO) then ! Deuterium 603 pente_MWL=8.0 604 ! fac_kcin=0.88 605 tdifrel(ixt)=1./0.9755 606 fac_kcin= (tdifrel(ixt)-1)/(tdifrel_O18-1) 607 tkcin0(ixt) = tkcin0_O18*fac_kcin 608 tkcin1(ixt) = tkcin1_O18*fac_kcin 609 tkcin2(ixt) = tkcin2_O18*fac_kcin 610 tnat(ixt)=155.76E-6 611 toce(ixt)=tnat(ixt)*(1.0+pente_MWL*deltaO18_oce/1000.0) 612 tcorr(ixt)=1.0+fac_enrichoce18*pente_MWL 613 talph1(ixt)=24844. 614 talph2(ixt)=-76.248 615 talph3(ixt)=52.612E-3 616 talps1(ixt)=16288. 617 talps2(ixt)=-0.0934 618 !ZXalpha_liq_sol=1.0192 ! Weston, Ralph, 1955 619 alpha_liq_sol(ixt)=1.0212 620 ! valeur de Lehmann & Siegenthaler, 1991, Journal of 621 ! Glaciology, vol 37, p 23 622 if (Rdefault_smow) then 623 Rdefault(ixt)=tnat(ixt)*((-6.0*pente_MWL+10.0)/1000.0+1.0) 624 else 625 Rdefault(ixt)=0.0 626 endif 627 Rmethox(ixt)=tnat(ixt)*(-25.0/1000.+1.) ! Zahn et al 2006 628 striso(ixt)='HDO' 629 write(*,*) 'isotopes_mod 548: ixt,striso(ixt)=',ixt,striso(ixt) 630 endif 631 632 ! write(*,*) 'iso_init 163: ZXalpha_liq_sol=',ZXalpha_liq_sol 633 if (ixt.eq.iso_eau) then ! Oxygene16 634 tkcin0(ixt) = 0.0 635 tkcin1(ixt) = 0.0 636 tkcin2(ixt) = 0.0 637 tnat(ixt)=1. 638 toce(ixt)=tnat(ixt) 639 tcorr(ixt)=1.0 640 tdifrel(ixt)=1. 641 talph1(ixt)=0. 642 talph2(ixt)=0. 643 talph3(ixt)=0. 644 talps1(ixt)=0. 645 talph3(ixt)=0. 646 alpha_liq_sol(ixt)=1. 647 if (Rdefault_smow) then 648 Rdefault(ixt)=tnat(ixt)*1.0 649 else 650 Rdefault(ixt)=1.0 651 endif 652 Rmethox(ixt)=1.0 653 striso(ixt)='eau' 654 endif 655 656 999 continue 657 658 ! test de sensibilité: 659 if (ok_nocinsfc) then ! no kinetic effect in sfc evaporation 660 do ixt=1,niso 661 tkcin0(ixt) = 0.0 662 tkcin1(ixt) = 0.0 663 tkcin2(ixt) = 0.0 664 enddo 665 endif 666 667 ! nom des isotopes 668 669 ! verif 670 write(*,*) 'iso_init 285: verif initialisation:' 671 672 do ixt=1,niso 673 write(*,*) '* striso(',ixt,')=<'//striso(ixt)//'>' 674 write(*,*) 'tnat(',ixt,')=',tnat(ixt) 675 ! write(*,*) 'alpha_liq_sol(',ixt,')=',alpha_liq_sol(ixt) 676 ! write(*,*) 'tkcin0(',ixt,')=',tkcin0(ixt) 677 ! write(*,*) 'tdifrel(',ixt,')=',tdifrel(ixt) 678 enddo 679 write(*,*) 'iso_init 69: lambda=',lambda_sursat 680 write(*,*) 'iso_init 69: thumxt1=',thumxt1 681 write(*,*) 'iso_init 69: h_land_ice=',h_land_ice 682 write(*,*) 'iso_init 69: P_veg=',P_veg 683 684 return 235 CALL msg(' 315: A_satlim='//real2str(A_satlim), modname, A_satlim > 1.0) 236 IF(A_satlim > 1.0) STOP 237 #endif 238 ! CALL get_in('slope_limiterxy', slope_limiterxy, 2.0) 239 ! CALL get_in('slope_limiterz', slope_limiterz, 2.0) 240 CALL get_in('modif_ratqs', modif_ratqs, 0) 241 CALL get_in('Pcrit_ratqs', Pcrit_ratqs, 500.0) 242 CALL get_in('ratqsbasnew', ratqsbasnew, 0.05) 243 CALL get_in('fac_modif_evaoce', fac_modif_evaoce, 1.0) 244 CALL get_in('ok_bidouille_wake', ok_bidouille_wake, 0) 245 ! si oui, la temperature de cond est celle de l'environnement, pour eviter 246 ! bugs quand temperature dans ascendances convs est mal calculee 247 CALL get_in('cond_temp_env', cond_temp_env, .FALSE.) 248 IF(ANY(isoName == 'HTO')) & 249 CALL get_in('ok_prod_nucl_tritium', ok_prod_nucl_tritium, .FALSE., .FALSE.) 250 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 283 284 CALL msg('iso_O18, iso_HDO, iso_eau = '//TRIM(strStack(int2str([iso_O18, iso_HDO, iso_eau]))), modname) 285 286 !-------------------------------------------------------------- 287 ! Parameters that depend on the nature of water isotopes: 288 !-------------------------------------------------------------- 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) 304 IF(.NOT.Rdefault_smow) Rdefault(:) = 0.0 305 306 !--- Sensitivity test: no kinetic effect in sfc evaporation 307 IF(ok_nocinsfc) THEN 308 tkcin0(1:niso) = 0.0 309 tkcin1(1:niso) = 0.0 310 tkcin2(1:niso) = 0.0 311 END IF 312 313 CALL msg('285: verif initialisation:', modname) 314 DO ixt=1,niso 315 sxt=int2str(ixt) 316 CALL msg(' * isoName('//TRIM(sxt)//') = <'//TRIM(isoName(ixt))//'>', modname) 317 CALL msg( ' tnat('//TRIM(sxt)//') = '//TRIM(real2str(tnat(ixt))), modname) 318 ! CALL msg(' alpha_liq_sol('//TRIM(sxt)//') = '//TRIM(real2str(alpha_liq_sol(ixt))), modname) 319 ! CALL msg( ' tkcin0('//TRIM(sxt)//') = '//TRIM(real2str(tkcin0(ixt))), modname) 320 ! CALL msg( ' tdifrel('//TRIM(sxt)//') = '//TRIM(real2str(tdifrel(ixt))), modname) 321 END DO 322 CALL msg('69: lambda = '//TRIM(real2str(lambda_sursat)), modname) 323 CALL msg('69: thumxt1 = '//TRIM(real2str(thumxt1)), modname) 324 CALL msg('69: h_land_ice = '//TRIM(real2str(h_land_ice)), modname) 325 CALL msg('69: P_veg = '//TRIM(real2str(P_veg)), modname) 326 685 327 END SUBROUTINE iso_init 686 328 329 330 SUBROUTINE getinp_s(nam, val, def, lDisp) 331 USE ioipsl_getincom, ONLY: getin 332 USE mod_phys_lmdz_mpi_data, ONLY : is_mpi_root 333 USE mod_phys_lmdz_omp_data, ONLY : is_omp_root 334 USE mod_phys_lmdz_transfert_para, ONLY : bcast 335 CHARACTER(LEN=*), INTENT(IN) :: nam 336 CHARACTER(LEN=*), INTENT(INOUT) :: val 337 CHARACTER(LEN=*), OPTIONAL, INTENT(IN) :: def 338 LOGICAL, OPTIONAL, INTENT(IN) :: lDisp 339 LOGICAL :: lD 340 !$OMP BARRIER 341 IF(is_mpi_root.AND.is_omp_root) THEN 342 IF(PRESENT(def)) val=def; CALL getin(nam,val) 343 lD=.TRUE.; IF(PRESENT(lDisp)) lD=lDisp 344 IF(lD) CALL msg(TRIM(nam)//' = '//TRIM(val)) 345 END IF 346 CALL bcast(val) 347 END SUBROUTINE getinp_s 348 349 SUBROUTINE getinp_i(nam, val, def, lDisp) 350 USE ioipsl_getincom, ONLY: getin 351 USE mod_phys_lmdz_mpi_data, ONLY : is_mpi_root 352 USE mod_phys_lmdz_omp_data, ONLY : is_omp_root 353 USE mod_phys_lmdz_transfert_para, ONLY : bcast 354 CHARACTER(LEN=*), INTENT(IN) :: nam 355 INTEGER, INTENT(INOUT) :: val 356 INTEGER, OPTIONAL, INTENT(IN) :: def 357 LOGICAL, OPTIONAL, INTENT(IN) :: lDisp 358 LOGICAL :: lD 359 !$OMP BARRIER 360 IF(is_mpi_root.AND.is_omp_root) THEN 361 IF(PRESENT(def)) val=def; CALL getin(nam,val) 362 lD=.TRUE.; IF(PRESENT(lDisp)) lD=lDisp 363 IF(lD) CALL msg(TRIM(nam)//' = '//TRIM(int2str(val))) 364 END IF 365 CALL bcast(val) 366 END SUBROUTINE getinp_i 367 368 SUBROUTINE getinp_r(nam, val, def, lDisp) 369 USE ioipsl_getincom, ONLY: getin 370 USE mod_phys_lmdz_mpi_data, ONLY : is_mpi_root 371 USE mod_phys_lmdz_omp_data, ONLY : is_omp_root 372 USE mod_phys_lmdz_transfert_para, ONLY : bcast 373 CHARACTER(LEN=*), INTENT(IN) :: nam 374 REAL, INTENT(INOUT) :: val 375 REAL, OPTIONAL, INTENT(IN) :: def 376 LOGICAL, OPTIONAL, INTENT(IN) :: lDisp 377 LOGICAL :: lD 378 !$OMP BARRIER 379 IF(is_mpi_root.AND.is_omp_root) THEN 380 IF(PRESENT(def)) val=def; CALL getin(nam,val) 381 lD=.TRUE.; IF(PRESENT(lDisp)) lD=lDisp 382 IF(lD) CALL msg(TRIM(nam)//' = '//TRIM(real2str(val))) 383 END IF 384 CALL bcast(val) 385 END SUBROUTINE getinp_r 386 387 SUBROUTINE getinp_l(nam, val, def, lDisp) 388 USE ioipsl_getincom, ONLY: getin 389 USE mod_phys_lmdz_mpi_data, ONLY : is_mpi_root 390 USE mod_phys_lmdz_omp_data, ONLY : is_omp_root 391 USE mod_phys_lmdz_transfert_para, ONLY : bcast 392 CHARACTER(LEN=*), INTENT(IN) :: nam 393 LOGICAL, INTENT(INOUT) :: val 394 LOGICAL, OPTIONAL, INTENT(IN) :: def 395 LOGICAL, OPTIONAL, INTENT(IN) :: lDisp 396 LOGICAL :: lD 397 !$OMP BARRIER 398 IF(is_mpi_root.AND.is_omp_root) THEN 399 IF(PRESENT(def)) val=def; CALL getin(nam,val) 400 lD=.TRUE.; IF(PRESENT(lDisp)) lD=lDisp 401 IF(lD) CALL msg(TRIM(nam)//' = '//TRIM(bool2str(val))) 402 END IF 403 CALL bcast(val) 404 END SUBROUTINE getinp_l 687 405 688 406 END MODULE isotopes_mod
Note: See TracChangeset
for help on using the changeset viewer.
