Changeset 5872 for LMDZ6/branches
- Timestamp:
- Nov 18, 2025, 3:41:00 PM (4 days ago)
- Location:
- LMDZ6/branches/lmdz-snow/libf
- Files:
-
- 2 edited
-
phylmd/lmdz_blowing_snow_sublim_sedim.f90 (modified) (3 diffs)
-
phylmdiso/physiq_mod.F90 (modified) (151 diffs)
Legend:
- Unmodified
- Added
- Removed
-
LMDZ6/branches/lmdz-snow/libf/phylmd/lmdz_blowing_snow_sublim_sedim.f90
r5871 r5872 25 25 !===== 26 26 27 integer, intent(in) :: ngrid,nlay 28 real, intent(in) :: dtime 29 real, intent(in), dimension(ngrid) :: ustar 30 real, intent(in), dimension(ngrid,nlay) :: temp 31 real, intent(in), dimension(ngrid,nlay) :: qv 32 real, intent(in), dimension(ngrid,nlay) :: qb 33 real, intent(in), dimension(ngrid,nlay) :: pplay 34 real, intent(in), dimension(ngrid,nlay+1) :: paprs 27 integer, intent(in) :: ngrid ! number of horizontal grid points 28 integer, intent(in) :: nlay ! number of vertical layers 29 real, intent(in) :: dtime ! physics time step [s] 30 real, intent(in), dimension(ngrid) :: ustar ! surface friction velocity [m/s] 31 real, intent(in), dimension(ngrid,nlay) :: temp ! temperature of the air [K] 32 real, intent(in), dimension(ngrid,nlay) :: qv ! specific content of water [kg/kg] 33 real, intent(in), dimension(ngrid,nlay) :: qb ! specific content of blowing snow [kg/kg] 34 real, intent(in), dimension(ngrid,nlay) :: pplay ! pressure at the middle of layers [Pa] 35 real, intent(in), dimension(ngrid,nlay+1) :: paprs ! pressure at the layer bottom interface [Pa] 35 36 36 37 … … 40 41 41 42 42 real, intent(out), dimension(ngrid,nlay) :: dtemp_bs 43 real, intent(out), dimension(ngrid,nlay) :: dqv_bs 44 real, intent(out), dimension(ngrid,nlay) :: dqb_bs 45 real, intent(out), dimension(ngrid,nlay+1) :: bsfl 46 real, intent(out), dimension(ngrid) :: precip_bs 43 real, intent(out), dimension(ngrid,nlay) :: dtemp_bs ! temperature tendency [K/s] 44 real, intent(out), dimension(ngrid,nlay) :: dqv_bs ! water vapor tendency [kg/kg/s] 45 real, intent(out), dimension(ngrid,nlay) :: dqb_bs ! blowing snow mass tendancy [kg/kg/s] 46 real, intent(out), dimension(ngrid,nlay+1) :: bsfl ! vertical profile of blowing snow vertical flux [kg/m2/s] 47 real, intent(out), dimension(ngrid) :: precip_bs ! surface sedimentation flux of blowing snow [kg/s/m2] 47 48 48 49 … … 293 294 294 295 ! 3D variables 296 ! tendencies of qb, qv and temperature 295 297 DO k=1,nlay 296 298 DO i=1,ngrid -
LMDZ6/branches/lmdz-snow/libf/phylmdiso/physiq_mod.F90
r5842 r5872 1 ! 1 2 2 ! $Id: physiq_mod.F90 3908 2021-05-20 07:11:13Z idelkadi $ 3 3 ! … … 18 18 ! For clarity, the "USE" section is now arranged in alphabetical order, 19 19 ! with a separate section for CPP keys 20 ! PLEASE try to follow this rule 20 ! PLEASE try to follow this rule 21 21 22 22 USE ACAMA_GWD_rando_m, only: ACAMA_GWD_rando … … 46 46 USE ioipsl_getin_p_mod, ONLY : getin_p 47 47 USE indice_sol_mod 48 USE infotrac_phy, ONLY: nqtot, nbtr, nqo, tracers, type_trac,addPhase, ivap, iliq, isol, ibs, icf, irvc48 USE infotrac_phy, ONLY: nqtot, nbtr, nqo, nqtke, tracers, type_trac, addPhase, ivap, iliq, isol, ibs, icf, irvc, itke 49 49 USE strings_mod, ONLY: strIdx 50 50 USE iophy … … 102 102 USE calltherm_mod, ONLY : calltherm 103 103 USE lmdz_thermcell_dtke, ONLY : thermcell_dtke 104 USE lmdz_blowing_snow_ini, ONLY : blowing_snow_ini , qbst_bs 104 USE lmdz_blowing_snow_ini, ONLY : blowing_snow_ini , qbst_bs 105 105 USE lmdz_lscp_ini, ONLY : lscp_ini 106 106 USE lmdz_ratqs_main, ONLY : ratqs_main … … 116 116 ! & fl_ebil, fl_cor_ebil 117 117 118 !!!!!!!!!!!!!!!!!! "USE" section for CPP keys !!!!!!!!!!!!!!!!!!!!!!!! 119 ! 120 ! 118 121 119 USE phytracr_spl_mod, ONLY: phytracr_spl, phytracr_spl_out_init 122 120 USE phys_output_write_spl_mod, ONLY: phys_output_write_spl … … 207 205 t_seri,q_seri,ql_seri,qs_seri,qbs_seri, & 208 206 u_seri,v_seri,cf_seri,rvc_seri,tr_seri, & 209 rhcl, & 207 rhcl, & 210 208 qx_seri, & ! CR 211 209 rhcl, & 212 210 ! Dynamic tendencies (diagnostics) 213 211 d_t_dyn,d_q_dyn,d_ql_dyn,d_qs_dyn,d_qbs_dyn, & 214 d_u_dyn,d_v_dyn,d_cf_dyn,d_rvc_dyn,d_t r_dyn, &212 d_u_dyn,d_v_dyn,d_cf_dyn,d_rvc_dyn,d_tke_dyn,d_tr_dyn, & 215 213 d_q_dyn2d,d_ql_dyn2d,d_qs_dyn2d,d_qbs_dyn2d, & 216 214 ! Physic tendencies 217 d_t_con,d_q_con,d_q_con_zmasse,d_u_con,d_v_con, & 215 d_t_con,d_q_con,d_u_con,d_v_con, & 216 d_t_con_zmasse,d_q_con_zmasse,d_u_con_zmasse,d_v_con_zmasse, & 218 217 d_tr, & !! to be removed?? (jyg) 219 218 d_t_wake,d_q_wake, & … … 286 285 topswai_aerop, solswai_aerop, & 287 286 topswad0_aerop, solswad0_aerop, & 288 topsw_aerop, topsw0_aerop, & 287 topsw_aerop, topsw0_aerop, & 289 288 solsw_aerop, solsw0_aerop, & 290 289 topswcf_aerop, solswcf_aerop, & … … 340 339 ! 341 340 wake_k, & 342 alp_wake, & 341 alp_wake, & 343 342 wake_h, wake_omg, & 344 343 ! tendencies of delta T and delta q: … … 352 351 !!! d_s_vdf, d_dens_a_vdf, d_dens_vdf, & ! due to vertical diffusion 353 352 !!! d_s_the, d_dens_a_the, d_dens_the, & ! due to thermals 354 ! 353 ! 355 354 ptconv, ratqsc, & 356 355 wbeff, convoccur, zmax_th, & … … 383 382 ! 384 383 rneblsvol, & 385 pfraclr, pfracld, cldfraliq, sigma2_icefracturb, mean_icefracturb, & 384 pfraclr, pfracld, & 385 cldfraliq, sigma2_icefracturb, mean_icefracturb, & 386 386 cldfraliqth, sigma2_icefracturbth, mean_icefracturbth, & 387 387 distcltop, temp_cltop, & … … 405 405 t2m, fluxlat, & 406 406 fsollw, evap_pot, & 407 fsolsw, wfbils, wfevap, & 407 fsolsw, wfbils, wfevap, & 408 408 prfl, psfl,bsfl, fraca, Vprecip, & 409 409 zw2, & … … 415 415 wwriteSTD, phiwriteSTD, & !pour calcul_STDlev.h 416 416 qwriteSTD, twriteSTD, rhwriteSTD, & !pour calcul_STDlev.h 417 ! 417 ! 418 418 beta_prec, & 419 419 rneb, & … … 462 462 USE conema3_mod_h 463 463 USE alpale_mod 464 USE yoethf_mod_h465 USE calcul_divers_mod_h, ONLY: calcul_divers466 USE compbl_mod_h467 USE nuage_params_mod_h468 USE dimpft_mod_h, ONLY: nvm_lmdz469 USE radepsi_mod_h470 USE radopt_mod_h471 USE regdim_mod_h464 USE yoethf_mod_h 465 USE calcul_divers_mod_h, ONLY: calcul_divers 466 USE compbl_mod_h 467 USE nuage_params_mod_h 468 USE dimpft_mod_h, ONLY: nvm_lmdz 469 USE radepsi_mod_h 470 USE radopt_mod_h 471 USE regdim_mod_h 472 472 IMPLICIT NONE 473 473 !>====================================================================== … … 479 479 !!AA - uniformisation des parametrisations ds phytrac 480 480 !!AA - stockage des moyennes des champs necessaires 481 !!AA en mode traceur off-line 481 !!AA en mode traceur off-line 482 482 !!====================================================================== 483 483 !! CLEFS CPP POUR LES IO … … 543 543 ! Clef iflag_cycle_diurne controlant l'activation du cycle diurne: 544 544 ! en attente du codage des cles par Fred 545 ! iflag_cycle_diurne est initialise par conf_phys et se trouve 545 ! iflag_cycle_diurne est initialise par conf_phys et se trouve 546 546 ! dans clesphys.h (IM) 547 547 !====================================================================== … … 568 568 !$OMP THREADPRIVATE(ok_instan) 569 569 ! 570 LOGICAL ok_LES ! sortir le fichier LES 571 SAVE ok_LES 572 !$OMP THREADPRIVATE(ok_LES) 573 ! 574 LOGICAL callstats ! sortir le fichier stats 575 SAVE callstats 576 !$OMP THREADPRIVATE(callstats) 570 LOGICAL ok_LES ! sortir le fichier LES 571 SAVE ok_LES 572 !$OMP THREADPRIVATE(ok_LES) 573 ! 574 LOGICAL callstats ! sortir le fichier stats 575 SAVE callstats 576 !$OMP THREADPRIVATE(callstats) 577 577 ! 578 578 LOGICAL ok_region ! sortir le fichier regional … … 582 582 SAVE seuil_inversion 583 583 !$OMP THREADPRIVATE(seuil_inversion) 584 585 586 584 585 586 587 587 real facteur 588 588 … … 643 643 !! real wght_cvfd(klon,klev) 644 644 !! ! Variables pour le lessivage convectif 645 !! ! RomP >>> 645 !! ! RomP >>> 646 646 !! real phi2(klon,klev,klev) 647 647 !! real d1a(klon,klev),dam(klon,klev) … … 659 659 INTEGER n 660 660 !ym INTEGER npoints 661 !ym PARAMETER(npoints=klon) 661 !ym PARAMETER(npoints=klon) 662 662 ! 663 663 INTEGER nregISCtot 664 PARAMETER(nregISCtot=1) 664 PARAMETER(nregISCtot=1) 665 665 ! 666 666 ! imin_debut, nbpti, jmin_debut, nbptj : parametres pour sorties … … 671 671 ! direction j (latitude) 672 672 !JLD INTEGER imin_debut, nbpti 673 !JLD INTEGER jmin_debut, nbptj 673 !JLD INTEGER jmin_debut, nbptj 674 674 !IM: region='3d' <==> sorties en global 675 675 CHARACTER*3 region … … 783 783 ! 784 784 !!! INTEGER, SAVE, DIMENSION(klon) :: wake_k 785 !!! !$OMP THREADPRIVATE(wake_k) 785 !!! !$OMP THREADPRIVATE(wake_k) 786 786 #ifdef ISO 787 787 REAL xt_x(ntraciso,klon,klev) … … 790 790 ! 791 791 !jyg< 792 !cc REAL wake_pe(klon) ! Wake potential energy - WAPE 792 !cc REAL wake_pe(klon) ! Wake potential energy - WAPE 793 793 !>jyg 794 794 … … 807 807 REAL d_q_adjwk(klon,klev) !jyg 808 808 LOGICAL,SAVE :: ok_adjwk=.FALSE. 809 !$OMP THREADPRIVATE(ok_adjwk) 809 !$OMP THREADPRIVATE(ok_adjwk) 810 810 INTEGER,SAVE :: iflag_adjwk=0 !jyg 811 811 !$OMP THREADPRIVATE(iflag_adjwk) !jyg 812 812 REAL,SAVE :: oliqmax=999.,oicemax=999. 813 !$OMP THREADPRIVATE(oliqmax,oicemax) 813 !$OMP THREADPRIVATE(oliqmax,oicemax) 814 814 REAL, SAVE :: alp_offset 815 815 !$OMP THREADPRIVATE(alp_offset) … … 830 830 REAL ztv(klon,klev),ztva(klon,klev) 831 831 REAL zpspsk(klon,klev) 832 REAL ztla(klon,klev),zqla(klon,klev) 832 REAL ztla(klon,klev),zqla(klon,klev) 833 833 REAL zthl(klon,klev) 834 834 … … 836 836 837 837 !--------Stochastic Boundary Layer Triggering: ALE_BL-------- 838 !---Propri\'et\'es du thermiques au LCL 838 !---Propri\'et\'es du thermiques au LCL 839 839 ! real zlcl_th(klon) ! Altitude du LCL calcul\'e 840 840 ! continument (pcon dans … … 844 844 real w_conv(klon) ! Vitesse verticale de grande \'echelle au LCL 845 845 real tke0(klon,klev+1) ! TKE au d\'ebut du pas de temps 846 real therm_tke_max0(klon) ! TKE dans les thermiques au LCL 847 real env_tke_max0(klon) ! TKE dans l'environnement au LCL 846 real therm_tke_max0(klon) ! TKE dans les thermiques au LCL 847 real env_tke_max0(klon) ! TKE dans l'environnement au LCL 848 848 INTEGER, SAVE :: iflag_thermcell_tke ! transtport TKE by thermals 849 849 !$OMP THREADPRIVATE(iflag_thermcell_tke) … … 853 853 854 854 REAL,SAVE :: random_notrig_max=1. 855 !$OMP THREADPRIVATE(random_notrig_max) 855 !$OMP THREADPRIVATE(random_notrig_max) 856 856 857 857 !--------Statistical Boundary Layer Closure: ALP_BL-------- … … 862 862 !-------Activer les tendances de TKE due a l'orograp??ie--------- 863 863 INTEGER, SAVE :: addtkeoro 864 !$OMP THREADPRIVATE(addtkeoro) 864 !$OMP THREADPRIVATE(addtkeoro) 865 865 REAL, SAVE :: alphatkeoro 866 !$OMP THREADPRIVATE(alphatkeoro) 866 !$OMP THREADPRIVATE(alphatkeoro) 867 867 LOGICAL, SAVE :: smallscales_tkeoro 868 !$OMP THREADPRIVATE(smallscales_tkeoro) 868 !$OMP THREADPRIVATE(smallscales_tkeoro) 869 869 870 870 … … 881 881 ! 882 882 !AA 883 !AA Pour phytrac 883 !AA Pour phytrac 884 884 REAL u1(klon) ! vents dans la premiere couche U 885 885 REAL v1(klon) ! vents dans la premiere couche V … … 890 890 REAL frac_impa(klon,klev) ! fractions d'aerosols lessivees (impaction) 891 891 REAL frac_nucl(klon,klev) ! idem (nucleation) 892 ! RomP >>> 892 ! RomP >>> 893 893 REAL beta_prec_fisrt(klon,klev) ! taux de conv de l'eau cond (fisrt) 894 894 ! RomP <<< … … 917 917 INTEGER lmt_pas 918 918 SAVE lmt_pas ! frequence de mise a jour 919 !$OMP THREADPRIVATE(lmt_pas) 920 real zmasse(klon, nbp_lev),exner(klon, nbp_lev) 919 !$OMP THREADPRIVATE(lmt_pas) 920 real zmasse(klon, nbp_lev),exner(klon, nbp_lev) 921 921 ! (column-density of mass of air in a cell, in kg m-2) 922 922 real, parameter:: dobson_u = 2.1415e-05 ! Dobson unit, in kg m-2 … … 963 963 REAL radocond(klon,klev) ! eau condensee nuageuse 964 964 ! 965 !XXX PB 965 !XXX PB 966 966 REAL fluxq(klon,klev, nbsrf) ! flux turbulent d'humidite 967 967 REAL fluxqbs(klon,klev, nbsrf) ! flux turbulent de neige soufflee … … 1028 1028 ! La convection n'est pas calculee tous les pas, il faut donc 1029 1029 ! sauvegarder les sorties de la convection 1030 !ym SAVE 1031 !ym SAVE 1032 !ym SAVE 1030 !ym SAVE 1031 !ym SAVE 1032 !ym SAVE 1033 1033 ! 1034 1034 INTEGER itapcv, itapwk … … 1076 1076 #endif 1077 1077 ! 1078 ! Flag pour pouvoir ne pas ajouter les tendances. 1078 ! Flag pour pouvoir ne pas ajouter les tendances. 1079 1079 ! Par defaut, les tendances doivente etre ajoutees et 1080 1080 ! flag_inhib_tend = 0 1081 ! flag_inhib_tend > 0 : tendances non ajoutees, avec un nombre 1081 ! flag_inhib_tend > 0 : tendances non ajoutees, avec un nombre 1082 1082 ! croissant de print quand la valeur du flag augmente 1083 1083 !!! attention, ce flag doit etre change avec prudence !!! … … 1085 1085 !! INTEGER :: flag_inhib_tend = 2 1086 1086 ! 1087 ! Logical switch to a bug : reseting to 0 convective variables at the 1087 ! Logical switch to a bug : reseting to 0 convective variables at the 1088 1088 ! begining of physiq. 1089 1089 LOGICAL, SAVE :: ok_bug_cv_trac = .TRUE. … … 1095 1095 !$OMP THREADPRIVATE(ok_bug_split_th) 1096 1096 1097 ! Logical switch to a bug : modifying directly wake_deltat by adding 1097 ! Logical switch to a bug : modifying directly wake_deltat by adding 1098 1098 ! the (w) dry adjustment tendency to wake_deltat 1099 1099 LOGICAL, SAVE :: ok_bug_ajs_cv = .TRUE. 1100 1100 !$OMP THREADPRIVATE(ok_bug_ajs_cv) 1101 1101 1102 ! 1102 1103 !******************************************************** … … 1175 1176 ! 1176 1177 REAL zx_tmp_fi2d(klon) ! variable temporaire grille physique 1177 REAL zx_tmp_fi3d(klon,klev) ! variable temporaire pour champs 3D 1178 REAL zx_tmp_fi3d(klon,klev) ! variable temporaire pour champs 3D 1178 1179 !JLD REAL zx_tmp_2d(nbp_lon,nbp_lat) 1179 1180 !JLD REAL zx_lon(nbp_lon,nbp_lat) … … 1206 1207 LOGICAL, SAVE :: ok_sync, ok_sync_omp 1207 1208 !$OMP THREADPRIVATE(ok_sync) 1209 ! ok_sync_omp should not be in a THREADPRIVATE statement 1208 1210 REAL date0 1209 1211 1210 1212 ! essai writephys 1211 1213 INTEGER fid_day, fid_mth, fid_ins 1212 PARAMETER (fid_ins = 1, fid_day = 2, fid_mth = 3) 1214 PARAMETER (fid_ins = 1, fid_day = 2, fid_mth = 3) 1213 1215 INTEGER prof2d_on, prof3d_on, prof2d_av, prof3d_av 1214 1216 PARAMETER (prof2d_on = 1, prof3d_on = 2, prof2d_av = 3, prof3d_av = 4) … … 1229 1231 INTEGER :: naero 1230 1232 ! Aerosol optical properties 1231 CHARACTER*4, DIMENSION(naero_grp) :: rfname 1233 CHARACTER*4, DIMENSION(naero_grp) :: rfname 1232 1234 REAL, DIMENSION(klon,klev) :: mass_solu_aero ! total mass 1233 1235 ! concentration … … 1247 1249 LOGICAL, SAVE :: aerosol_couple ! true : calcul des aerosols dans INCA 1248 1250 ! false : lecture des aerosol dans un fichier 1249 !$OMP THREADPRIVATE(aerosol_couple) 1250 LOGICAL, SAVE :: chemistry_couple ! true : use INCA chemistry O3 1251 !$OMP THREADPRIVATE(aerosol_couple) 1252 LOGICAL, SAVE :: chemistry_couple ! true : use INCA chemistry O3 1251 1253 ! false : use offline chemistry O3 1252 !$OMP THREADPRIVATE(chemistry_couple) 1253 INTEGER, SAVE :: flag_aerosol 1254 !$OMP THREADPRIVATE(flag_aerosol) 1254 !$OMP THREADPRIVATE(chemistry_couple) 1255 INTEGER, SAVE :: flag_aerosol 1256 !$OMP THREADPRIVATE(flag_aerosol) 1255 1257 LOGICAL, SAVE :: flag_bc_internal_mixture 1256 1258 !$OMP THREADPRIVATE(flag_bc_internal_mixture) … … 1284 1286 REAL, ALLOCATABLE, SAVE :: press_edg_climoz(:) ! Edges of pressure intervals 1285 1287 REAL, ALLOCATABLE, SAVE :: time_climoz(:) ! Time vector 1288 1286 1289 CHARACTER(LEN=13), PARAMETER :: vars_climoz(2) & 1287 1290 = ["tro3 ","tro3_daylight"] … … 1292 1295 1293 1296 include "FCTTRE.h" 1294 1297 !IM 100106 END : pouvoir sortir les ctes de la physique 1298 ! 1295 1299 ! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 1296 1300 ! Declarations pour Simulateur COSP … … 1306 1310 INTEGER, ALLOCATABLE, SAVE :: tabCFMIP(:) 1307 1311 REAL, ALLOCATABLE, SAVE :: lonCFMIP(:), latCFMIP(:) 1308 !$OMP THREADPRIVATE(tabCFMIP, lonCFMIP, latCFMIP) 1312 !$OMP THREADPRIVATE(tabCFMIP, lonCFMIP, latCFMIP) 1309 1313 INTEGER, ALLOCATABLE, SAVE :: tabijGCM(:) 1310 1314 REAL, ALLOCATABLE, SAVE :: lonGCM(:), latGCM(:) … … 1345 1349 1346 1350 !lwoff=y : offset LW CRE for radiation code and other schemes 1347 REAL, SAVE :: betalwoff 1348 ! OMP THREADPRIVATE(betalwoff)1351 REAL, SAVE :: betalwoff 1352 !$OMP THREADPRIVATE(betalwoff) 1349 1353 ! 1350 1354 INTEGER :: nbtr_tmp ! Number of tracer inside concvl … … 1388 1392 CHARACTER(len=512) :: namelist_ecrad_file 1389 1393 1394 1390 1395 ! Subgrid scale wind : 1391 1396 ! Need to be allocatable/save because the number of bin is not known (provided by surf_wind_ini) 1392 integer, save :: nsurfwind=1 1397 integer, save :: nsurfwind=1 1393 1398 real, dimension(:,:), allocatable, save :: surf_wind_value, surf_wind_proba ! module and probability of sugrdi wind wind sample 1394 1399 !$OMP THREADPRIVATE(nsurfwind,surf_wind_value, surf_wind_proba) 1400 1395 1401 1396 1402 !======================================================================! … … 1403 1409 call getin_p('iflag_physiq', iflag_physiq) ! 1404 1410 endif ! 1405 if ( iflag_physiq == 2 ) then 1411 if ( iflag_physiq == 2 ) then ! 1406 1412 #ifdef ISO 1407 1413 abort_message='isotopes pas encore dans physiqex' … … 1423 1429 ! set-up call to alerte function 1424 1430 call_alert = (alert_first_call .AND. is_master) 1425 1431 1426 1432 ! Ehouarn: set value of jjmp1 since it is no longer a "fixed parameter" 1427 1433 jjmp1=nbp_lat … … 1436 1442 1437 1443 IF (using_xios) THEN 1438 ! switch to XIOS LMDZ physics context 1444 ! switch to XIOS LMDZ physics context 1439 1445 IF (.NOT. debut .AND. is_omp_master) THEN 1440 1446 CALL wxios_set_context() … … 1456 1462 'nlon,klev,nqtot,debut,lafin, jD_cur, jH_cur,pdtphys' 1457 1463 write(lunout,*) & 1458 nlon,klev,nqtot,debut,lafin, jD_cur, jH_cur,pdtphys 1464 nlon,klev,nqtot,debut,lafin, jD_cur, jH_cur,pdtphys 1459 1465 1460 1466 write(lunout,*) 'paprs, play, phi, u, v, t' … … 1473 1479 "physiq_mod paprs bad order") 1474 1480 1475 IF (first) THEN 1481 IF (first) THEN 1476 1482 ! CALL init_etat0_limit_unstruct 1477 1483 IF (.NOT. create_etat0_limit) CALL init_limit_read(days_elapsed) … … 1502 1508 CALL phys_state_var_init(read_climoz) 1503 1509 CALL phys_output_var_init 1504 IF (read_climoz>=1 .AND. create_etat0_limit .AND. grid_type==unstructured) & 1510 IF (read_climoz>=1 .AND. create_etat0_limit .AND. grid_type==unstructured) & 1505 1511 CALL regr_horiz_time_climoz(read_climoz,ok_daily_climoz) 1506 1512 … … 1553 1559 1554 1560 IF (ok_bs) THEN 1555 #ifdef ISO 1556 abort_message='blowing snow cannot be activated with water isotopes yet' 1557 CALL abort_physic(modname,abort_message, 1) 1558 #endif 1559 IF ((ok_ice_supersat.AND.nqo .LT.6).OR.(.NOT.ok_ice_supersat.AND.nqo.LT.5)) THEN 1561 IF ((ok_ice_supersat.AND.nqo .LT.6).OR.(.NOT.ok_ice_supersat.AND.nqo.LT.4)) THEN 1560 1562 WRITE (lunout, *) 'activation of blowing snow needs a specific H2O tracer', & 1561 1563 'but nqo=', nqo … … 1564 1566 ENDIF 1565 1567 ENDIF 1568 1569 IF (ok_advtke) THEN 1570 IF (nqtke .LT. 1) THEN 1571 WRITE (lunout, *) 'activation of TKE advection need a specific TKE tracer', & 1572 'but nqtke=', nqtke 1573 abort_message='see above' 1574 CALL abort_physic(modname,abort_message, 1) 1575 ENDIF 1576 ENDIF 1577 1566 1578 1567 1579 Ncvpaseq1 = 0 … … 1579 1591 ENDIF ! first 1580 1592 1581 !ym => necessaire pour iflag_con != 2 1593 !ym => necessaire pour iflag_con != 2 1582 1594 pmfd(:,:) = 0. 1583 1595 pen_u(:,:) = 0. … … 1609 1621 call getin_p('iflag_thermcell_tke', iflag_thermcell_tke) ! 1610 1622 1611 CALL getin_p('iflag_alp_wk_cond', iflag_alp_wk_cond) 1623 CALL getin_p('iflag_alp_wk_cond', iflag_alp_wk_cond) 1612 1624 CALL getin_p('random_notrig_max',random_notrig_max) 1613 CALL getin_p('ok_adjwk',ok_adjwk) 1625 CALL getin_p('ok_adjwk',ok_adjwk) 1614 1626 IF (ok_adjwk) iflag_adjwk=2 ! for compatibility with older versions 1615 1627 ! iflag_adjwk: ! 0 = Default: no convective adjustment of w-region … … 1652 1664 WRITE(lunout,*) 'ok_adjwk=', ok_adjwk 1653 1665 WRITE(lunout,*) 'iflag_adjwk=', iflag_adjwk 1654 WRITE(lunout,*) 'qtcon_multistep_max=',dtcon_multistep_max 1655 WRITE(lunout,*) 'qdcon_multistep_max=',dqcon_multistep_max 1666 WRITE(lunout,*) 'qtcon_multistep_max=',dtcon_multistep_max 1667 WRITE(lunout,*) 'qdcon_multistep_max=',dqcon_multistep_max 1656 1668 WRITE(lunout,*) 'ratqsp0=', ratqsp0 1657 WRITE(lunout,*) 'ratqsdp=', ratqsdp 1669 WRITE(lunout,*) 'ratqsdp=', ratqsdp 1658 1670 WRITE(lunout,*) 'iflag_wake_tend=', iflag_wake_tend 1659 WRITE(lunout,*) 'ok_bad_ecmwf_thermo=',ok_bad_ecmwf_thermo 1660 WRITE(lunout,*) 'ok_bug_cv_trac=', ok_bug_cv_trac 1671 WRITE(lunout,*) 'ok_bad_ecmwf_thermo=',ok_bad_ecmwf_thermo 1672 WRITE(lunout,*) 'ok_bug_cv_trac=', ok_bug_cv_trac 1661 1673 WRITE(lunout,*) 'ok_bug_split_th=', ok_bug_split_th 1662 1674 WRITE(lunout,*) 'fl_ebil=', fl_ebil … … 1694 1706 !des caracteristiques du thermique 1695 1707 wght_th(:,:)=1. 1696 lalim_conv(:)=1 1708 lalim_conv(:)=1 1697 1709 !RC 1698 1710 ustar(:,:)=0. … … 1720 1732 CALL getin_p('config_inca',config_inca) 1721 1733 1722 ELSE 1734 ELSE 1723 1735 config_inca='none' ! default 1724 1736 ENDIF … … 1800 1812 1801 1813 IF (iflag_pbl>1) THEN 1802 PRINT*, "Using method MELLOR&YAMADA" 1814 PRINT*, "Using method MELLOR&YAMADA" 1803 1815 ENDIF 1804 1816 … … 1818 1830 IF (MOD(NINT(86400./phys_tstep),nbapp_rad).EQ.0) THEN 1819 1831 radpas = NINT( 86400./phys_tstep)/nbapp_rad 1820 ELSE 1832 ELSE 1821 1833 WRITE(lunout,*) 'le nombre de pas de temps physique doit etre un ', & 1822 1834 'multiple de nbapp_rad' … … 1834 1846 cvpas = cvpas_0 1835 1847 print *,'physiq, cvpas ',cvpas 1836 ELSE 1848 ELSE 1837 1849 WRITE(lunout,*) 'le nombre de pas de temps physique doit etre un ', & 1838 1850 'multiple de nbapp_cv' … … 1846 1858 wkpas = NINT( 86400./phys_tstep)/nbapp_wk 1847 1859 ! print *,'physiq, wkpas ',wkpas 1848 ELSE 1860 ELSE 1849 1861 WRITE(lunout,*) 'le nombre de pas de temps physique doit etre un ', & 1850 1862 'multiple de nbapp_wk' … … 1882 1894 ALLOCATE(iGCM(nCFMIP), jGCM(nCFMIP)) 1883 1895 ! 1884 ! lecture des nCFMIP stations CFMIP, de leur numero 1896 ! lecture des nCFMIP stations CFMIP, de leur numero 1885 1897 ! et des coordonnees geographiques lonCFMIP, latCFMIP 1886 1898 ! … … 1930 1942 ! WRITE(lunout,*)'OK freq_outNMC(3)=',freq_outNMC(3) 1931 1943 1932 1933 1944 #ifndef CPP_XIOS 1934 1945 CALL ini_paramLMDZ_phy(phys_tstep,nid_ctesGCM) … … 1939 1950 1940 1951 !XXXPB Positionner date0 pour initialisation de ORCHIDEE 1941 date0 = jD_ref 1952 date0 = jD_ref 1942 1953 WRITE(*,*) 'physiq date0 : ',date0 1943 1954 ! … … 2000 2011 CALL iniradia(klon,klev,paprs(1,1:klev+1)) 2001 2012 2002 2013 2003 2014 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 2004 2015 CALL surf_wind_ini(klon,lunout) 2005 2016 CALL getin_p('nsurfwind',nsurfwind) 2006 2017 allocate(surf_wind_value(klon,nsurfwind),surf_wind_proba(klon,nsurfwind)) 2007 2018 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 2019 ! initialisation des variables tkeoro 2020 addtkeoro=0 2021 CALL getin_p('addtkeoro',addtkeoro) 2022 2023 alphatkeoro=1. 2024 CALL getin_p('alphatkeoro',alphatkeoro) 2025 alphatkeoro=min(max(0.,alphatkeoro),1.) 2026 2027 smallscales_tkeoro=.FALSE. 2028 CALL getin_p('smallscales_tkeoro',smallscales_tkeoro) 2029 2008 2030 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 2009 2031 CALL wake_ini(iflag_wake,rg,rd,rv,prt_level) … … 2025 2047 IF (ok_cdnc.AND.NRADLP.NE.3) THEN 2026 2048 abort_message='RRTM choix incoherent NRADLP doit etre egal a 3 ' & 2027 // 'pour ok_cdnc' 2049 // 'pour ok_cdnc' 2028 2050 CALL abort_physic(modname,abort_message,1) 2029 2051 ENDIF … … 2034 2056 #endif 2035 2057 ENDIF 2036 ENDIF 2058 ENDIF 2037 2059 CALL cloud_optics_prop_ini(klon, klev, prt_level, lunout, flag_aerosol, & 2038 2060 & ok_cdnc, bl95_b0, & … … 2059 2081 flag_aerosol, flag_aerosol_strat, ok_cdnc) 2060 2082 ELSE 2061 ! phys_output_write écrit des variables traceurs seulement si iflag_phytrac == 1 2083 ! phys_output_write écrit des variables traceurs seulement si iflag_phytrac == 1 2062 2084 ! donc seulement dans ce cas on doit appeler phytrac_init() 2063 2085 IF (iflag_phytrac == 1 ) THEN … … 2076 2098 IF (is_omp_master) CALL xios_update_calendar(1) 2077 2099 ENDIF 2078 2100 2079 2101 IF(read_climoz>=1 .AND. create_etat0_limit) CALL regr_horiz_time_climoz(read_climoz,ok_daily_climoz) 2080 2102 CALL create_etat0_limit_unstruct … … 2138 2160 ENDIF 2139 2161 ! 2140 IF (phys_tstep*REAL(radpas).GT.21600..AND.iflag_cycle_diurne.GE.1) THEN 2162 IF (phys_tstep*REAL(radpas).GT.21600..AND.iflag_cycle_diurne.GE.1) THEN 2141 2163 WRITE(lunout,*)'Nbre d appels au rayonnement insuffisant' 2142 2164 WRITE(lunout,*)"Au minimum 4 appels par jour si cycle diurne" … … 2267 2289 !============================================================= 2268 2290 2269 IF (using_xios) THEN 2291 IF (using_xios) THEN 2270 2292 ! Get "missing_val" value from XML files (from temperature variable) 2271 2293 IF (is_omp_master) CALL xios_get_field_attr("temp",default_value=missing_val) … … 2273 2295 ENDIF 2274 2296 2275 IF (using_xios) THEN 2297 IF (using_xios) THEN 2276 2298 ! Need to put this initialisation after phyetat0 as in the coupled model the XIOS context is only 2277 2299 ! initialised at that moment … … 2279 2301 IF (is_omp_master) CALL xios_get_field_attr("temp",default_value=missing_val) 2280 2302 CALL bcast_omp(missing_val) 2281 ! 2303 ! 2282 2304 ! Now we activate some double radiation call flags only if some 2283 2305 ! diagnostics are requested, otherwise there is no point in doing this 2284 2306 IF (is_master) THEN 2285 !--setting up swaero_diag to TRUE in XIOS case 2286 IF (xios_field_is_active("topswad").OR.xios_field_is_active("topswad0").OR. & 2287 xios_field_is_active("solswad").OR.xios_field_is_active("solswad0").OR. & 2288 xios_field_is_active("topswai").OR.xios_field_is_active("solswai").OR. & 2307 !--setting up swaero_diag to TRUE in XIOS case 2308 IF (xios_field_is_active("topswad").OR.xios_field_is_active("topswad0").OR. & 2309 xios_field_is_active("solswad").OR.xios_field_is_active("solswad0").OR. & 2310 xios_field_is_active("topswai").OR.xios_field_is_active("solswai").OR. & 2289 2311 (iflag_rrtm==1.AND.(xios_field_is_active("toplwad").OR.xios_field_is_active("toplwad0").OR. & 2290 2312 xios_field_is_active("sollwad").OR.xios_field_is_active("sollwad0")))) & 2291 !!!--for now these fields are not in the XML files so they are omitted 2292 !!! xios_field_is_active("toplwai").OR.xios_field_is_active("sollwai") !))) & 2293 swaero_diag=.TRUE. 2294 2313 !!!--for now these fields are not in the XML files so they are omitted 2314 !!! xios_field_is_active("toplwai").OR.xios_field_is_active("sollwai") !))) & 2315 swaero_diag=.TRUE. 2295 2316 !--setting up swaerofree_diag to TRUE in XIOS case 2296 2317 IF (xios_field_is_active("SWdnSFCcleanclr").OR.xios_field_is_active("SWupSFCcleanclr").OR. & … … 2298 2319 xios_field_is_active("rsdcsaf") .OR. xios_field_is_active("LWdnSFCcleanclr").OR. & 2299 2320 xios_field_is_active("LWupTOAcleanclr")) & 2300 swaerofree_diag=.TRUE. 2301 2302 !--setting up dryaod_diag to TRUE in XIOS case 2321 swaerofree_diag=.TRUE. 2322 2323 !--setting up dryaod_diag to TRUE in XIOS case 2303 2324 DO naero = 1, naero_tot-1 2304 IF (xios_field_is_active("dryod550_"//name_aero_tau(naero))) dryaod_diag=.TRUE. 2325 IF (xios_field_is_active("dryod550_"//name_aero_tau(naero))) dryaod_diag=.TRUE. 2305 2326 ENDDO 2306 2327 ! 2307 !--setting up ok_4xCO2atm to TRUE in XIOS case 2308 IF (xios_field_is_active("rsut4co2").OR.xios_field_is_active("rlut4co2").OR. & 2328 !--setting up ok_4xCO2atm to TRUE in XIOS case 2329 IF (xios_field_is_active("rsut4co2").OR.xios_field_is_active("rlut4co2").OR. & 2309 2330 xios_field_is_active("rsutcs4co2").OR.xios_field_is_active("rlutcs4co2").OR. & 2310 2331 xios_field_is_active("rsu4co2").OR.xios_field_is_active("rsucs4co2").OR. & … … 2312 2333 xios_field_is_active("rlu4co2").OR.xios_field_is_active("rlucs4co2").OR. & 2313 2334 xios_field_is_active("rld4co2").OR.xios_field_is_active("rldcs4co2")) & 2314 ok_4xCO2atm=.TRUE. 2335 ok_4xCO2atm=.TRUE. 2315 2336 ENDIF 2316 2337 !$OMP BARRIER … … 2494 2515 ENDIF 2495 2516 WRITE(*,*)'ok_lwoff=',ok_lwoff 2496 ! 2517 ! 2497 2518 !lwoff=y to begin only sollw and sollwdown are set up to CS values 2498 2519 sollw = sollw + betalwoff * (sollw0 - sollw) … … 2522 2543 ! 2523 2544 ! 2524 ! Update fraction of the sub-surfaces (pctsrf) and 2525 ! initialize, where a new fraction has appeared, all variables depending 2545 ! Update fraction of the sub-surfaces (pctsrf) and 2546 ! initialize, where a new fraction has appeared, all variables depending 2526 2547 ! on the surface fraction. 2527 2548 ! … … 2583 2604 beta_prec(:,:)=0. 2584 2605 ! 2585 ! Output variables from the convective scheme should not be set to 0 2606 ! Output variables from the convective scheme should not be set to 0 2586 2607 ! since convection is not always called at every time step. 2587 2608 IF (ok_bug_cv_trac) THEN … … 2631 2652 ENDDO 2632 2653 2654 ! in case of TKE advection, we interpolate vertically 2655 ! since TKE is defined at the bottom interface of layers but 2656 ! it is interpolated onto middle layers for advection 2657 2658 IF (ok_advtke) THEN 2659 DO k=2,klev 2660 DO i=1,klon 2661 pbl_tke(i,k,:)=(qx(i,k-1,itke)*zmasse(i,k-1)+qx(i,k,itke)*zmasse(i,k))/(zmasse(i,k-1)+zmasse(i,k)) 2662 ENDDO 2663 ENDDO 2664 ENDIF 2665 2666 tke0(:,:)=pbl_tke(:,:,is_ave) 2667 2668 2633 2669 ! C Risi: dispatcher les isotopes dans les xt_seri 2634 2670 #ifdef ISO … … 2695 2731 !--fin mass fixer 2696 2732 2697 tke0(:,:)=pbl_tke(:,:,is_ave) 2698 IF (nqtot > nqo) THEN 2733 IF (nqtot > (nqo+nqtke)) THEN 2699 2734 ! water isotopes are not included in tr_seri 2700 2735 itr = 0 … … 2721 2756 IF(.NOT.tracers(iq)%isInPhysics) CYCLE 2722 2757 itr = itr+1 2723 tr_ancien(:,:,itr)=tr_seri(:,:,itr) 2758 tr_ancien(:,:,itr)=tr_seri(:,:,itr) 2724 2759 enddo 2725 2760 ENDIF … … 2794 2829 d_cf_dyn(:,:) = (cf_seri(:,:)-cf_ancien(:,:))/phys_tstep 2795 2830 d_rvc_dyn(:,:)= (rvc_seri(:,:)-rvc_ancien(:,:))/phys_tstep 2831 d_tke_dyn(:,:)= (pbl_tke(:,:,is_ave)-tke_ancien(:,:))/phys_tstep 2796 2832 CALL water_int(klon,klev,q_seri,zmasse,zx_tmp_fi2d) 2797 2833 d_q_dyn2d(:)=(zx_tmp_fi2d(:)-prw_ancien(:))/phys_tstep … … 2803 2839 d_qbs_dyn2d(:)=(zx_tmp_fi2d(:)-prbsw_ancien(:))/phys_tstep 2804 2840 ! !! RomP >>> td dyn traceur 2805 IF (nqtot > nqo) d_tr_dyn(:,:,:)=(tr_seri(:,:,:)-tr_ancien(:,:,:))/phys_tstep2841 IF (nqtot > (nqo+nqtke)) d_tr_dyn(:,:,:)=(tr_seri(:,:,:)-tr_ancien(:,:,:))/phys_tstep 2806 2842 ! !! RomP <<< 2807 2843 … … 2889 2925 d_cf_dyn(:,:) = 0.0 2890 2926 d_rvc_dyn(:,:)= 0.0 2927 d_tke_dyn(:,:)= 0.0 2891 2928 d_q_dyn2d(:) = 0.0 2892 2929 d_ql_dyn2d(:) = 0.0 … … 2913 2950 2914 2951 ! !! RomP >>> td dyn traceur 2915 IF (nqtot > nqo) d_tr_dyn(:,:,:)= 0.02952 IF (nqtot > (nqo+nqtke)) d_tr_dyn(:,:,:)= 0.0 2916 2953 ! !! RomP <<< 2917 2954 ancien_ok = .TRUE. … … 3145 3182 day_since_equinox = (jD_cur + jH_cur) - jD_eq 3146 3183 ! 3147 ! choix entre calcul de la longitude solaire vraie ou valeur fixee a 3184 ! choix entre calcul de la longitude solaire vraie ou valeur fixee a 3148 3185 ! solarlong0 3149 3186 IF (solarlong0<-999.) THEN … … 3157 3194 ELSE 3158 3195 zlongi=solarlong0 ! longitude solaire vraie 3159 dist=1. ! distance au soleil / moyenne 3196 dist=1. ! distance au soleil / moyenne 3160 3197 ENDIF 3161 3198 … … 3179 3216 ! recode par Olivier Boucher en sept 2015 3180 3217 SELECT CASE (iflag_cycle_diurne) 3181 CASE(0) 3218 CASE(0) 3182 3219 ! Sans cycle diurne 3183 3220 CALL angle(zlongi, latitude_deg, fract, rmu0) … … 3185 3222 JrNt = 1.0 3186 3223 zrmu0 = rmu0 3187 CASE(1) 3224 CASE(1) 3188 3225 ! Avec cycle diurne sans application des poids 3189 3226 ! bit comparable a l ancienne formulation cycle_diurne=true … … 3197 3234 JrNt = 0.0 3198 3235 WHERE (fract.GT.0.0) JrNt = 1.0 3199 CASE(2) 3236 CASE(2) 3200 3237 ! Avec cycle diurne sans application des poids 3201 3238 ! On integre entre gmtime-pdtphys et gmtime+pdtphys*(radpas-1) … … 3206 3243 ! premier pas de temps de la physique pendant lequel 3207 3244 ! itaprad=0 3208 zdtime1=phys_tstep*REAL(-MOD(itaprad,radpas)-1) 3209 zdtime2=phys_tstep*REAL(radpas-MOD(itaprad,radpas)-1) 3245 zdtime1=phys_tstep*REAL(-MOD(itaprad,radpas)-1) 3246 zdtime2=phys_tstep*REAL(radpas-MOD(itaprad,radpas)-1) 3210 3247 CALL zenang(zlongi,jH_cur,zdtime1,zdtime2, & 3211 3248 latitude_deg,longitude_deg,rmu0,fract) … … 3222 3259 ! Calcul du flag jour-nuit 3223 3260 JrNt = 0.0 3224 WHERE (zfract.GT.0.0) JrNt = 1.0 3261 WHERE (zfract.GT.0.0) JrNt = 1.0 3225 3262 END SELECT 3226 3263 ENDIF … … 3238 3275 ! Cela implique tous les interactions des sous-surfaces et la 3239 3276 ! partie diffusion turbulent du couche limit. 3240 ! 3241 ! Certains varibales de sorties de pbl_surface sont utiliser que pour 3277 ! 3278 ! Certains varibales de sorties de pbl_surface sont utiliser que pour 3242 3279 ! ecriture des fihiers hist_XXXX.nc, ces sont : 3243 3280 ! qsol, zq2m, s_pblh, s_lcl, … … 3251 3288 ! wfbils, fluxt, fluxu, fluxv, 3252 3289 ! 3253 ! Certains ne sont pas utiliser du tout : 3290 ! Certains ne sont pas utiliser du tout : 3254 3291 ! dsens, devap, zxsnow, zxfluxt, zxfluxq, q2m, fluxq 3255 3292 ! … … 3297 3334 !ym : Warning gustiness non inialized for iflag_gusts=2 & iflag_gusts=3 3298 3335 gustiness=0 !ym missing init 3299 3336 3300 3337 IF (iflag_gusts==0) THEN 3301 3338 gustiness(1:klon)=0 … … 3344 3381 cdragh, cdragm, u1, v1, & 3345 3382 beta_aridity, & 3346 !albedo SB >>> 3347 ! albsol1, albsol2, sens, evap, & 3348 albsol_dir, albsol_dif, sens, evap, snowerosion, icesub_lic, & 3349 !albedo SB <<< 3383 albsol_dir, albsol_dif, sens, evap, snowerosion, icesub_lic, & 3350 3384 albsol3_lic,runoff, snowhgt, qsnow, to_ice, sissnow, & 3351 3385 zxtsol, zxfluxlat, zt2m, qsat2m, zn2mout, & … … 3373 3407 z0m, z0h, agesno, fsollw, fsolsw, & 3374 3408 d_ts, fevap, fluxlat, t2m, & 3375 wfbils, wfevap, & 3409 wfbils, wfevap, & 3376 3410 fluxt, fluxu, fluxv, & 3377 3411 dsens, devap, zxsnow, & … … 3575 3609 IF (ok_bs) THEN 3576 3610 3577 CALL call_blowing_snow_sublim_sedim(klon,klev,phys_tstep, t_seri,q_seri,qbs_seri,pplay,paprs, &3611 CALL call_blowing_snow_sublim_sedim(klon,klev,phys_tstep,zustar,t_seri,q_seri,qbs_seri,pplay,paprs, & 3578 3612 d_t_bsss,d_q_bsss,d_qbs_bsss,bsfl,bs_fall) 3579 3613 3580 3614 CALL add_phys_tend & 3581 3615 (du0,dv0,d_t_bsss,d_q_bsss,dql0,dqi0,d_qbs_bsss,paprs,& 3582 'bs ',abortphy,flag_inhib_tend,itap,0 &3616 'bsss',abortphy,flag_inhib_tend,itap,0 & 3583 3617 #ifdef ISO 3584 3618 & ,dxt0,dxtl0,dxti0 & … … 3891 3925 !>jyg 3892 3926 ! 3893 3927 3894 3928 !! print *,'physiq. q_w(1,k), q_x(1,k) ', & 3895 3929 !! (k, q_w(1,k), q_x(1,k),k=1,25) … … 4093 4127 ! 4094 4128 !jyg< 4095 ! If convective tendencies are too large, then call convection 4129 ! If convective tendencies are too large, then call convection 4096 4130 ! every time step 4097 4131 cvpas = cvpas_0 … … 4248 4282 ENDIF 4249 4283 4250 !--saving d_ q_con * zmass for next timestep if convection is not called every timestep4251 IF (ok_ conserv_d_q_con) THEN4284 !--saving d_X_con * zmass for next timestep if convection is not called every timestep 4285 IF (ok_mass_dqcon) THEN 4252 4286 d_q_con_zmasse(:,:) = d_q_con(:,:) * zmasse(:,:) 4253 4287 ENDIF 4288 4289 IF (ok_mass_dtcon) THEN 4290 d_t_con_zmasse(:,:) = d_t_con(:,:) * zmasse(:,:) 4291 ENDIF 4292 4293 IF (ok_mass_duvcon) THEN 4294 d_u_con_zmasse(:,:) = d_u_con(:,:) * zmasse(:,:) 4295 d_v_con_zmasse(:,:) = d_v_con(:,:) * zmasse(:,:) 4296 ENDIF 4297 4254 4298 4255 4299 ! CALL homogene(paprs, q_seri, d_q_con, u_seri,v_seri, … … 4271 4315 ENDIF 4272 4316 4273 !!!jyg Appel diagnostique a add_phys_tend pour tester la conservation de 4317 !!!jyg Appel diagnostique a add_phys_tend pour tester la conservation de 4274 4318 !!! l'energie dans les courants satures. 4275 4319 !! d_t_con_sat(:,:) = d_t_con(:,:) - ftd(:,:)*dtime … … 4277 4321 !! dql_sat(:,:) = (wdtrainA(:,:)+wdtrainM(:,:))*dtime/zmasse(:,:) 4278 4322 !! CALL add_phys_tend(d_u_con, d_v_con, d_t_con_sat, d_q_con_sat, dql_sat, & 4279 !! dqi0, paprs, 'convection_sat', abortphy, flag_inhib_tend,& 4323 !! dqi0, paprs, 'convection_sat', abortphy, flag_inhib_tend,& 4280 4324 !! itap, 1) 4281 4325 !! call prt_enerbil('convection_sat',itap) … … 4283 4327 !! 4284 4328 4285 !--recompute d_ q_con with zmasse from new timestep4286 IF (ok_ conserv_d_q_con) THEN4329 !--recompute d_X_con with zmasse from new timestep 4330 IF (ok_mass_dqcon) THEN 4287 4331 d_q_con(:,:)=d_q_con_zmasse(:,:)/zmasse(:,:) 4288 4332 ENDIF 4333 4334 IF (ok_mass_dtcon) THEN 4335 d_t_con(:,:)=d_t_con_zmasse(:,:)/zmasse(:,:) 4336 ENDIF 4337 4338 IF (ok_mass_duvcon) THEN 4339 d_u_con(:,:)=d_u_con_zmasse(:,:)/zmasse(:,:) 4340 d_v_con(:,:)=d_v_con_zmasse(:,:)/zmasse(:,:) 4341 ENDIF 4342 4343 4344 4289 4345 4290 4346 CALL add_phys_tend(d_u_con, d_v_con, d_t_con, d_q_con, dql0, dqi0, dqbs0, paprs, & … … 4394 4450 ! 4395 4451 !========================================================================== 4396 !RR:Evolution de la poche froide: on ne fait pas de separation wake/env 4452 !RR:Evolution de la poche froide: on ne fait pas de separation wake/env 4397 4453 !pour la couche limite diffuse pour l instant 4398 4454 ! … … 4411 4467 DO k=1,klev 4412 4468 DO i=1,klon 4413 dt_dwn(i,k) = ftd(i,k) 4414 dq_dwn(i,k) = fqd(i,k) 4469 dt_dwn(i,k) = ftd(i,k) 4470 dq_dwn(i,k) = fqd(i,k) 4415 4471 M_dwn(i,k) = dnwd0(i,k) 4416 4472 M_up(i,k) = upwd(i,k) 4417 dt_a(i,k) = d_t_con(i,k)/phys_tstep - ftd(i,k) 4473 dt_a(i,k) = d_t_con(i,k)/phys_tstep - ftd(i,k) 4418 4474 dq_a(i,k) = d_q_con(i,k)/phys_tstep - fqd(i,k) 4419 4475 #ifdef ISO … … 4425 4481 ENDDO 4426 4482 ENDDO 4427 4483 4428 4484 IF (iflag_wake==2) THEN 4429 4485 ok_wk_lsp(:)=max(sign(1.,wake_s(:)-wake_s_min_lsp),0.) … … 4462 4518 ENDDO 4463 4519 ENDIF 4464 4520 4465 4521 #ifdef ISOVERIF 4466 4522 write(*,*) 'physiq 3977: verif des inputs de calwake' … … 4733 4789 !jyg< 4734 4790 !! IF (mod(iflag_pbl_split/2,2) .EQ. 1) THEN 4735 IF (mod(iflag_pbl_split/10,10) .GE. 1) THEN4791 IF (mod(iflag_pbl_split/10,10) .GE. 1) THEN 4736 4792 ! Appel des thermiques avec les profils exterieurs aux poches 4737 4793 DO k=1,klev … … 4823 4879 DO k=1,klev 4824 4880 DO i=1,klon 4825 !4826 4881 d_deltat_the(i,k) = - d_t_ajs(i,k) 4827 4882 d_deltaq_the(i,k) = - d_q_ajs(i,k) 4828 !4829 4883 d_u_ajs(i,k) = d_u_ajs(i,k)*(1.-wake_s(i)) 4830 4884 d_v_ajs(i,k) = d_v_ajs(i,k)*(1.-wake_s(i)) … … 4887 4941 ! zmax_th(i)=pphi(i,lmax_th(i))/rg 4888 4942 !CR:04/05/12:correction calcul zmax 4889 zmax_th(i)=zmax0(i) 4943 zmax_th(i)=zmax0(i) 4890 4944 ENDDO 4891 4945 … … 4910 4964 ! Attention : le call ajsec_convV2 n'est maintenu que momentanneement 4911 4965 ! pour des test de convergence numerique. 4912 ! Le nouveau ajsec est a priori mieux, meme pour le cas 4966 ! Le nouveau ajsec est a priori mieux, meme pour le cas 4913 4967 ! iflag_thermals = 0 (l'ancienne version peut faire des tendances 4914 4968 ! non nulles numeriquement pour des mailles non concernees. … … 4994 5048 4995 5049 ENDIF 4996 4997 5050 ! 4998 5051 !=================================================================== … … 5000 5053 ! Developed for dust lifting. Could be extended to coupling with ocean and others 5001 5054 ! by default : 1 bin equal to the mean wind 5002 5055 5003 5056 call surf_wind(klon,nsurfwind,zu10m,zv10m,wake_s,wake_Cstar,zustar,ale_bl,surf_wind_value,surf_wind_proba) 5004 5005 5006 ! 5057 5007 5058 !=================================================================== 5008 ! Computation of ratqs, the width (normalized) of the subrid scale 5059 ! Computation of ratqs, the width (normalized) of the subrid scale 5009 5060 ! water distribution 5010 5061 … … 5274 5325 radocond(i,k)=radocond(i,k)+qbs_seri(i,k) 5275 5326 picefra(i,k)=(radocond(i,k)*picefra(i,k)+qbs_seri(i,k))/(radocond(i,k)) 5276 qbsfra=min(qbs_seri(i,k)/qbst_bs,1.0) 5327 qbsfra=min(qbs_seri(i,k)/qbst_bs,1.0) 5277 5328 cldfra(i,k)=max(cldfra(i,k),qbsfra) 5278 5329 ENDDO … … 5509 5560 ! profonde. 5510 5561 5511 !IM/FH: 2011/02/23 5562 !IM/FH: 2011/02/23 5512 5563 ! definition des points sur lesquels ls thermiques sont actifs 5513 5564 … … 5659 5710 ENDDO 5660 5711 5661 !IM Calcul temp.potentielle a 2m (tpot) et temp. potentielle 5712 !IM Calcul temp.potentielle a 2m (tpot) et temp. potentielle 5662 5713 ! equivalente a 2m (tpote) pour diagnostique 5663 5714 ! … … 5795 5846 5796 5847 ! 5797 ELSE IF (NSW.EQ.2) THEN 5848 ELSE IF (NSW.EQ.2) THEN 5798 5849 !--for now we use the old aerosol properties 5799 5850 ! … … 5895 5946 ENDIF 5896 5947 ELSE 5897 tausum_aero(:,:,id_STRAT_phy) = 0. 5948 tausum_aero(:,:,id_STRAT_phy) = 0. 5898 5949 ENDIF 5899 5950 ! … … 5911 5962 #endif 5912 5963 !--fin STRAT AEROSOL 5913 ! 5964 ! 5914 5965 5915 5966 ! Calculer les parametres optiques des nuages et quelques 5916 5967 ! parametres pour diagnostiques: 5917 5968 ! 5918 IF (aerosol_couple.AND.config_inca=='aero') THEN 5919 mass_solu_aero(:,:) = ccm(:,:,1) 5920 mass_solu_aero_pi(:,:) = ccm(:,:,2) 5969 IF (aerosol_couple.AND.config_inca=='aero') THEN 5970 mass_solu_aero(:,:) = ccm(:,:,1) 5971 mass_solu_aero_pi(:,:) = ccm(:,:,2) 5921 5972 ENDIF 5922 5973 … … 5929 5980 cldtaupi, distcltop, temp_cltop, re, fl, ref_liq, ref_ice, & 5930 5981 ref_liq_pi, ref_ice_pi, scdnc, cldncl, reffclwtop, lcc, reffclws, & 5931 reffclwc, cldnvi, lcc3d, lcc3dcon, lcc3dstra, icc3dcon, icc3dstra, & 5982 reffclwc, cldnvi, lcc3d, lcc3dcon, lcc3dstra, icc3dcon, icc3dstra, & 5932 5983 zfice, dNovrN, ptconv, rnebcon, clwcon) 5933 5984 … … 5996 6047 ENDIF 5997 6048 5998 !lecture de la chlorophylle pour le nouvel albedo de Sunghye Baek 6049 !lecture de la chlorophylle pour le nouvel albedo de Sunghye Baek 5999 6050 IF (ok_chlorophyll) THEN 6000 6051 print*,"-- reading chlorophyll" … … 6002 6053 ENDIF 6003 6054 6004 !--if ok_suntime_rrtm we use ancillay data for RSUN 6055 !--if ok_suntime_rrtm we use ancillay data for RSUN 6005 6056 !--previous values are therefore overwritten 6006 6057 !--this is needed for CMIP6 runs 6007 6058 !--and only possible for new radiation scheme 6008 IF (iflag_rrtm.EQ.1.AND.ok_suntime_rrtm) THEN 6059 IF (iflag_rrtm.EQ.1.AND.ok_suntime_rrtm) THEN 6009 6060 #ifdef CPP_RRTM 6010 6061 CALL read_rsun_rrtm(debut) … … 6028 6079 ENDIF 6029 6080 6030 IF (aerosol_couple.AND.config_inca=='aero') THEN 6081 IF (aerosol_couple.AND.config_inca=='aero') THEN 6031 6082 IF (CPPKEY_INCA) THEN 6032 6083 CALL radlwsw_inca & … … 6059 6110 RCFC11 = RCFC11_act 6060 6111 RCFC12 = RCFC12_act 6061 ! 6112 ! 6062 6113 !--interactive CO2 in ppm from carbon cycle 6063 6114 IF (carbon_cycle_rad) RCO2=RCO2_glo … … 6081 6132 flag_aerosol, flag_aerosol_strat, flag_aer_feedback, & 6082 6133 tau_aero, piz_aero, cg_aero, & 6083 tau_aero_sw_rrtm, piz_aero_sw_rrtm, cg_aero_sw_rrtm, & 6134 tau_aero_sw_rrtm, piz_aero_sw_rrtm, cg_aero_sw_rrtm, & 6084 6135 ! Rajoute par OB pour RRTM 6085 tau_aero_lw_rrtm, & 6136 tau_aero_lw_rrtm, & 6086 6137 cldtaupirad, m_allaer, & 6087 6138 ! zqsat, flwcrad, fiwcrad, & … … 6122 6173 sollwdown(:)) 6123 6174 cool = cool + betalwoff * (cool0 - cool) 6124 6175 6125 6176 IF (.NOT. using_xios) THEN 6126 6177 ! … … 6133 6184 RN2O_per.NE.RN2O_act.OR. & 6134 6185 RCFC11_per.NE.RCFC11_act.OR. & 6135 RCFC12_per.NE.RCFC12_act) ok_4xCO2atm =.TRUE. 6186 RCFC12_per.NE.RCFC12_act) ok_4xCO2atm =.TRUE. 6136 6187 ENDIF 6137 ! 6188 ! 6138 6189 IF (ok_4xCO2atm) THEN 6139 6190 ! … … 6155 6206 !albedo SB >>> 6156 6207 ! paprs, pplay,zxtsol,albsol1, albsol2, & 6157 paprs, pplay,zxtsol,SFRWL,albsol_dir, albsol_dif, & 6208 paprs, pplay,zxtsol,SFRWL,albsol_dir, albsol_dif, & 6158 6209 !albedo SB <<< 6159 6210 t_seri,q_seri,wo, & … … 6201 6252 #ifdef CPP_ECRAD 6202 6253 IF (ok_3Deffect) then 6203 ! print*,'ok_3Deffect = ',ok_3Deffect 6254 ! print*,'ok_3Deffect = ',ok_3Deffect 6204 6255 namelist_ecrad_file='namelist_ecrad_s2' 6205 6256 CALL radlwsw & … … 6217 6268 ref_liq, ref_ice, ref_liq_pi, ref_ice_pi, & 6218 6269 namelist_ecrad_file, & 6219 ! A modifier 6270 ! A modifier 6220 6271 heat_s2,heat0_s2,cool_s2,cool0_s2,albpla_s2, & 6221 6272 heat_volc,cool_volc, & … … 6338 6389 ! a l'echelle sous-maille: 6339 6390 ! 6340 6341 ! calculation of nm_oro 6391 6392 ! calculation of nm_oro 6342 6393 DO i=1,klon 6343 6394 ! nm_oro is a proxy for the number of subgrid scale mountains … … 6347 6398 ! nm_oro_t=0. 6348 6399 nm_oro(i)=zsig(i)*sqrt(cell_area(i)*(pctsrf(i,is_ter)+pctsrf(i,is_lic)))/(4.*MAX(zstd(i),1.e-8))-1. 6349 END 6400 ENDDO 6350 6401 6351 6402 IF (prt_level .GE.10) THEN 6352 6403 print *,' call orography ? ', ok_orodr 6353 6404 ENDIF 6354 6405 ! 6355 6406 IF (ok_orodr) THEN 6356 6407 ! … … 6546 6597 !IM calcul composantes axiales du moment angulaire et couple des montagnes 6547 6598 ! 6548 IF (is_sequential .and. ok_orodr) THEN 6599 IF (is_sequential .and. ok_orodr) THEN 6549 6600 CALL aaam_bud (27,klon,klev,jD_cur-jD_ref,jH_cur, & 6550 6601 ra,rg,romega, & … … 6639 6690 ! Additional tendency of TKE due to orography 6640 6691 !=============================================================== 6641 !6642 ! Inititialization6643 !------------------6644 6645 addtkeoro=06646 CALL getin_p('addtkeoro',addtkeoro)6647 6648 IF (prt_level.ge.5) &6649 print*,'addtkeoro', addtkeoro6650 6651 alphatkeoro=1.6652 CALL getin_p('alphatkeoro',alphatkeoro)6653 alphatkeoro=min(max(0.,alphatkeoro),1.)6654 6655 smallscales_tkeoro=.FALSE.6656 CALL getin_p('smallscales_tkeoro',smallscales_tkeoro)6657 6658 6692 6659 6693 dtadd(:,:)=0. … … 6662 6696 6663 6697 ! Choices for addtkeoro: 6664 ! ** 0 no TKE tendency from orography 6698 ! ** 0 no TKE tendency from orography 6665 6699 ! ** 1 we include a fraction alphatkeoro of the whole tendency duoro 6666 6700 ! ** 2 we include a fraction alphatkeoro of the gravity wave part of duoro … … 6687 6721 ! Etienne: ici je prends en compte plus de relief que la routine drag_noro_strato 6688 6722 ! car on peut s'attendre a ce que les petites echelles produisent aussi de la TKE 6689 ! Mais attention, cela ne va pas dans le sens de la conservation de l'energie! 6723 ! Mais attention, cela ne va pas dans le sens de la conservation de l'energie! 6690 6724 IF ((zstd(i).GT.1.0) .AND.(nm_oro(i).GT.nm_oro_t)) THEN 6691 6725 itest(i)=1 … … 6695 6729 ENDDO 6696 6730 6697 ELSE 6731 ELSE 6698 6732 6699 6733 igwd=0 … … 6791 6825 6792 6826 IF (ok_cosp) THEN 6793 ! adeclarer 6827 ! adeclarer 6794 6828 IF (CPPKEY_COSP) THEN 6795 6829 IF (itap.eq.1.or.MOD(itap,NINT(freq_cosp/phys_tstep)).EQ.0) THEN … … 6959 6993 ELSE 6960 6994 sh_in(:,:) = qx(:,:,ivap) 6961 IF (nqo >= 3) THEN 6995 IF (nqo >= 3) THEN 6962 6996 ch_in(:,:) = qx(:,:,iliq) + qx(:,:,isol) 6963 ELSE 6997 ELSE 6964 6998 ch_in(:,:) = qx(:,:,iliq) 6965 6999 ENDIF … … 7119 7153 ENDDO 7120 7154 ENDIF 7121 7155 7122 7156 DO i = 1, klon 7123 7157 !--compute ratio of what q+ql should be with conservation to what it is 7124 7158 IF (ok_bs) THEN 7125 corrqql=(qql1(i)+(evap(i)- rain_fall(i)-snow_fall(i)-bs_fall(i))*pdtphys)/qql2(i)7159 corrqql=(qql1(i)+(evap(i)-snowerosion(i)-rain_fall(i)-snow_fall(i)-bs_fall(i))*pdtphys)/qql2(i) 7126 7160 ELSE 7127 7161 corrqql=(qql1(i)+(evap(i)-rain_fall(i)-snow_fall(i))*pdtphys)/qql2(i) … … 7142 7176 7143 7177 !cc prw = eau precipitable 7144 ! prlw = colonne eau liquide 7178 ! prlw = colonne eau liquide 7145 7179 ! prlw = colonne eau solide 7146 7180 ! prbsw = colonne neige soufflee … … 7255 7289 ENDDO 7256 7290 ENDDO 7291 7292 ! in case of advection of TKE 7293 IF (ok_advtke) THEN 7294 DO k=1,klev 7295 DO i=1,klon 7296 d_qx(i,k,itke)=((pbl_tke(i,k,is_ave)+pbl_tke(i,k+1,is_ave))/2. - qx(i,k,itke)) / phys_tstep 7297 ENDDO 7298 ENDDO 7299 ENDIF 7300 ! 7257 7301 7258 7302 ! Lea Raillard qs_ini for cloud phase param. … … 7303 7347 cf_ancien(:,:) = cf_seri(:,:) 7304 7348 rvc_ancien(:,:)= rvc_seri(:,:) 7305 #ifdef ISO 7349 tke_ancien(:,:)= pbl_tke(:,:,is_ave) 7350 7351 7352 #ifdef ISO 7306 7353 xt_ancien(:,:,:)=xt_seri(:,:,:) 7307 7354 xtl_ancien(:,:,:)=xtl_seri(:,:,:) … … 7389 7436 7390 7437 ! Recupere des varibles calcule dans differents modules 7391 ! pour ecriture dans histxxx.nc 7438 ! pour ecriture dans histxxx.nc 7392 7439 7393 7440 ! Get some variables from module fonte_neige_mod … … 7440 7487 CALL phys_output_write(itap, pdtphys, paprs, pphis, & 7441 7488 pplay, lmax_th, aerosol_couple, & 7442 ok_ade, ok_aie, ok_volcan, ivap, iliq, isol, ibs, & 7489 ok_ade, ok_aie, ok_volcan, ivap, iliq, isol, ibs, & 7443 7490 ok_sync, ptconv, read_climoz, clevSTD, & 7444 7491 ptconvth, d_u, d_t, qx, d_qx, zmasse, & … … 7498 7545 alert_first_call = .FALSE. 7499 7546 7500 7547 7501 7548 IF (lafin) THEN 7502 7549 itau_phy = itau_phy + itap … … 7505 7552 ! write(97) u_seri,v_seri,t_seri,q_seri 7506 7553 ! close(97) 7507 7554 7508 7555 IF (is_omp_master) THEN 7509 7556 7510 7557 IF (read_climoz >= 1) THEN 7511 7558 IF (is_mpi_root) CALL nf95_close(ncid_climoz) … … 7513 7560 DEALLOCATE(press_cen_climoz) 7514 7561 ENDIF 7515 7562 7516 7563 ENDIF 7517 7564 … … 7531 7578 7532 7579 WRITE(lunout,*) ' physiq fin, nombre de steps ou cvpas = 1 : ', Ncvpaseq1 7533 7580 7534 7581 ENDIF 7535 7582
Note: See TracChangeset
for help on using the changeset viewer.
