Changeset 3940 for LMDZ6/trunk/libf/phylmdiso/pbl_surface_mod.F90
- Timestamp:
- Jun 15, 2021, 1:18:14 PM (3 years ago)
- File:
-
- 1 edited
-
LMDZ6/trunk/libf/phylmdiso/pbl_surface_mod.F90 (modified) (84 diffs)
Legend:
- Unmodified
- Added
- Removed
-
LMDZ6/trunk/libf/phylmdiso/pbl_surface_mod.F90
r3927 r3940 1 1 ! 2 ! $Id: pbl_surface_mod.F90 3 435 2019-01-22 15:21:59Z fairhead$2 ! $Id: pbl_surface_mod.F90 3906 2021-05-19 10:35:18Z jyg $ 3 3 ! 4 4 MODULE pbl_surface_mod … … 23 23 USE climb_wind_mod, ONLY : climb_wind_down, climb_wind_up 24 24 USE coef_diff_turb_mod, ONLY : coef_diff_turb 25 USE wx_pbl_mod, ONLY : wx_pbl_init, wx_pbl_final, & 26 !! wx_pbl_fuse_no_dts, wx_pbl_split_no_dts, & 27 !! wx_pbl_fuse, wx_pbl_split 28 wx_pbl0_fuse, wx_pbl0_split 25 USE ioipsl_getin_p_mod, ONLY : getin_p 26 USE cdrag_mod 27 USE stdlevvar_mod 28 USE wx_pbl_var_mod, ONLY : wx_pbl_init, wx_pbl_final, & 29 wx_pbl_prelim_0, wx_pbl_prelim_beta 30 USE wx_pbl_mod, ONLY : wx_pbl0_merge, wx_pbl_split, wx_pbl_dts_merge, & 31 wx_pbl_check, wx_pbl_dts_check, wx_evappot 32 use config_ocean_skin_m, only: activate_ocean_skin 29 33 30 34 IMPLICIT NONE … … 33 37 REAL, ALLOCATABLE, DIMENSION(:), PRIVATE, SAVE :: fder ! flux drift 34 38 !$OMP THREADPRIVATE(fder) 35 REAL, ALLOCATABLE, DIMENSION(:,:), PUBLIC, SAVE :: snow ! snow at surface39 REAL, ALLOCATABLE, DIMENSION(:,:), PUBLIC, SAVE :: snow ! snow at surface 36 40 !$OMP THREADPRIVATE(snow) 37 41 REAL, ALLOCATABLE, DIMENSION(:,:), PRIVATE, SAVE :: qsurf ! humidity at surface 38 42 !$OMP THREADPRIVATE(qsurf) 39 REAL, ALLOCATABLE, DIMENSION(:,:,:), SAVE :: ftsoil ! soil temperature43 REAL, ALLOCATABLE, DIMENSION(:,:,:), SAVE :: ftsoil ! soil temperature 40 44 !$OMP THREADPRIVATE(ftsoil) 45 REAL, ALLOCATABLE, DIMENSION(:), SAVE :: ydTs0, ydqs0 46 ! nul forced temperature and humidity differences 47 !$OMP THREADPRIVATE(ydTs0, ydqs0) 41 48 42 49 #ifdef ISO … … 51 58 INTEGER, SAVE :: iflag_pbl_surface_t2m_bug 52 59 !$OMP THREADPRIVATE(iflag_pbl_surface_t2m_bug) 60 INTEGER, SAVE :: iflag_new_t2mq2m 61 !$OMP THREADPRIVATE(iflag_new_t2mq2m) 62 53 63 !FC 54 64 ! integer, save :: iflag_frein … … 78 88 REAL, DIMENSION(klon, nbsrf), INTENT(IN) :: qsurf_rst 79 89 REAL, DIMENSION(klon, nsoilmx, nbsrf), INTENT(IN) :: ftsoil_rst 80 81 90 82 91 ! Local variables … … 102 111 IF (ierr /= 0) CALL abort_physic('pbl_surface_init', 'pb in allocation',1) 103 112 113 ALLOCATE(ydTs0(klon), stat=ierr) 114 IF (ierr /= 0) CALL abort_physic('pbl_surface_init', 'pb in allocation',1) 115 116 ALLOCATE(ydqs0(klon), stat=ierr) 117 IF (ierr /= 0) CALL abort_physic('pbl_surface_init', 'pb in allocation',1) 118 104 119 fder(:) = fder_rst(:) 105 120 snow(:,:) = snow_rst(:,:) 106 121 qsurf(:,:) = qsurf_rst(:,:) 107 122 ftsoil(:,:,:) = ftsoil_rst(:,:,:) 108 123 ydTs0(:) = 0. 124 ydqs0(:) = 0. 109 125 110 126 !**************************************************************************************** … … 152 168 iflag_pbl_surface_t2m_bug=0 153 169 CALL getin_p('iflag_pbl_surface_t2m_bug',iflag_pbl_surface_t2m_bug) 170 WRITE(lunout,*) 'iflag_pbl_surface_t2m_bug=',iflag_pbl_surface_t2m_bug 154 171 !FC 155 172 ! iflag_frein = 0 … … 240 257 debut, lafin, & 241 258 rlon, rlat, rugoro, rmu0, & 242 zsig, lwdown_m, pphi,cldt, &243 rain_f, snow_f, solsw_m, sol lw_m, &259 lwdown_m, cldt, & 260 rain_f, snow_f, solsw_m, solswfdiff_m, sollw_m, & 244 261 gustiness, & 245 262 t, q, u, v, & … … 252 269 ts,SFRWL, alb_dir, alb_dif,ustar, u10m, v10m,wstar, & 253 270 cdragh, cdragm, zu1, zv1, & 271 !jyg< (26/09/2019) 272 beta, & 273 !>jyg 254 274 alb_dir_m, alb_dif_m, zxsens, zxevap, & 255 275 alb3_lic, runoff, snowhgt, qsnow, to_ice, sissnow, & 256 zxtsol, zxfluxlat, zt2m, qsat2m, 276 zxtsol, zxfluxlat, zt2m, qsat2m, zn2mout, & 257 277 d_t, d_q, d_u, d_v, d_t_diss, & 258 278 !!! nrlmd+jyg le 02/05/2011 et le 20/02/2012 … … 275 295 s_therm, s_trmb1, s_trmb2, s_trmb3, & 276 296 zustar,zu10m, zv10m, fder_print, & 277 zxqsurf, rh2m, zxfluxu, zxfluxv, & 297 zxqsurf, delta_qsurf, & 298 rh2m, zxfluxu, zxfluxv, & 278 299 z0m, z0h, agesno, sollw, solsw, & 279 300 d_ts, evap, fluxlat, t2m, & … … 283 304 !jyg< 284 305 !! zxfluxt, zxfluxq, q2m, flux_q, tke, & 285 zxfluxt, zxfluxq, q2m, flux_q, tke_x, 306 zxfluxt, zxfluxq, q2m, flux_q, tke_x, & 286 307 !>jyg 287 308 !!! nrlmd+jyg le 02/05/2011 et le 20/02/2012 … … 338 359 ! z0m, z0h ----input-R- longeur de rugosite (en m) 339 360 ! Martin 340 ! zsig-----input-R- slope341 361 ! cldt-----input-R- total cloud fraction 342 ! pphi-----input-R- geopotentiel de chaque couche (g z) (reference sol)343 362 ! Martin 344 363 ! … … 370 389 USE carbon_cycle_mod, ONLY : carbon_cycle_cpl, carbon_cycle_tr, level_coupling_esm 371 390 USE carbon_cycle_mod, ONLY : co2_send, nbcf_out, fields_out, yfields_out, cfname_out 391 use hbtm_mod, only: hbtm 372 392 USE indice_sol_mod 373 393 USE time_phylmdz_mod, ONLY : day_ini,annee_ref,itau_phy … … 385 405 #endif 386 406 USE ioipsl_getin_p_mod, ONLY : getin_p 407 use phys_state_var_mod, only: ds_ns, dt_ns, delta_sst, delta_sal, zsig, zmea 408 use phys_output_var_mod, only: dter, dser, tkt, tks, taur, sss 409 #ifdef CPP_XIOS 410 USE wxios, ONLY: missing_val 411 #else 412 use netcdf, only: missing_val => nf90_fill_real 413 #endif 414 415 416 387 417 388 418 IMPLICIT NONE … … 412 442 REAL, DIMENSION(klon), INTENT(IN) :: snow_f ! snow fall 413 443 REAL, DIMENSION(klon), INTENT(IN) :: solsw_m ! net shortwave radiation at mean surface 444 REAL, DIMENSION(klon), INTENT(IN) :: solswfdiff_m ! diffuse fraction fordownward shortwave radiation at mean surface 414 445 REAL, DIMENSION(klon), INTENT(IN) :: sollw_m ! net longwave radiation at mean surface 415 446 REAL, DIMENSION(klon,klev), INTENT(IN) :: t ! temperature (K) … … 421 452 REAL, DIMENSION(klon, nbsrf), INTENT(IN) :: pctsrf ! sub-surface fraction 422 453 ! Martin 423 REAL, DIMENSION(klon), INTENT(IN) :: zsig ! slope424 454 REAL, DIMENSION(klon), INTENT(IN) :: lwdown_m ! downward longwave radiation at mean s 425 455 REAL, DIMENSION(klon), INTENT(IN) :: gustiness ! gustiness 426 456 427 457 REAL, DIMENSION(klon), INTENT(IN) :: cldt ! total cloud fraction 428 REAL, DIMENSION(klon,klev), INTENT(IN) :: pphi ! geopotential (m2/s2)429 ! Martin430 458 431 459 #ifdef ISO … … 451 479 ! Input/Output variables 452 480 !**************************************************************************************** 481 !jyg< 482 REAL, DIMENSION(klon, nbsrf), INTENT(INOUT) :: beta ! Aridity factor 483 !>jyg 453 484 REAL, DIMENSION(klon, nbsrf), INTENT(INOUT) :: ts ! temperature at surface (K) 454 485 REAL, DIMENSION(klon, nbsrf), INTENT(INOUT) :: delta_tsurf !surface temperature difference between … … 496 527 REAL, DIMENSION(klon), INTENT(OUT) :: zxfluxlat ! latent flux, mean for each grid point 497 528 REAL, DIMENSION(klon), INTENT(OUT) :: zt2m ! temperature at 2m, mean for each grid point 529 INTEGER, DIMENSION(klon, 6), INTENT(OUT) :: zn2mout ! number of times the 2m temperature is out of the [tsol,temp] 498 530 REAL, DIMENSION(klon), INTENT(OUT) :: qsat2m 499 531 REAL, DIMENSION(klon, klev), INTENT(OUT) :: d_t ! change in temperature … … 560 592 REAL, DIMENSION(klon), INTENT(OUT) :: fder_print ! fder for printing (=fder(i) + dflux_t(i) + dflux_q(i)) 561 593 REAL, DIMENSION(klon), INTENT(OUT) :: zxqsurf ! humidity at surface, mean for each grid point 594 REAL, DIMENSION(klon), INTENT(OUT) :: delta_qsurf! humidity difference at surface, mean for each grid point 562 595 REAL, DIMENSION(klon), INTENT(OUT) :: rh2m ! relative humidity at 2m 563 596 REAL, DIMENSION(klon, klev), INTENT(OUT) :: zxfluxu ! u wind tension, mean for each grid point 564 597 REAL, DIMENSION(klon, klev), INTENT(OUT) :: zxfluxv ! v wind tension, mean for each grid point 565 REAL, DIMENSION(klon, nbsrf+1), INTENT(INOUT) 566 REAL, DIMENSION(klon, nbsrf), INTENT(INOUT) 598 REAL, DIMENSION(klon, nbsrf+1), INTENT(INOUT) :: z0m,z0h ! rugosity length (m) 599 REAL, DIMENSION(klon, nbsrf), INTENT(INOUT) :: agesno ! age of snow at surface 567 600 REAL, DIMENSION(klon, nbsrf), INTENT(OUT) :: solsw ! net shortwave radiation at surface 568 601 REAL, DIMENSION(klon, nbsrf), INTENT(OUT) :: sollw ! net longwave radiation at surface 569 602 REAL, DIMENSION(klon, nbsrf), INTENT(OUT) :: d_ts ! change in temperature at surface 570 REAL, DIMENSION(klon, nbsrf), INTENT(INOUT) 603 REAL, DIMENSION(klon, nbsrf), INTENT(INOUT) :: evap ! evaporation at surface 571 604 REAL, DIMENSION(klon, nbsrf), INTENT(OUT) :: fluxlat ! latent flux 572 605 REAL, DIMENSION(klon, nbsrf), INTENT(OUT) :: t2m ! temperature at 2 meter height … … 605 638 606 639 ! Martin 607 ! sisvat640 ! inlandsis 608 641 REAL, DIMENSION(klon), INTENT(OUT) :: qsnow ! snow water content 609 642 REAL, DIMENSION(klon), INTENT(OUT) :: snowhgt ! snow height … … 632 665 INTEGER :: n 633 666 ! << PC 634 INTEGER :: iflag_split 667 INTEGER :: iflag_split, iflag_split_ref 635 668 INTEGER :: i, k, nsrf 636 669 INTEGER :: knon, j … … 643 676 REAL, DIMENSION(klon) :: r_co2_ppm ! taux CO2 atmosphere 644 677 REAL, DIMENSION(klon) :: yts, yz0m, yz0h, ypct 678 REAL, DIMENSION(klon) :: yz0h_old 645 679 !albedo SB >>> 646 680 REAL, DIMENSION(klon) :: yalb,yalb_vis 647 681 !albedo SB <<< 648 682 REAL, DIMENSION(klon) :: yt1, yq1, yu1, yv1 683 REAL, DIMENSION(klon) :: yqa 649 684 REAL, DIMENSION(klon) :: ysnow, yqsurf, yagesno, yqsol 650 685 REAL, DIMENSION(klon) :: yrain_f, ysnow_f … … 670 705 REAL, DIMENSION(klon) :: y_flux_u1, y_flux_v1 671 706 REAL, DIMENSION(klon) :: yt2m, yq2m, yu10m 707 INTEGER, DIMENSION(klon, nbsrf, 6) :: yn2mout, yn2mout_x, yn2mout_w 708 INTEGER, DIMENSION(klon, nbsrf, 6) :: n2mout, n2mout_x, n2mout_w 672 709 REAL, DIMENSION(klon) :: yustar 673 710 REAL, DIMENSION(klon) :: ywstar … … 690 727 REAL, DIMENSION(klon) :: yz0h_oupas 691 728 REAL, DIMENSION(klon) :: yfluxsens 729 REAL, DIMENSION(klon) :: AcoefH_0, AcoefQ_0, BcoefH_0, BcoefQ_0 692 730 REAL, DIMENSION(klon) :: AcoefH, AcoefQ, BcoefH, BcoefQ 693 731 #ifdef ISO … … 696 734 REAL, DIMENSION(klon) :: AcoefU, AcoefV, BcoefU, BcoefV 697 735 REAL, DIMENSION(klon) :: ypsref 698 REAL, DIMENSION(klon) :: yevap, y tsurf_new, yalb3_new736 REAL, DIMENSION(klon) :: yevap, yevap_pot, ytsurf_new, yalb3_new 699 737 !albedo SB >>> 700 738 REAL, DIMENSION(klon,nsw) :: yalb_dir_new, yalb_dif_new … … 708 746 REAL, DIMENSION(klon,klev) :: y_flux_u, y_flux_v 709 747 REAL, DIMENSION(klon,klev) :: ycoefh, ycoefm,ycoefq 710 REAL, DIMENSION(klon) :: ycdragh, ycdrag m748 REAL, DIMENSION(klon) :: ycdragh, ycdragq, ycdragm 711 749 REAL, DIMENSION(klon,klev) :: yu, yv 712 750 REAL, DIMENSION(klon,klev) :: yt, yq … … 746 784 REAL, DIMENSION(klon,klev) :: ycoefh_x, ycoefm_x, ycoefh_w, ycoefm_w 747 785 REAL, DIMENSION(klon,klev) :: ycoefq_x, ycoefq_w 748 REAL, DIMENSION(klon) :: ycdragh_x, ycdragm_x, ycdragh_w, ycdragm_w 786 REAL, DIMENSION(klon) :: ycdragh_x, ycdragh_w, ycdragq_x, ycdragq_w 787 REAL, DIMENSION(klon) :: ycdragm_x, ycdragm_w 749 788 REAL, DIMENSION(klon) :: AcoefH_x, AcoefQ_x, BcoefH_x, BcoefQ_x 750 789 REAL, DIMENSION(klon) :: AcoefH_w, AcoefQ_w, BcoefH_w, BcoefQ_w … … 767 806 REAL :: zx_qs_surf, zcor_surf, zdelta_surf 768 807 REAL, DIMENSION(klon) :: ytsurf_th, yqsatsurf 808 !jyg< 769 809 REAL, DIMENSION(klon) :: ybeta 810 REAL, DIMENSION(klon) :: ybeta_prev 811 !>jyg 770 812 REAL, DIMENSION(klon, klev) :: d_u_x 771 813 REAL, DIMENSION(klon, klev) :: d_u_w … … 915 957 !!! nrlmd le 13/06/2011 916 958 REAL, DIMENSION(klon) :: y_delta_flux_t1, y_delta_flux_q1, y_delta_flux_u1, y_delta_flux_v1 917 REAL, DIMENSION(klon) :: y_delta_tsurf,delta_coef,tau_eq 959 REAL, DIMENSION(klon) :: y_delta_tsurf, y_delta_tsurf_new 960 REAL, DIMENSION(klon) :: delta_coef, tau_eq 961 REAL, DIMENSION(klon) :: HTphiT_b, dd_HTphiT, HTphiQ_b, dd_HTphiQ, HTRn_b, dd_HTRn 962 REAL, DIMENSION(klon) :: phiT0_b, dphiT0, phiQ0_b, dphiQ0, Rn0_b, dRn0 963 REAL, DIMENSION(klon) :: y_delta_qsurf 964 REAL, DIMENSION(klon) :: y_delta_qsats 965 REAL, DIMENSION(klon) :: yg_T, yg_Q 966 REAL, DIMENSION(klon) :: yGamma_dTs_phiT, yGamma_dQs_phiQ 967 REAL, DIMENSION(klon) :: ydTs_ins, ydqs_ins 968 ! 918 969 REAL, PARAMETER :: facteur=2./sqrt(3.14) 919 970 REAL, PARAMETER :: inertia=2000. … … 928 979 REAL, DIMENSION(klon) :: Kech_m 929 980 REAL, DIMENSION(klon) :: Kech_m_x, Kech_m_w 930 REAL, DIMENSION(klon) :: yts_x,yts_w 981 REAL, DIMENSION(klon) :: yts_x, yts_w 982 REAL, DIMENSION(klon) :: yqsatsrf0_x, yqsatsrf0_w 983 REAL, DIMENSION(klon) :: yqsurf_x, yqsurf_w 931 984 !jyg< 932 985 !! REAL, DIMENSION(klon) :: Kech_Hp, Kech_H_xp, Kech_H_wp … … 935 988 !! REAL, DIMENSION(klon) :: Kech_Vp, Kech_V_xp, Kech_V_wp 936 989 !>jyg 937 !jyg< 938 REAL , DIMENSION(klon) :: ah, bh ! coefficients of the delta_Tsurf equation939 !>jyg 990 991 REAL :: fact_cdrag 992 REAL :: z1lay 940 993 941 994 REAL :: vent … … 971 1024 REAL, DIMENSION(klon) :: ytoice 972 1025 REAL, DIMENSION(klon) :: ysnowhgt, yqsnow, ysissnow, yrunoff 1026 REAL, DIMENSION(klon) :: yzmea 973 1027 REAL, DIMENSION(klon) :: yzsig 974 REAL, DIMENSION(klon,klev) :: ypphi975 1028 REAL, DIMENSION(klon) :: ycldt 976 1029 REAL, DIMENSION(klon) :: yrmu0 977 1030 ! Martin 1031 1032 REAL, DIMENSION(klon):: ydelta_sst, ydelta_sal, yds_ns, ydt_ns, ydter, ydser, & 1033 ytkt, ytks, ytaur, ysss 1034 ! compression of delta_sst, delta_sal, ds_ns, dt_ns, dter, dser, tkt, tks, 1035 ! taur, sss on ocean points 1036 978 1037 #ifdef ISO 979 1038 REAL, DIMENSION(klon) :: h1 … … 991 1050 ! 992 1051 !!jyg iflag_split = mod(iflag_pbl_split,2) 993 iflag_split = mod(iflag_pbl_split,10) 1052 !!jyg iflag_split = mod(iflag_pbl_split,10) 1053 iflag_split_ref = mod(iflag_pbl_split,10) 994 1054 995 1055 #ifdef ISO … … 1038 1098 1039 1099 IF (first_call) THEN 1100 1101 iflag_new_t2mq2m=1 1102 CALL getin_p('iflag_new_t2mq2m',iflag_new_t2mq2m) 1103 WRITE(lunout,*) 'pbl_iflag_new_t2mq2m=',iflag_new_t2mq2m 1104 1040 1105 print*,'PBL SURFACE AVEC GUSTINESS' 1041 1106 first_call=.FALSE. 1042 1107 1043 1108 ! Initialize ok_flux_surf (for 1D model) 1044 if (klon_glo>1) ok_flux_surf=.FALSE. 1109 IF (klon_glo>1) ok_flux_surf=.FALSE. 1110 IF (klon_glo>1) ok_forc_tsurf=.FALSE. 1045 1111 1046 1112 ! intialize beta_land … … 1113 1179 zxfluxlat(:)=0. 1114 1180 zt2m(:)=0. ; zq2m(:)=0. ; qsat2m(:)=0. ; rh2m(:)=0. 1181 zn2mout(:,:)=0 ; 1115 1182 d_t(:,:)=0. ; d_t_diss(:,:)=0. ; d_q(:,:)=0. ; d_u(:,:)=0. ; d_v(:,:)=0. 1116 1183 zcoefh(:,:,:)=0. ; zcoefm(:,:,:)=0. … … 1128 1195 fder_print(:)=0. 1129 1196 zxqsurf(:)=0. 1197 delta_qsurf(:) = 0. 1130 1198 zxfluxu(:,:)=0. ; zxfluxv(:,:)=0. 1131 1199 solsw(:,:)=0. ; sollw(:,:)=0. … … 1164 1232 !! tke(:,:,is_ave)=0. 1165 1233 tke_x(:,:,is_ave)=0. 1234 1166 1235 wake_dltke(:,:,is_ave)=0. 1167 1236 !>jyg … … 1184 1253 yqsurf = 0.0 ; yalb = 0.0 ; yalb_vis = 0.0 1185 1254 !albedo SB <<< 1186 yrain_f = 0.0 ; ysnow_f = 0.0 ; yfder = 0.0 ; ysolsw = 0.0 1255 yrain_f = 0.0 ; ysnow_f = 0.0 ; yfder = 0.0 ; ysolsw = 0.0 1187 1256 ysollw = 0.0 ; yz0m = 0.0 ; yz0h = 0.0 ; yu1 = 0.0 1188 1257 yv1 = 0.0 ; ypaprs = 0.0 ; ypplay = 0.0 … … 1195 1264 !! d_t_diss= 0.0 ;d_u = 0.0 ; d_v = 0.0 1196 1265 yqsol = 0.0 1197 ytherm = 0.0 ; ytke=0. 1266 1267 ytke=0. 1198 1268 !FC 1199 1269 y_treedrg=0. … … 1202 1272 ysnowhgt = 0.0; yqsnow = 0.0 ; yrunoff = 0.0 ; ytoice =0.0 1203 1273 yalb3_new = 0.0 ; ysissnow = 0.0 1204 y pphi = 0.0 ; ycldt = 0.0 ; yrmu0 = 0.01274 ycldt = 0.0 ; yrmu0 = 0.0 1205 1275 ! Martin 1206 1276 … … 1218 1288 y_delta_flux_t1=0. 1219 1289 ydtsurf_th=0. 1220 yts_x=0. ; yts_w=0. 1221 y_delta_tsurf=0. 1290 yts_x(:)=0. ; yts_w(:)=0. 1291 y_delta_tsurf(:)=0. ; y_delta_qsurf(:)=0. 1292 yqsurf_x(:)=0. ; yqsurf_w(:)=0. 1293 yg_T(:) = 0. ; yg_Q(:) = 0. 1294 yGamma_dTs_phiT(:) = 0. ; yGamma_dQs_phiQ(:) = 0. 1295 ydTs_ins(:) = 0. ; ydqs_ins(:) = 0. 1296 1222 1297 !!! 1223 1298 ytsoil = 999999. … … 1410 1485 DO i = 1, klon 1411 1486 sollw(i,nsrf) = sollw_m(i) + 4.0*RSIGMA*ztsol(i)**3 * (ztsol(i)-ts(i,nsrf)) 1412 1413 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 1414 ! ! Martin 1415 ! Apparently introduced for sisvat but not used 1416 ! sollwd(i,nsrf)= sollwd_m(i) 1417 ! ! Martin 1418 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 1419 1487 !--OB this line is not satisfactory because alb is the direct albedo not total albedo 1420 1488 solsw(i,nsrf) = solsw_m(i) * (1.-alb(i,nsrf)) / (1.-alb_m(i)) 1421 1489 ENDDO … … 1440 1508 !>al1 1441 1509 1510 !--OB add diffuse fraction of SW down 1511 DO n=1,nbcf_out 1512 IF (cfname_out(n) == "swdownfdiff" ) fields_out(:,n) = solswfdiff_m(:) 1513 ENDDO 1442 1514 ! >> PC 1443 1515 IF (carbon_cycle_cpl .AND. carbon_cycle_tr .AND. nbcf_out.GT.0 ) THEN … … 1461 1533 ! 1462 1534 !**************************************************************************************** 1463 1464 loop_nbsrf: DO nsrf = 1, nbsrf 1535 !<<<<<<<<<<<<< 1536 loop_nbsrf: DO nsrf = 1, nbsrf !<<<<<<<<<<<<< 1537 !<<<<<<<<<<<<< 1465 1538 IF (prt_level >=10) print *,' Loop nsrf ',nsrf 1539 ! 1540 IF (iflag_split_ref == 3) THEN 1541 IF (nsrf == is_oce) THEN 1542 iflag_split = 1 1543 ELSE 1544 iflag_split=0 1545 ENDIF !! (nsrf == is_oce) 1546 ELSE 1547 iflag_split = iflag_split_ref 1548 ENDIF !! (iflag_split_ref == 3) 1466 1549 1467 1550 ! Search for index(ni) and size(knon) of domaine to treat … … 1499 1582 !**************************************************************************************** 1500 1583 1584 ! 1585 !jyg< (20190926) 1586 ! Provisional : set ybeta to standard values 1587 IF (nsrf .NE. is_ter) THEN 1588 ybeta(:) = 1. 1589 ELSE 1590 IF (iflag_split .EQ. 0) THEN 1591 ybeta(:) = 1. 1592 ELSE 1593 DO j = 1, knon 1594 i = ni(j) 1595 ybeta(j) = beta(i,nsrf) 1596 ENDDO 1597 ENDIF ! (iflag_split .LE.1) 1598 ENDIF ! (nsrf .NE. is_ter) 1599 !>jyg 1600 ! 1501 1601 DO j = 1, knon 1502 1602 i = ni(j) … … 1531 1631 ywindsp(j) = windsp(i,nsrf) 1532 1632 !>jyg 1533 ! Martin 1633 ! Martin and Etienne 1634 yzmea(j) = zmea(i) 1534 1635 yzsig(j) = zsig(i) 1535 1636 ycldt(j) = cldt(i) … … 1656 1757 ytke_w(j,k) = tke_x(i,k,nsrf)+wake_dltke(i,k,nsrf) 1657 1758 ywake_dltke(j,k) = wake_dltke(i,k,nsrf) 1759 1658 1760 !>jyg 1659 1761 ENDDO … … 1695 1797 ENDDO 1696 1798 ENDIF 1799 1800 if (nsrf == is_oce .and. activate_ocean_skin >= 1) then 1801 if (activate_ocean_skin == 2 .and. type_ocean == "couple") then 1802 ydelta_sal(:knon) = delta_sal(ni(:knon)) 1803 ydelta_sst(:knon) = delta_sst(ni(:knon)) 1804 end if 1805 1806 yds_ns(:knon) = ds_ns(ni(:knon)) 1807 ydt_ns(:knon) = dt_ns(ni(:knon)) 1808 end if 1697 1809 1698 1810 !**************************************************************************************** … … 1735 1847 ENDDO 1736 1848 ENDIF 1849 1737 1850 IF (prt_level >=10) print *,'clcdrag -> ycdragh ', ycdragh 1738 1851 ELSE !(iflag_split .eq.0) … … 1749 1862 speed_x(i) = SQRT(yu_x(i,1)**2+yv_x(i,1)**2) 1750 1863 ENDDO 1751 CALL cdrag(knon, nsrf, & 1864 1865 1866 CALL cdrag(knon, nsrf, & 1752 1867 speed_x, yt_x(:,1), yq_x(:,1), zgeo1_x, ypaprs(:,1),& 1753 yts_x, yqsurf , yz0m, yz0h, &1868 yts_x, yqsurf_x, yz0m, yz0h, & 1754 1869 ycdragm_x, ycdragh_x, zri1_x, pref_x ) 1755 1870 … … 1778 1893 CALL cdrag(knon, nsrf, & 1779 1894 speed_w, yt_w(:,1), yq_w(:,1), zgeo1_w, ypaprs(:,1),& 1780 yts_w, yqsurf , yz0m, yz0h, &1895 yts_w, yqsurf_w, yz0m, yz0h, & 1781 1896 ycdragm_w, ycdragh_w, zri1_w, pref_w ) 1782 1897 ! … … 1818 1933 print *,' args coef_diff_turb: ycdragh ', ycdragh 1819 1934 print *,' args coef_diff_turb: ytke ', ytke 1935 1820 1936 ENDIF 1821 1937 CALL coef_diff_turb(dtime, nsrf, knon, ni, & … … 1847 1963 print *,' args coef_diff_turb: ycdragh_x ', ycdragh_x 1848 1964 print *,' args coef_diff_turb: ytke_x ', ytke_x 1965 1849 1966 ENDIF 1850 1967 CALL coef_diff_turb(dtime, nsrf, knon, ni, & 1851 ypaprs, ypplay, yu_x, yv_x, yq_x, yt_x, yts_x, yqsurf , ycdragm_x, &1968 ypaprs, ypplay, yu_x, yv_x, yq_x, yt_x, yts_x, yqsurf_x, ycdragm_x, & 1852 1969 ycoefm_x, ycoefh_x, ytke_x,y_treedrg) 1853 1970 ! ycoefm_x, ycoefh_x, ytke_x) … … 1877 1994 ENDIF 1878 1995 CALL coef_diff_turb(dtime, nsrf, knon, ni, & 1879 ypaprs, ypplay, yu_w, yv_w, yq_w, yt_w, yts_w, yqsurf , ycdragm_w, &1996 ypaprs, ypplay, yu_w, yv_w, yq_w, yt_w, yts_w, yqsurf_w, ycdragm_w, & 1880 1997 ycoefm_w, ycoefh_w, ytke_w,y_treedrg) 1881 1998 ! ycoefm_w, ycoefh_w, ytke_w) … … 2029 2146 yxt1(:,:) = yxt(:,:,1) 2030 2147 #endif 2031 !! ELSE IF (iflag_split .eq. 1) THEN 2032 !!! 2033 !jyg< 2034 !! CALL wx_pbl_fuse_no_dts(knon, dtime, ypplay, ywake_s, & 2035 !! yt_x, yt_w, yq_x, yq_w, & 2036 !! yu_x, yu_w, yv_x, yv_w, & 2037 !! ycdragh_x, ycdragh_w, ycdragm_x, ycdragm_w, & 2038 !! AcoefH_x, AcoefH_w, AcoefQ_x, AcoefQ_w, & 2039 !! AcoefU_x, AcoefU_w, AcoefV_x, AcoefV_w, & 2040 !! BcoefH_x, BcoefH_w, BcoefQ_x, BcoefQ_w, & 2041 !! BcoefU_x, BcoefU_w, BcoefV_x, BcoefV_w, & 2042 !! AcoefH, AcoefQ, AcoefU, AcoefV, & 2043 !! BcoefH, BcoefQ, BcoefU, BcoefV, & 2044 !! ycdragh, ycdragm, & 2045 !! yt1, yq1, yu1, yv1 & 2046 !! ) 2148 2047 2149 ELSE IF (iflag_split .ge. 1) THEN 2048 CALL wx_pbl0_fuse(knon, dtime, ypplay, ywake_s, & 2150 #ifdef ISO 2151 call abort_gcm('pbl_surface_mod 2149','isos pas encore dans iflag_split=1',1) 2152 #endif 2153 2154 ! 2155 ! Cdragq computation 2156 ! ------------------ 2157 !****************************************************************************** 2158 ! Cdragq computed from cdrag 2159 ! The difference comes only from a factor (f_z0qh_oce) on z0, so that 2160 ! it can be computed inside wx_pbl0_merge 2161 ! More complicated appraches may require the propagation through 2162 ! pbl_surface of an independant cdragq variable. 2163 !****************************************************************************** 2164 ! 2165 IF ( f_z0qh_oce .ne. 1. .and. nsrf .eq.is_oce) THEN 2166 ! Si on suit les formulations par exemple de Tessel, on 2167 ! a z0h=0.4*nu/u*, z0q=0.62*nu/u*, d'ou f_z0qh_oce=0.62/0.4=1.55 2168 !! ycdragq_x(1:knon)=ycdragh_x(1:knon)* & 2169 !! log(z1lay(1:knon)/yz0h(1:knon))/log(z1lay(1:knon)/(f_z0qh_oce*yz0h(1:knon))) 2170 !! ycdragq_w(1:knon)=ycdragh_w(1:knon)* & 2171 !! log(z1lay(1:knon)/yz0h(1:knon))/log(z1lay(1:knon)/(f_z0qh_oce*yz0h(1:knon))) 2172 ! 2173 DO j = 1,knon 2174 z1lay = zgeo1(j)/RG 2175 fact_cdrag = log(z1lay/yz0h(j))/log(z1lay/(f_z0qh_oce*yz0h(j))) 2176 ycdragq_x(j)=ycdragh_x(j)*fact_cdrag 2177 ycdragq_w(j)=ycdragh_w(j)*fact_cdrag 2178 !! Print *,'YYYYpbl0: fact_cdrag ', fact_cdrag 2179 ENDDO ! j = 1,knon 2180 ! 2181 !! Print *,'YYYYpbl0: z1lay, yz0h, f_z0qh_oce, ycdragh_w, ycdragq_w ', & 2182 !! z1lay, yz0h(1:knon), f_z0qh_oce, ycdragh_w(1:knon), ycdragq_w(1:knon) 2183 ELSE 2184 ycdragq_x(1:knon)=ycdragh_x(1:knon) 2185 ycdragq_w(1:knon)=ycdragh_w(1:knon) 2186 ENDIF ! ( f_z0qh_oce .ne. 1. .and. nsrf .eq.is_oce) 2187 ! 2188 CALL wx_pbl_prelim_0(knon, nsrf, dtime, ypplay, ypaprs, ywake_s, & 2189 yts, y_delta_tsurf, ygustiness, & 2049 2190 yt_x, yt_w, yq_x, yq_w, & 2050 2191 yu_x, yu_w, yv_x, yv_w, & 2051 ycdragh_x, ycdragh_w, ycdragm_x, ycdragm_w, & 2192 ycdragh_x, ycdragh_w, ycdragq_x, ycdragq_w, & 2193 ycdragm_x, ycdragm_w, & 2194 AcoefH_x, AcoefH_w, AcoefQ_x, AcoefQ_w, & 2195 AcoefU_x, AcoefU_w, AcoefV_x, AcoefV_w, & 2196 BcoefH_x, BcoefH_w, BcoefQ_x, BcoefQ_w, & 2197 BcoefU_x, BcoefU_w, BcoefV_x, BcoefV_w & 2198 ) 2199 CALL wx_pbl_prelim_beta(knon, dtime, ywake_s, ybeta, & 2200 BcoefQ_x, BcoefQ_w & 2201 ) 2202 CALL wx_pbl0_merge(knon, ypplay, ypaprs, & 2203 ywake_s, ydTs0, ydqs0, & 2204 yt_x, yt_w, yq_x, yq_w, & 2205 yu_x, yu_w, yv_x, yv_w, & 2206 ycdragh_x, ycdragh_w, ycdragq_x, ycdragq_w, & 2207 ycdragm_x, ycdragm_w, & 2052 2208 AcoefH_x, AcoefH_w, AcoefQ_x, AcoefQ_w, & 2053 2209 AcoefU_x, AcoefU_w, AcoefV_x, AcoefV_w, & 2054 2210 BcoefH_x, BcoefH_w, BcoefQ_x, BcoefQ_w, & 2055 2211 BcoefU_x, BcoefU_w, BcoefV_x, BcoefV_w, & 2056 AcoefH , AcoefQ, AcoefU, AcoefV, &2057 BcoefH , BcoefQ, BcoefU, BcoefV, &2058 ycdragh, ycdrag m, &2212 AcoefH_0, AcoefQ_0, AcoefU, AcoefV, & 2213 BcoefH_0, BcoefQ_0, BcoefU, BcoefV, & 2214 ycdragh, ycdragq, ycdragm, & 2059 2215 yt1, yq1, yu1, yv1 & 2060 #ifdef ISO2061 ,yxt_x,yxt_w,yxt1 &2062 #endif2063 2216 ) 2064 !! ELSE IF (iflag_split .ge.2) THEN 2065 !!! Provisoire 2066 !! ah(:) = 0. 2067 !! bh(:) = 0. 2068 !! IF (nsrf == is_oce) THEN 2069 !! ybeta(:) = 1. 2070 !! ELSE 2071 !! ybeta(:) = beta_land 2072 !! ENDIF 2073 !! ycdragh(:) = ywake_s(:)*ycdragh_w(:) + (1.-ywake_s(:))*ycdragh_x(:) 2074 !! CALL wx_dts(knon, nsrf, ywake_cstar, ywake_s, ywake_dens, & 2075 !! yts, ypplay(:,1), ybeta, ycdragh , ypaprs(:,1), & 2076 !! yq(:,1), yt(:,1), yu(:,1), yv(:,1), ygustiness, & 2077 !! ah, bh & 2078 !! ) 2079 !!! 2080 !! CALL wx_pbl_fuse(knon, dtime, ypplay, ywake_s, & 2081 !! yt_x, yt_w, yq_x, yq_w, & 2082 !! yu_x, yu_w, yv_x, yv_w, & 2083 !! ycdragh_x, ycdragh_w, ycdragm_x, ycdragm_w, & 2084 !! AcoefH_x, AcoefH_w, AcoefQ_x, AcoefQ_w, & 2085 !! AcoefU_x, AcoefU_w, AcoefV_x, AcoefV_w, & 2086 !! BcoefH_x, BcoefH_w, BcoefQ_x, BcoefQ_w, & 2087 !! BcoefU_x, BcoefU_w, BcoefV_x, BcoefV_w, & 2088 !! ah, bh, & 2089 !! AcoefH, AcoefQ, AcoefU, AcoefV, & 2090 !! BcoefH, BcoefQ, BcoefU, BcoefV, & 2091 !! ycdragh, ycdragm, & 2092 !! yt1, yq1, yu1, yv1 & 2093 !! ) 2094 !>jyg 2095 !!! 2096 ENDIF ! (iflag_split .eq.0) 2217 IF (iflag_split .eq. 2 .AND. nsrf .ne. is_oce) THEN 2218 CALL wx_pbl_dts_merge(knon, dtime, ypplay, ypaprs, & 2219 ywake_s, ybeta, ywake_cstar, ywake_dens, & 2220 AcoefH_x, AcoefH_w, & 2221 BcoefH_x, BcoefH_w, & 2222 AcoefH_0, AcoefQ_0, BcoefH_0, BcoefQ_0, & 2223 AcoefH, AcoefQ, BcoefH, BcoefQ, & 2224 HTphiT_b, dd_HTphiT, HTphiQ_b, dd_HTphiQ, HTRn_b, dd_HTRn, & 2225 phiT0_b, dphiT0, phiQ0_b, dphiQ0, Rn0_b, dRn0, & 2226 yg_T, yg_Q, & 2227 yGamma_dTs_phiT, yGamma_dQs_phiQ, & 2228 ydTs_ins, ydqs_ins & 2229 ) 2230 ELSE ! 2231 AcoefH(:) = AcoefH_0(:) 2232 AcoefQ(:) = AcoefQ_0(:) 2233 BcoefH(:) = BcoefH_0(:) 2234 BcoefQ(:) = BcoefQ_0(:) 2235 yg_T(:) = 0. 2236 yg_Q(:) = 0. 2237 yGamma_dTs_phiT(:) = 0. 2238 yGamma_dQs_phiQ(:) = 0. 2239 ydTs_ins(:) = 0. 2240 ydqs_ins(:) = 0. 2241 ENDIF ! (iflag_split .eq. 2) 2242 ENDIF ! (iflag_split .eq.0) 2097 2243 !!! 2098 2244 IF (prt_level >=10) THEN 2099 PRINT *,'pbl_surface (fuse->): yt(1,:) ',yt(1,:) 2100 PRINT *,'pbl_surface (fuse->): yq(1,:) ',yq(1,:) 2101 PRINT *,'pbl_surface (fuse->): yu(1,:) ',yu(1,:) 2102 PRINT *,'pbl_surface (fuse->): yv(1,:) ',yv(1,:) 2103 PRINT *,'pbl_surface (fuse->): AcoefH(1) ',AcoefH(1) 2104 PRINT *,'pbl_surface (fuse->): BcoefH(1) ',BcoefH(1) 2245 PRINT *,'pbl_surface (merge->): yt(1,:) ',yt(1,:) 2246 PRINT *,'pbl_surface (merge->): yq(1,:) ',yq(1,:) 2247 PRINT *,'pbl_surface (merge->): yu(1,:) ',yu(1,:) 2248 PRINT *,'pbl_surface (merge->): yv(1,:) ',yv(1,:) 2249 PRINT *,'pbl_surface (merge->): AcoefH(1), AcoefQ(1), AcoefU(1), AcoefV(1) ', & 2250 AcoefH(1), AcoefQ(1), AcoefU(1), AcoefV(1) 2251 PRINT *,'pbl_surface (merge->): BcoefH(1), BcoefQ(1), BcoefU(1), BcoefV(1) ', & 2252 BcoefH(1), BcoefQ(1), BcoefU(1), BcoefV(1) 2253 2105 2254 ENDIF 2106 2255 2256 ! Save initial value of z0h for use in evappot (z0h wiil be computed again in the surface models) 2257 yz0h_old(1:knon) = yz0h(1:knon) 2258 ! 2107 2259 !**************************************************************************************** 2108 2260 ! … … 2119 2271 2120 2272 ! Calculate the temperature et relative humidity at 2m and the wind at 10m 2273 IF (iflag_new_t2mq2m==1) THEN 2274 CALL stdlevvarn(klon, knon, is_ter, zxli, & 2275 yu(:,1), yv(:,1), yt(:,1), yq(:,1), zgeo1, & 2276 yts, yqsurf, yz0m, yz0h, ypaprs(:,1), ypplay(:,1), & 2277 yt2m, yq2m, yt10m, yq10m, yu10m, yustar, & 2278 yn2mout(:, nsrf, :)) 2279 ELSE 2121 2280 CALL stdlevvar(klon, knon, is_ter, zxli, & 2122 2281 yu(:,1), yv(:,1), yt(:,1), yq(:,1), zgeo1, & 2123 2282 yts, yqsurf, yz0m, yz0h, ypaprs(:,1), ypplay(:,1), & 2124 2283 yt2m, yq2m, yt10m, yq10m, yu10m, yustar) 2284 ENDIF 2125 2285 2126 2286 ENDIF … … 2211 2371 CALL surf_landice(itap, dtime, knon, ni, & 2212 2372 rlon, rlat, debut, lafin, & 2213 yrmu0, ylwdown, yalb, ypphi(:,1), &2373 yrmu0, ylwdown, yalb, zgeo1, & 2214 2374 ysolsw, ysollw, yts, ypplay(:,1), & 2215 2375 !!jyg ycdragh, ycdragm, yrain_f, ysnow_f, yt(:,1), yq(:,1),& … … 2221 2381 ytsoil, yz0m, yz0h, SFRWL, yalb_dir_new, yalb_dif_new, yevap,yfluxsens,yfluxlat, & 2222 2382 ytsurf_new, y_dflux_t, y_dflux_q, & 2223 yz sig, ycldt, &2383 yzmea, yzsig, ycldt, & 2224 2384 ysnowhgt, yqsnow, ytoice, ysissnow, & 2225 2385 yalb3_new, yrunoff, & … … 2284 2444 yz0m, yz0h, SFRWL,yalb_dir_new, yalb_dif_new, yevap, yfluxsens,yfluxlat,& 2285 2445 ytsurf_new, y_dflux_t, y_dflux_q, slab_wfbils, & 2286 y_flux_u1, y_flux_v1 & 2446 y_flux_u1, y_flux_v1, ydelta_sst(:knon), ydelta_sal(:knon), & 2447 yds_ns(:knon), ydt_ns(:knon), ydter(:knon), ydser(:knon), & 2448 ytkt(:knon), ytks(:knon), ytaur(:knon), ysss & 2287 2449 #ifdef ISO 2288 2450 & ,yxtrain_f, yxtsnow_f,yxt1,Roce, & … … 2392 2554 ! 2393 2555 !**************************************************************************************** 2394 2395 !!! 2396 !!! jyg le 10/04/2013 2556 !! 2557 !!! 2558 !!! jyg le 10/04/2013 et EV 10/2020 2559 2560 IF (ok_forc_tsurf) THEN 2561 DO j=1,knon 2562 ytsurf_new(j)=tg 2563 y_d_ts(j) = ytsurf_new(j) - yts(j) 2564 ENDDO 2565 ENDIF ! ok_forc_tsurf 2566 2397 2567 !!! 2398 2568 IF (ok_flux_surf) THEN … … 2429 2599 #endif 2430 2600 ENDDO 2431 ENDIF 2432 2433 IF (prt_level >=10) THEN 2434 DO j=1,knon 2435 print*,'y_flux_t1,yfluxlat,wakes' & 2436 & , y_flux_t1(j), yfluxlat(j), ywake_s(j) 2437 print*,'beta,ytsurf_new', ybeta(j), ytsurf_new(j) 2438 print*,'inertia,facteur,cstar', inertia, facteur,wake_cstar(j) 2439 ENDDO 2440 ENDIF 2441 2442 !!! jyg le 07/02/2012 puis le 10/04/2013 2443 !! IF (iflag_split .eq.1) THEN 2444 !!!!! 2445 !!!jyg< 2446 !! CALL wx_pbl_split_no_dts(knon, ywake_s, & 2447 !! AcoefH_x, AcoefH_w, & 2448 !! AcoefQ_x, AcoefQ_w, & 2449 !! AcoefU_x, AcoefU_w, & 2450 !! AcoefV_x, AcoefV_w, & 2451 !! y_flux_t1, y_flux_q1, y_flux_u1, y_flux_v1, & 2452 !! y_flux_t1_x, y_flux_t1_w, & 2453 !! y_flux_q1_x, y_flux_q1_w, & 2454 !! y_flux_u1_x, y_flux_u1_w, & 2455 !! y_flux_v1_x, y_flux_v1_w, & 2456 !! yfluxlat_x, yfluxlat_w & 2457 !! ) 2458 !! ELSE IF (iflag_split .ge. 2) THEN 2601 ENDIF ! (ok_flux_surf) 2602 ! 2603 ! ------------------------------------------------------------------------------ 2604 ! 12a) Splitting 2605 ! ------------------------------------------------------------------------------ 2606 2459 2607 IF (iflag_split .GE. 1) THEN 2460 CALL wx_pbl0_split(knon, dtime, ywake_s, & 2608 #ifdef ISO 2609 call abort_gcm('pbl_surface_mod 2607','isos pas encore dans iflag_split=1',1) 2610 #endif 2611 ! 2612 IF (nsrf .ne. is_oce) THEN 2613 ! 2614 ! Compute potential evaporation and aridity factor (jyg, 20200328) 2615 ybeta_prev(:) = ybeta(:) 2616 DO j = 1, knon 2617 yqa(j) = AcoefQ(j) - BcoefQ(j)*yevap(j)*dtime 2618 ENDDO 2619 ! 2620 CALL wx_evappot(knon, yqa, yTsurf_new, yevap_pot) 2621 ! 2622 ybeta(1:knon) = min(yevap(1:knon)/yevap_pot(1:knon), 1.) 2623 2624 IF (prt_level >=10) THEN 2625 DO j=1,knon 2626 print*,'y_flux_t1,yfluxlat,wakes' & 2627 & , y_flux_t1(j), yfluxlat(j), ywake_s(j) 2628 print*,'beta_prev, beta, ytsurf_new', ybeta_prev(j), ybeta(j), ytsurf_new(j) 2629 print*,'inertia,facteur,cstar', inertia, facteur,wake_cstar(j) 2630 ENDDO 2631 ENDIF ! (prt_level >=10) 2632 ! 2633 ! Second call to wx_pbl0_merge and wx_pbl_dts_merge in order to take into account 2634 ! the update of the aridity coeficient beta. 2635 ! 2636 CALL wx_pbl_prelim_beta(knon, dtime, ywake_s, ybeta, & 2637 BcoefQ_x, BcoefQ_w & 2638 ) 2639 CALL wx_pbl0_merge(knon, ypplay, ypaprs, & 2640 ywake_s, ydTs0, ydqs0, & 2641 yt_x, yt_w, yq_x, yq_w, & 2642 yu_x, yu_w, yv_x, yv_w, & 2643 ycdragh_x, ycdragh_w, ycdragq_x, ycdragq_w, & 2644 ycdragm_x, ycdragm_w, & 2645 AcoefH_x, AcoefH_w, AcoefQ_x, AcoefQ_w, & 2646 AcoefU_x, AcoefU_w, AcoefV_x, AcoefV_w, & 2647 BcoefH_x, BcoefH_w, BcoefQ_x, BcoefQ_w, & 2648 BcoefU_x, BcoefU_w, BcoefV_x, BcoefV_w, & 2649 AcoefH_0, AcoefQ_0, AcoefU, AcoefV, & 2650 BcoefH_0, BcoefQ_0, BcoefU, BcoefV, & 2651 ycdragh, ycdragq, ycdragm, & 2652 yt1, yq1, yu1, yv1 & 2653 ) 2654 IF (iflag_split .eq. 2) THEN 2655 CALL wx_pbl_dts_merge(knon, dtime, ypplay, ypaprs, & 2656 ywake_s, ybeta, ywake_cstar, ywake_dens, & 2657 AcoefH_x, AcoefH_w, & 2658 BcoefH_x, BcoefH_w, & 2659 AcoefH_0, AcoefQ_0, BcoefH_0, BcoefQ_0, & 2660 AcoefH, AcoefQ, BcoefH, BcoefQ, & 2661 HTphiT_b, dd_HTphiT, HTphiQ_b, dd_HTphiQ, HTRn_b, dd_HTRn, & 2662 phiT0_b, dphiT0, phiQ0_b, dphiQ0, Rn0_b, dRn0, & 2663 yg_T, yg_Q, & 2664 yGamma_dTs_phiT, yGamma_dQs_phiQ, & 2665 ydTs_ins, ydqs_ins & 2666 ) 2667 ELSE ! 2668 AcoefH(:) = AcoefH_0(:) 2669 AcoefQ(:) = AcoefQ_0(:) 2670 BcoefH(:) = BcoefH_0(:) 2671 BcoefQ(:) = BcoefQ_0(:) 2672 yg_T(:) = 0. 2673 yg_Q(:) = 0. 2674 yGamma_dTs_phiT(:) = 0. 2675 yGamma_dQs_phiQ(:) = 0. 2676 ydTs_ins(:) = 0. 2677 ydqs_ins(:) = 0. 2678 ENDIF ! (iflag_split .eq. 2) 2679 ! 2680 ELSE ! (nsrf .ne. is_oce) 2681 ybeta(1:knon) = 1. 2682 yevap_pot(1:knon) = yevap(1:knon) 2683 AcoefH(:) = AcoefH_0(:) 2684 AcoefQ(:) = AcoefQ_0(:) 2685 BcoefH(:) = BcoefH_0(:) 2686 BcoefQ(:) = BcoefQ_0(:) 2687 yg_T(:) = 0. 2688 yg_Q(:) = 0. 2689 yGamma_dTs_phiT(:) = 0. 2690 yGamma_dQs_phiQ(:) = 0. 2691 ydTs_ins(:) = 0. 2692 ydqs_ins(:) = 0. 2693 ENDIF ! (nsrf .ne. is_oce) 2694 ! 2695 CALL wx_pbl_split(knon, nsrf, dtime, ywake_s, ybeta, iflag_split, & 2696 yg_T, yg_Q, & 2697 yGamma_dTs_phiT, yGamma_dQs_phiQ, & 2698 ydTs_ins, ydqs_ins, & 2461 2699 y_flux_t1, y_flux_q1, y_flux_u1, y_flux_v1, & 2700 !!!! HTRn_b, dd_HTRn, HTphiT_b, dd_HTphiT, & 2701 phiQ0_b, phiT0_b, & 2462 2702 y_flux_t1_x, y_flux_t1_w, & 2463 2703 y_flux_q1_x, y_flux_q1_w, & … … 2465 2705 y_flux_v1_x, y_flux_v1_w, & 2466 2706 yfluxlat_x, yfluxlat_w, & 2467 y_delta_tsurf & 2707 y_delta_qsats, & 2708 y_delta_tsurf_new, y_delta_qsurf & 2468 2709 ) 2710 ! 2711 CALL wx_pbl_check(knon, dtime, ypplay, ypaprs, ywake_s, ybeta, iflag_split, & 2712 yTs, y_delta_tsurf, & 2713 yqsurf, yTsurf_new, & 2714 y_delta_tsurf_new, y_delta_qsats, & 2715 AcoefH_x, AcoefH_w, & 2716 BcoefH_x, BcoefH_w, & 2717 AcoefH_0, AcoefQ_0, BcoefH_0, BcoefQ_0, & 2718 AcoefH, AcoefQ, BcoefH, BcoefQ, & 2719 y_flux_t1, y_flux_q1, & 2720 y_flux_t1_x, y_flux_t1_w, & 2721 y_flux_q1_x, y_flux_q1_w) 2722 ! 2723 IF (nsrf .ne. is_oce) THEN 2724 CALL wx_pbl_dts_check(knon, dtime, ypplay, ypaprs, ywake_s, ybeta, iflag_split, & 2725 yTs, y_delta_tsurf, & 2726 yqsurf, yTsurf_new, & 2727 y_delta_qsats, y_delta_tsurf_new, y_delta_qsurf, & 2728 AcoefH_x, AcoefH_w, & 2729 BcoefH_x, BcoefH_w, & 2730 AcoefH_0, AcoefQ_0, BcoefH_0, BcoefQ_0, & 2731 AcoefH, AcoefQ, BcoefH, BcoefQ, & 2732 HTphiT_b, dd_HTphiT, HTphiQ_b, dd_HTphiQ, HTRn_b, dd_HTRn, & 2733 phiT0_b, dphiT0, phiQ0_b, dphiQ0, Rn0_b, dRn0, & 2734 yg_T, yg_Q, & 2735 yGamma_dTs_phiT, yGamma_dQs_phiQ, & 2736 ydTs_ins, ydqs_ins, & 2737 y_flux_t1, y_flux_q1, & 2738 y_flux_t1_x, y_flux_t1_w, & 2739 y_flux_q1_x, y_flux_q1_w ) 2740 ENDIF ! (nsrf .ne. is_oce) 2741 ! 2742 ELSE ! (iflag_split .ge. 1) 2743 ybeta(1:knon) = 1. 2744 yevap_pot(1:knon) = yevap(1:knon) 2469 2745 ENDIF ! (iflag_split .ge. 1) 2746 ! 2747 IF (prt_level >= 10) THEN 2748 print *,'pbl_surface, ybeta , yevap, yevap_pot ', & 2749 ybeta , yevap, yevap_pot 2750 ENDIF ! (prt_level >= 10) 2751 ! 2470 2752 !>jyg 2471 2753 ! … … 2648 2930 ENDIF ! (iflag_split .eq.0) 2649 2931 !!! 2650 2651 DO j = 1, knon2652 y_dflux_t(j) = y_dflux_t(j) * ypct(j)2653 y_dflux_q(j) = y_dflux_q(j) * ypct(j)2654 ENDDO2655 2932 !! 2933 !! DO j = 1, knon 2934 !! y_dflux_t(j) = y_dflux_t(j) * ypct(j) 2935 !! y_dflux_q(j) = y_dflux_q(j) * ypct(j) 2936 !! ENDDO 2937 !! 2656 2938 !**************************************************************************************** 2657 2939 ! 13) Transform variables for output format : … … 2821 3103 i = ni(j) 2822 3104 evap(i,nsrf) = - flux_q(i,1,nsrf) !jyg 3105 beta(i,nsrf) = ybeta(j) !jyg 2823 3106 d_ts(i,nsrf) = y_d_ts(j) 2824 3107 !albedo SB >>> … … 2836 3119 cdragh(i) = cdragh(i) + ycdragh(j)*ypct(j) 2837 3120 cdragm(i) = cdragm(i) + ycdragm(j)*ypct(j) 2838 dflux_t(i) = dflux_t(i) + y_dflux_t(j) 2839 dflux_q(i) = dflux_q(i) + y_dflux_q(j) 3121 dflux_t(i) = dflux_t(i) + y_dflux_t(j)*ypct(j) 3122 dflux_q(i) = dflux_q(i) + y_dflux_q(j)*ypct(j) 2840 3123 #ifdef ISO 2841 3124 do ixt=1,niso … … 2844 3127 do ixt=1,ntraciso 2845 3128 xtevap(ixt,i,nsrf) = - flux_xt(ixt,i,1,nsrf) 2846 dflux_xt(ixt,i) = dflux_xt(ixt,i) + y_dflux_xt(ixt,j) 3129 dflux_xt(ixt,i) = dflux_xt(ixt,i) + y_dflux_xt(ixt,j)*ypct(j) 2847 3130 enddo 2848 3131 IF (nsrf == is_lic) THEN … … 2874 3157 !!! nrlmd le 13/06/2011 2875 3158 !!jyg20170131 delta_tsurf(i,nsrf)=y_delta_tsurf(j)*ypct(j) 2876 delta_tsurf(i,nsrf)=y_delta_tsurf(j) 3159 !!jyg20210118 delta_tsurf(i,nsrf)=y_delta_tsurf(j) 3160 delta_tsurf(i,nsrf)=y_delta_tsurf_new(j) 3161 ! 3162 delta_qsurf(i) = delta_qsurf(i) + y_delta_qsurf(j)*ypct(j) 2877 3163 ! 2878 3164 cdragh_x(i) = cdragh_x(i) + ycdragh_x(j)*ypct(j) … … 2918 3204 tke_x(i,k,nsrf) = ytke(j,k) 2919 3205 tke_x(i,k,is_ave) = tke_x(i,k,is_ave) + ytke(j,k)*ypct(j) 3206 2920 3207 !>jyg 2921 3208 ENDDO … … 2931 3218 !! tke(i,k,is_ave) = tke(i,k,is_ave) + tke(i,k,nsrf)*ypct(j) 2932 3219 tke_x(i,k,nsrf) = ytke_x(j,k) 2933 tke_x(i,k,is_ave) = tke_x(i,k,is_ave) + tke_x(i,k,nsrf)*ypct(j) 3220 tke_x(i,k,is_ave) = tke_x(i,k,is_ave) + tke_x(i,k,nsrf)*ypct(j) 2934 3221 wake_dltke(i,k,is_ave) = wake_dltke(i,k,is_ave) + wake_dltke(i,k,nsrf)*ypct(j) 3222 2935 3223 2936 3224 !>jyg … … 3069 3357 d_t_w(:,1), d_t_x(:,1), d_t(:,1) 3070 3358 ENDIF 3359 3360 if (nsrf == is_oce .and. activate_ocean_skin >= 1) then 3361 delta_sal = missing_val 3362 ds_ns = missing_val 3363 dt_ns = missing_val 3364 delta_sst = missing_val 3365 dter = missing_val 3366 dser = missing_val 3367 tkt = missing_val 3368 tks = missing_val 3369 taur = missing_val 3370 sss = missing_val 3371 3372 delta_sal(ni(:knon)) = ydelta_sal(:knon) 3373 ds_ns(ni(:knon)) = yds_ns(:knon) 3374 dt_ns(ni(:knon)) = ydt_ns(:knon) 3375 delta_sst(ni(:knon)) = ydelta_sst(:knon) 3376 dter(ni(:knon)) = ydter(:knon) 3377 dser(ni(:knon)) = ydser(:knon) 3378 tkt(ni(:knon)) = ytkt(:knon) 3379 tks(ni(:knon)) = ytks(:knon) 3380 taur(ni(:knon)) = ytaur(:knon) 3381 sss(ni(:knon)) = ysss(:knon) 3382 end if 3383 3384 3071 3385 3072 3386 !**************************************************************************************** … … 3106 3420 * (ypaprs(j,1)-ypplay(j,1)) 3107 3421 tairsol(j) = yts(j) + y_d_ts(j) 3108 tairsol_x(j) = tairsol(j) - ywake_s(j)*y_delta_tsurf(j) 3422 !! tairsol_x(j) = tairsol(j) - ywake_s(j)*y_delta_tsurf(j) 3423 tairsol_x(j) = tairsol(j) - ywake_s(j)*y_delta_tsurf_new(j) 3109 3424 qairsol(j) = yqsurf(j) 3110 3425 ENDDO … … 3145 3460 !!! jyg le 07/02/2012 3146 3461 IF (iflag_split .eq.0) THEN 3462 IF (iflag_new_t2mq2m==1) THEN 3463 CALL stdlevvarn(klon, knon, nsrf, zxli, & 3464 uzon, vmer, tair1, qair1, zgeo1, & 3465 tairsol, qairsol, yz0m, yz0h_oupas, psfce, patm, & 3466 yt2m, yq2m, yt10m, yq10m, yu10m, yustar, & 3467 yn2mout(:, nsrf, :)) 3468 ELSE 3147 3469 CALL stdlevvar(klon, knon, nsrf, zxli, & 3148 3470 uzon, vmer, tair1, qair1, zgeo1, & 3149 3471 tairsol, qairsol, yz0m, yz0h_oupas, psfce, patm, & 3150 3472 yt2m, yq2m, yt10m, yq10m, yu10m, yustar) 3473 ENDIF 3151 3474 ELSE !(iflag_split .eq.0) 3475 IF (iflag_new_t2mq2m==1) THEN 3476 CALL stdlevvarn(klon, knon, nsrf, zxli, & 3477 uzon_x, vmer_x, tair1_x, qair1_x, zgeo1_x, & 3478 tairsol_x, qairsol, yz0m, yz0h_oupas, psfce, patm, & 3479 yt2m_x, yq2m_x, yt10m_x, yq10m_x, yu10m_x, yustar_x, & 3480 yn2mout_x(:, nsrf, :)) 3481 CALL stdlevvarn(klon, knon, nsrf, zxli, & 3482 uzon_w, vmer_w, tair1_w, qair1_w, zgeo1_w, & 3483 tairsol_w, qairsol, yz0m, yz0h_oupas, psfce, patm, & 3484 yt2m_w, yq2m_w, yt10m_w, yq10m_w, yu10m_w, yustar_w, & 3485 yn2mout_w(:, nsrf, :)) 3486 ELSE 3152 3487 CALL stdlevvar(klon, knon, nsrf, zxli, & 3153 3488 uzon_x, vmer_x, tair1_x, qair1_x, zgeo1_x, & … … 3158 3493 tairsol_w, qairsol, yz0m, yz0h_oupas, psfce, patm, & 3159 3494 yt2m_w, yq2m_w, yt10m_w, yq10m_w, yu10m_w, yustar_w) 3495 ENDIF 3160 3496 !!! 3161 3497 ENDIF ! (iflag_split .eq.0) … … 3171 3507 u10m(i,nsrf)=(yu10m(j) * uzon(j))/SQRT(uzon(j)**2+vmer(j)**2) 3172 3508 v10m(i,nsrf)=(yu10m(j) * vmer(j))/SQRT(uzon(j)**2+vmer(j)**2) 3509 ! 3510 DO k = 1, 6 3511 n2mout(i,nsrf,k) = yn2mout(j,nsrf,k) 3512 END DO 3513 ! 3173 3514 ENDDO 3174 3515 ELSE !(iflag_split .eq.0) … … 3181 3522 u10m_x(i,nsrf)=(yu10m_x(j) * uzon_x(j))/SQRT(uzon_x(j)**2+vmer_x(j)**2) 3182 3523 v10m_x(i,nsrf)=(yu10m_x(j) * vmer_x(j))/SQRT(uzon_x(j)**2+vmer_x(j)**2) 3524 ! 3525 DO k = 1, 6 3526 n2mout_x(i,nsrf,k) = yn2mout_x(j,nsrf,k) 3527 END DO 3528 ! 3183 3529 ENDDO 3184 3530 DO j=1, knon … … 3194 3540 u10m(i,nsrf) = u10m_x(i,nsrf) + wake_s(i)*(u10m_w(i,nsrf)-u10m_x(i,nsrf)) 3195 3541 v10m(i,nsrf) = v10m_x(i,nsrf) + wake_s(i)*(v10m_w(i,nsrf)-v10m_x(i,nsrf)) 3542 ! 3543 DO k = 1, 6 3544 n2mout_w(i,nsrf,k) = yn2mout_w(j,nsrf,k) 3545 END DO 3546 ! 3196 3547 ENDDO 3197 3548 !!! … … 3479 3830 #endif 3480 3831 !!! 3481 3832 3482 3833 ! 3483 3834 ! Incrementer la temperature du sol 3484 3835 ! 3485 3836 zxtsol(:) = 0.0 ; zxfluxlat(:) = 0.0 3486 zt2m(:) = 0.0 ; zq2m(:) = 0.0 3837 zt2m(:) = 0.0 ; zq2m(:) = 0.0 ; zn2mout(:,:) = 0 3487 3838 zustar(:)=0.0 ; zu10m(:) = 0.0 ; zv10m(:) = 0.0 3488 3839 s_pblh(:) = 0.0 ; s_plcl(:) = 0.0 … … 3537 3888 zt2m(i) = zt2m(i) + t2m(i,nsrf) * pctsrf(i,nsrf) 3538 3889 zq2m(i) = zq2m(i) + q2m(i,nsrf) * pctsrf(i,nsrf) 3890 ! 3891 DO k = 1, 6 3892 zn2mout(i,k) = zn2mout(i,k) + n2mout(i,nsrf,k) * pctsrf(i,nsrf) 3893 ENDDO 3894 ! 3539 3895 zustar(i) = zustar(i) + ustar(i,nsrf) * pctsrf(i,nsrf) 3540 3896 wstar(i,is_ave)=wstar(i,is_ave)+wstar(i,nsrf)*pctsrf(i,nsrf) … … 3568 3924 zt2m(i) = zt2m(i) + (t2m_x(i,nsrf)+wake_s(i)*(t2m_w(i,nsrf)-t2m_x(i,nsrf))) * pctsrf(i,nsrf) 3569 3925 zq2m(i) = zq2m(i) + q2m_x(i,nsrf) * pctsrf(i,nsrf) 3926 ! 3927 DO k = 1, 6 3928 zn2mout(i,k) = zn2mout(i,k) + n2mout_x(i,nsrf,k) * pctsrf(i,nsrf) 3929 ENDDO 3930 ! 3570 3931 zustar(i) = zustar(i) + ustar_x(i,nsrf) * pctsrf(i,nsrf) 3571 3932 wstar(i,is_ave)=wstar(i,is_ave)+wstar_x(i,nsrf)*pctsrf(i,nsrf) … … 3649 4010 DO nsrf = 1, nbsrf 3650 4011 DO i = 1, klon 3651 zxqsurf(i) = zxqsurf(i) + qsurf(i,nsrf) * pctsrf(i,nsrf)4012 zxqsurf(i) = zxqsurf(i) + MAX(qsurf(i,nsrf),0.0) * pctsrf(i,nsrf) 3652 4013 zxsnow(i) = zxsnow(i) + snow(i,nsrf) * pctsrf(i,nsrf) 3653 4014 #ifdef ISO … … 3718 4079 IF (ALLOCATED(qsurf)) DEALLOCATE(qsurf) 3719 4080 IF (ALLOCATED(ftsoil)) DEALLOCATE(ftsoil) 4081 IF (ALLOCATED(ydTs0)) DEALLOCATE(ydTs0) 4082 IF (ALLOCATED(ydqs0)) DEALLOCATE(ydqs0) 3720 4083 #ifdef ISO 3721 4084 IF (ALLOCATED(xtsnow)) DEALLOCATE(xtsnow) … … 3739 4102 3740 4103 !albedo SB >>> 3741 SUBROUTINE pbl_surface_newfrac(itime, pctsrf_new, pctsrf_old, &3742 evap, z0m, z0h, agesno, &3743 tsurf,alb_dir,alb_dif, ustar, u10m, v10m, tke &4104 SUBROUTINE pbl_surface_newfrac(itime, pctsrf_new, pctsrf_old, & 4105 evap, z0m, z0h, agesno, & 4106 tsurf,alb_dir,alb_dif, ustar, u10m, v10m, tke & 3744 4107 #ifdef ISO 3745 4108 ,xtevap & … … 3750 4113 3751 4114 USE indice_sol_mod 4115 use phys_state_var_mod, only: delta_sal, ds_ns, dt_ns, delta_sst 4116 use config_ocean_skin_m, only: activate_ocean_skin 3752 4117 #ifdef ISO 3753 4118 USE infotrac_phy, ONLY: ntraciso … … 3852 4217 alb_dif(i,k,nsrf) = 0.06 3853 4218 ENDDO 4219 if (activate_ocean_skin >= 1) then 4220 if (activate_ocean_skin == 2 & 4221 .and. type_ocean == "couple") then 4222 delta_sal(i) = 0. 4223 delta_sst(i) = 0. 4224 end if 4225 4226 ds_ns(i) = 0. 4227 dt_ns(i) = 0. 4228 end if 3854 4229 ELSE IF (nsrf.EQ.is_sic) THEN 3855 4230 tsurf(i,nsrf) = 271.35
Note: See TracChangeset
for help on using the changeset viewer.
