Changeset 996 for LMDZ4/trunk/libf/phylmd/pbl_surface_mod.F90
- Timestamp:
- Sep 9, 2008, 3:22:23 PM (16 years ago)
- File:
-
- 1 edited
-
LMDZ4/trunk/libf/phylmd/pbl_surface_mod.F90 (modified) (33 diffs)
Legend:
- Unmodified
- Added
- Removed
-
LMDZ4/trunk/libf/phylmd/pbl_surface_mod.F90
r987 r996 13 13 USE mod_phys_lmdz_para, ONLY : mpi_size 14 14 USE ioipsl 15 USE surface_data, ONLY : ocean, ok_veget15 USE surface_data, ONLY : type_ocean, ok_veget 16 16 USE surf_land_mod, ONLY : surf_land 17 17 USE surf_landice_mod, ONLY : surf_landice … … 151 151 !**************************************************************************************** 152 152 153 IF ( ocean /= 'slab ' .AND. ocean /= 'force ' .AND.ocean /= 'couple') THEN153 IF (type_ocean /= 'slab ' .AND. type_ocean /= 'force ' .AND. type_ocean /= 'couple') THEN 154 154 WRITE(lunout,*)' *** Warning ***' 155 WRITE(lunout,*)'Option couplage pour l''ocean = ', ocean155 WRITE(lunout,*)'Option couplage pour l''ocean = ', type_ocean 156 156 abort_message='option pour l''ocean non valable' 157 157 CALL abort_gcm(modname,abort_message,1) … … 187 187 zxtsol, zxfluxlat, zt2m, qsat2m, & 188 188 d_t, d_q, d_u, d_v, & 189 zcoefh, pctsrf_new,&189 zcoefh, slab_wfbils, & 190 190 qsol_d, zq2m, s_pblh, s_plcl, & 191 191 s_capCL, s_oliqCL, s_cteiCL, s_pblT, & … … 299 299 REAL, DIMENSION(klon, nbsrf), INTENT(INOUT) :: u10m ! u speed at 10m 300 300 REAL, DIMENSION(klon, nbsrf), INTENT(INOUT) :: v10m ! v speed at 10m 301 REAL, DIMENSION(klon, klev+1, nbsrf), INTENT(INOUT) :: tke 301 302 302 303 ! Output variables … … 321 322 REAL, DIMENSION(klon, klev), INTENT(OUT) :: d_v ! change in v speed 322 323 REAL, DIMENSION(klon, klev), INTENT(OUT) :: zcoefh ! coef for turbulent diffusion of T and Q, mean for each grid point 323 REAL, DIMENSION(klon, nbsrf), INTENT(OUT) :: pctsrf_new ! new sub-surface fraction324 324 325 325 ! Output only for diagnostics 326 REAL, DIMENSION(klon), INTENT(OUT) :: slab_wfbils! heat balance at surface only for slab at ocean points 326 327 REAL, DIMENSION(klon), INTENT(OUT) :: qsol_d ! water height in the soil (mm) 327 328 REAL, DIMENSION(klon), INTENT(OUT) :: zq2m ! water vapour at 2m, mean for each grid point … … 367 368 REAL, DIMENSION(klon, nbsrf),INTENT(OUT) :: q2m ! water vapour at 2 meter height 368 369 REAL, DIMENSION(klon, klev, nbsrf), INTENT(OUT) :: flux_q ! water vapour flux(latent flux) (kg/m**2/s) 369 370 ! Input/output371 REAL, DIMENSION(klon, klev+1, nbsrf), INTENT(INOUT) :: tke372 370 373 371 … … 431 429 REAL, DIMENSION(klon) :: ypsref 432 430 REAL, DIMENSION(klon) :: yevap, ytsurf_new, yalb1_new, yalb2_new 433 REAL, DIMENSION(klon) :: pctsrf_nsrf434 431 REAL, DIMENSION(klon) :: ztsol 435 432 REAL, DIMENSION(klon) :: alb_m ! mean albedo for whole SW interval … … 446 443 REAL, DIMENSION(klon,klev+1) :: ytke 447 444 REAL, DIMENSION(klon,nsoilmx) :: ytsoil 448 REAL, DIMENSION(klon,nbsrf) :: pctsrf_pot449 445 CHARACTER(len=80) :: abort_message 450 446 CHARACTER(len=20) :: modname = 'pbl_surface' … … 470 466 REAL, DIMENSION(klon,nbsrf) :: trmb2 ! inhibition 471 467 REAL, DIMENSION(klon,nbsrf) :: trmb3 ! point Omega 472 REAL, DIMENSION(klon,nbsrf) :: zx_rh2m, zx_qsat2m 473 REAL, DIMENSION(klon,nbsrf) :: zx_qs1, zx_t1 474 REAL, DIMENSION(klon,nbsrf) :: zdelta1, zcor1 468 REAL, DIMENSION(klon,nbsrf) :: zx_t1 475 469 REAL, DIMENSION(klon, nbsrf) :: alb ! mean albedo for whole SW interval 476 470 REAL, DIMENSION(klon) :: ylwdown ! jg : temporary (ysollwdown) 477 471 478 479 !jg+ temporary test 480 REAL, DIMENSION(klon,klev) :: y_flux_u_old, y_flux_v_old 481 REAL, DIMENSION(klon,klev) :: y_d_u_old, y_d_v_old 482 !jg- 483 472 REAL :: zx_qs1, zcor1, zdelta1 473 484 474 !**************************************************************************************** 485 475 ! Declarations specifiques pour le 1D. A reprendre … … 558 548 yv1 = 0.0 ; ypaprs = 0.0 ; ypplay = 0.0 559 549 ydelp = 0.0 ; yu = 0.0 ; yv = 0.0 ; yt = 0.0 560 yq = 0.0 ; pctsrf_new = 0.0 ; y_dflux_t = 0.0; y_dflux_q = 0.0550 yq = 0.0 ; y_dflux_t = 0.0 ; y_dflux_q = 0.0 561 551 yrugoro = 0.0 ; yu10mx = 0.0 ; yu10my = 0.0 ; ywindsp = 0.0 562 552 d_ts = 0.0 ; yfluxlat=0.0 ; flux_t = 0.0 ; flux_q = 0.0 … … 661 651 ! 4) Loop over different surfaces 662 652 ! 663 ! All points with a possibility of the current surface are used. This is done 664 ! to allow the sea-ice to appear or disappear. It is considered here that the 665 ! entier domaine of the ocean possibly can contain sea-ice. 653 ! Only points containing a fraction of the sub surface will be threated. 666 654 ! 667 655 !**************************************************************************************** 668 669 pctsrf_pot = pctsrf 670 pctsrf_pot(:,is_oce) = 1. - zmasq(:) 671 pctsrf_pot(:,is_sic) = 1. - zmasq(:) 672 656 673 657 loop_nbsrf: DO nsrf = 1, nbsrf 674 658 … … 677 661 knon = 0 678 662 DO i = 1, klon 679 IF (pctsrf _pot(i,nsrf).GT.epsfra) THEN663 IF (pctsrf(i,nsrf) > 0.) THEN 680 664 knon = knon + 1 681 665 ni(knon) = i … … 730 714 ypaprs(j,k) = paprs(i,k) 731 715 ypplay(j,k) = pplay(i,k) 732 ydelp(j,k) = delp(i,k)733 ytke(j,k) =tke(i,k,nsrf)716 ydelp(j,k) = delp(i,k) 717 ytke(j,k) = tke(i,k,nsrf) 734 718 yu(j,k) = u(i,k) 735 719 yv(j,k) = v(i,k) … … 762 746 CALL coef_diff_turb(dtime, nsrf, knon, ni, & 763 747 ypaprs, ypplay, yu, yv, yq, yt, yts, yrugos, yqsurf, & 764 ycoefm, ycoefh, ytke)748 ycoefm, ycoefh, ytke) 765 749 766 750 !**************************************************************************************** … … 817 801 ysnow, yqsol, yagesno, ytsoil, & 818 802 yz0_new, yalb1_new, yalb2_new, yevap, yfluxsens, yfluxlat, & 819 yqsurf, ytsurf_new, y_dflux_t, y_dflux_q, pctsrf_nsrf,&803 yqsurf, ytsurf_new, y_dflux_t, y_dflux_q, & 820 804 ylwdown) 821 805 … … 828 812 ysnow, yqsurf, yqsol, yagesno, & 829 813 ytsoil, yz0_new, yalb1_new, yalb2_new, yevap, yfluxsens, yfluxlat, & 830 ytsurf_new, y_dflux_t, y_dflux_q , pctsrf_nsrf)814 ytsurf_new, y_dflux_t, y_dflux_q) 831 815 832 816 CASE(is_oce) 833 817 CALL surf_ocean(rlon, rlat, ysolsw, ysollw, yalb1, & 834 yrugos, ywindsp, rmu0, yfder, &818 yrugos, ywindsp, rmu0, yfder, yts, & 835 819 itap, dtime, jour, knon, ni, & 836 820 debut, & … … 840 824 ysnow, yqsurf, yagesno, & 841 825 yz0_new, yalb1_new, yalb2_new, yevap, yfluxsens, yfluxlat, & 842 ytsurf_new, y_dflux_t, y_dflux_q, pctsrf_nsrf)826 ytsurf_new, y_dflux_t, y_dflux_q, slab_wfbils) 843 827 844 828 CASE(is_sic) … … 852 836 ysnow, yqsurf, yqsol, yagesno, ytsoil, & 853 837 yz0_new, yalb1_new, yalb2_new, yevap, yfluxsens, yfluxlat, & 854 ytsurf_new, y_dflux_t, y_dflux_q , pctsrf_nsrf)838 ytsurf_new, y_dflux_t, y_dflux_q) 855 839 856 840 … … 861 845 END SELECT 862 846 863 !****************************************************************************************864 ! Save the fraction of this sub-surface865 !866 !****************************************************************************************867 pctsrf_new(:,nsrf) = pctsrf_nsrf(:)868 847 869 848 !**************************************************************************************** … … 882 861 !**************************************************************************************** 883 862 ! H and Q 884 ! print *,'pbl_surface: ok_flux_surf=',ok_flux_surf885 ! print *,'pbl_surface: fsens flat RLVTT=',fsens,flat,RLVTT863 ! print *,'pbl_surface: ok_flux_surf=',ok_flux_surf 864 ! print *,'pbl_surface: fsens flat RLVTT=',fsens,flat,RLVTT 886 865 if (ok_flux_surf) then 887 866 y_flux_t1(:) = fsens … … 916 895 !**************************************************************************************** 917 896 918 tke(:,:,nsrf) =0.897 tke(:,:,nsrf) = 0. 919 898 DO k = 1, klev 920 899 DO j = 1, knon … … 922 901 ycoefh(j,k) = ycoefh(j,k) * ypct(j) 923 902 ycoefm(j,k) = ycoefm(j,k) * ypct(j) 924 y_d_t(j,k) = y_d_t(j,k) * ypct(j)925 y_d_q(j,k) = y_d_q(j,k) * ypct(j)926 y_d_u(j,k) = y_d_u(j,k) * ypct(j)927 y_d_v(j,k) = y_d_v(j,k) * ypct(j)903 y_d_t(j,k) = y_d_t(j,k) * ypct(j) 904 y_d_q(j,k) = y_d_q(j,k) * ypct(j) 905 y_d_u(j,k) = y_d_u(j,k) * ypct(j) 906 y_d_v(j,k) = y_d_v(j,k) * ypct(j) 928 907 929 908 flux_t(i,k,nsrf) = y_flux_t(j,k) … … 932 911 flux_v(i,k,nsrf) = y_flux_v(j,k) 933 912 934 tke(i,k,nsrf) =ytke(j,k)913 tke(i,k,nsrf) = ytke(j,k) 935 914 936 915 ENDDO … … 1021 1000 trmb2(:,nsrf) = 0. ! inhibition 1022 1001 trmb3(:,nsrf) = 0. ! Point Omega 1002 rh2m(:) = 0. 1003 qsat2m(:) = 0. 1023 1004 1024 1005 #undef T2m 1025 1006 #define T2m 1026 1007 #ifdef T2m 1027 ! diagnostic t,q a 2m et u, v a10m1008 ! Calculations of diagnostic t,q at 2m and u, v at 10m 1028 1009 1029 1010 DO j=1, knon … … 1055 1036 i = ni(j) 1056 1037 t2m(i,nsrf)=yt2m(j) 1038 q2m(i,nsrf)=yq2m(j) 1057 1039 1058 q2m(i,nsrf)=yq2m(j) 1059 1060 ! u10m, v10m : composantes du vent a 10m sans spirale de Ekman 1040 ! u10m, v10m : composantes du vent a 10m sans spirale de Ekman 1061 1041 u10m(i,nsrf)=(yu10m(j) * uzon(j))/SQRT(uzon(j)**2+vmer(j)**2) 1062 1042 v10m(i,nsrf)=(yu10m(j) * vmer(j))/SQRT(uzon(j)**2+vmer(j)**2) 1063 1064 1043 END DO 1065 1044 1045 !IM Calcule de l'humidite relative a 2m (rh2m) pour diagnostique 1046 !IM Ajoute dependance type surface 1047 IF (thermcep) THEN 1048 DO j = 1, knon 1049 i=ni(j) 1050 zdelta1 = MAX(0.,SIGN(1., rtt-yt2m(j) )) 1051 zx_qs1 = r2es * FOEEW(yt2m(j),zdelta1)/paprs(i,1) 1052 zx_qs1 = MIN(0.5,zx_qs1) 1053 zcor1 = 1./(1.-RETV*zx_qs1) 1054 zx_qs1 = zx_qs1*zcor1 1055 1056 rh2m(i) = rh2m(i) + yq2m(j)/zx_qs1 * pctsrf(i,nsrf) 1057 qsat2m(i) = qsat2m(i) + zx_qs1 * pctsrf(i,nsrf) 1058 END DO 1059 END IF 1066 1060 1067 1061 CALL HBTM(knon, ypaprs, ypplay, & … … 1086 1080 1087 1081 #else 1088 ! not defined T2m1082 ! T2m not defined 1089 1083 ! No calculation 1090 ! Set output variables to zero 1091 DO j = 1, knon 1092 i = ni(j) 1093 pblh(i,nsrf) = 0. 1094 plcl(i,nsrf) = 0. 1095 capCL(i,nsrf) = 0. 1096 oliqCL(i,nsrf) = 0. 1097 cteiCL(i,nsrf) = 0. 1098 pblT(i,nsrf) = 0. 1099 therm(i,nsrf) = 0. 1100 trmb1(i,nsrf) = 0. 1101 trmb2(i,nsrf) = 0. 1102 trmb3(i,nsrf) = 0. 1103 END DO 1104 DO j = 1, knon 1105 i = ni(j) 1106 t2m(i,nsrf)=0. 1107 q2m(i,nsrf)=0. 1108 u10m(i,nsrf)=0. 1109 v10m(i,nsrf)=0. 1110 END DO 1084 PRINT*,' Warning !!! No T2m calculation. Output is set to zero.' 1111 1085 #endif 1112 1086 … … 1120 1094 ! 16) Calculate the mean value over all sub-surfaces for som variables 1121 1095 ! 1122 ! NB!!! jg : Pour garder la convergence numerique j'utilise pctsrf_new comme c'etait1123 ! fait dans physiq.F mais ca devrait plutot etre pctsrf???!!!!! A verifier!1124 1096 !**************************************************************************************** 1125 1097 … … 1129 1101 DO k = 1, klev 1130 1102 DO i = 1, klon 1131 zxfluxt(i,k) = zxfluxt(i,k) + flux_t(i,k,nsrf) * pctsrf _new(i,nsrf)1132 zxfluxq(i,k) = zxfluxq(i,k) + flux_q(i,k,nsrf) * pctsrf _new(i,nsrf)1133 zxfluxu(i,k) = zxfluxu(i,k) + flux_u(i,k,nsrf) * pctsrf _new(i,nsrf)1134 zxfluxv(i,k) = zxfluxv(i,k) + flux_v(i,k,nsrf) * pctsrf _new(i,nsrf)1103 zxfluxt(i,k) = zxfluxt(i,k) + flux_t(i,k,nsrf) * pctsrf(i,nsrf) 1104 zxfluxq(i,k) = zxfluxq(i,k) + flux_q(i,k,nsrf) * pctsrf(i,nsrf) 1105 zxfluxu(i,k) = zxfluxu(i,k) + flux_u(i,k,nsrf) * pctsrf(i,nsrf) 1106 zxfluxv(i,k) = zxfluxv(i,k) + flux_v(i,k,nsrf) * pctsrf(i,nsrf) 1135 1107 END DO 1136 1108 END DO … … 1143 1115 ENDDO 1144 1116 1145 1146 DO i = 1, klon1147 IF ( ABS( pctsrf_new(i, is_ter) + pctsrf_new(i, is_lic) + &1148 pctsrf_new(i, is_oce) + pctsrf_new(i, is_sic) - 1.) .GT. EPSFRA) &1149 THEN1150 WRITE(*,*) 'physiq : pb sous surface au point ', i, &1151 pctsrf_new(i, 1 : nbsrf)1152 ENDIF1153 ENDDO1154 1155 1117 ! 1156 1118 ! Incrementer la temperature du sol … … 1171 1133 1172 1134 wfbils(i,nsrf) = ( solsw(i,nsrf) + sollw(i,nsrf) & 1173 + flux_t(i,1,nsrf) + fluxlat(i,nsrf) ) * pctsrf _new(i,nsrf)1135 + flux_t(i,1,nsrf) + fluxlat(i,nsrf) ) * pctsrf(i,nsrf) 1174 1136 wfbilo(i,nsrf) = (evap(i,nsrf) - (rain_f(i) + snow_f(i))) * & 1175 pctsrf _new(i,nsrf)1176 1177 zxtsol(i) = zxtsol(i) + ts(i,nsrf) * pctsrf _new(i,nsrf)1178 zxfluxlat(i) = zxfluxlat(i) + fluxlat(i,nsrf) * pctsrf _new(i,nsrf)1137 pctsrf(i,nsrf) 1138 1139 zxtsol(i) = zxtsol(i) + ts(i,nsrf) * pctsrf(i,nsrf) 1140 zxfluxlat(i) = zxfluxlat(i) + fluxlat(i,nsrf) * pctsrf(i,nsrf) 1179 1141 1180 zt2m(i) = zt2m(i) + t2m(i,nsrf) * pctsrf _new(i,nsrf)1181 zq2m(i) = zq2m(i) + q2m(i,nsrf) * pctsrf _new(i,nsrf)1182 zu10m(i) = zu10m(i) + u10m(i,nsrf) * pctsrf _new(i,nsrf)1183 zv10m(i) = zv10m(i) + v10m(i,nsrf) * pctsrf _new(i,nsrf)1184 1185 s_pblh(i) = s_pblh(i) + pblh(i,nsrf) * pctsrf _new(i,nsrf)1186 s_plcl(i) = s_plcl(i) + plcl(i,nsrf) * pctsrf _new(i,nsrf)1187 s_capCL(i) = s_capCL(i) + capCL(i,nsrf) * pctsrf _new(i,nsrf)1188 s_oliqCL(i) = s_oliqCL(i) + oliqCL(i,nsrf)* pctsrf _new(i,nsrf)1189 s_cteiCL(i) = s_cteiCL(i) + cteiCL(i,nsrf)* pctsrf _new(i,nsrf)1190 s_pblT(i) = s_pblT(i) + pblT(i,nsrf) * pctsrf _new(i,nsrf)1191 s_therm(i) = s_therm(i) + therm(i,nsrf) * pctsrf _new(i,nsrf)1192 s_trmb1(i) = s_trmb1(i) + trmb1(i,nsrf) * pctsrf _new(i,nsrf)1193 s_trmb2(i) = s_trmb2(i) + trmb2(i,nsrf) * pctsrf _new(i,nsrf)1194 s_trmb3(i) = s_trmb3(i) + trmb3(i,nsrf) * pctsrf _new(i,nsrf)1142 zt2m(i) = zt2m(i) + t2m(i,nsrf) * pctsrf(i,nsrf) 1143 zq2m(i) = zq2m(i) + q2m(i,nsrf) * pctsrf(i,nsrf) 1144 zu10m(i) = zu10m(i) + u10m(i,nsrf) * pctsrf(i,nsrf) 1145 zv10m(i) = zv10m(i) + v10m(i,nsrf) * pctsrf(i,nsrf) 1146 1147 s_pblh(i) = s_pblh(i) + pblh(i,nsrf) * pctsrf(i,nsrf) 1148 s_plcl(i) = s_plcl(i) + plcl(i,nsrf) * pctsrf(i,nsrf) 1149 s_capCL(i) = s_capCL(i) + capCL(i,nsrf) * pctsrf(i,nsrf) 1150 s_oliqCL(i) = s_oliqCL(i) + oliqCL(i,nsrf)* pctsrf(i,nsrf) 1151 s_cteiCL(i) = s_cteiCL(i) + cteiCL(i,nsrf)* pctsrf(i,nsrf) 1152 s_pblT(i) = s_pblT(i) + pblT(i,nsrf) * pctsrf(i,nsrf) 1153 s_therm(i) = s_therm(i) + therm(i,nsrf) * pctsrf(i,nsrf) 1154 s_trmb1(i) = s_trmb1(i) + trmb1(i,nsrf) * pctsrf(i,nsrf) 1155 s_trmb2(i) = s_trmb2(i) + trmb2(i,nsrf) * pctsrf(i,nsrf) 1156 s_trmb3(i) = s_trmb3(i) + trmb3(i,nsrf) * pctsrf(i,nsrf) 1195 1157 END DO 1196 1158 END DO … … 1205 1167 PRINT*,' debut apres d_ts min max ftsol(ts)',itap,amn,amx 1206 1168 ENDIF 1207 ! 1208 ! If a sub-surface does not exsist for a grid point, the mean value for all1209 ! sub-surfaces is distributed.1210 ! 1211 DO nsrf = 1, nbsrf1212 DO i = 1, klon1213 IF ((pctsrf_new(i,nsrf) .LT. epsfra) .OR. (t2m(i,nsrf).EQ.0.)) THEN1214 ts(i,nsrf) = zxtsol(i)1215 t2m(i,nsrf) = zt2m(i)1216 q2m(i,nsrf) = zq2m(i)1217 u10m(i,nsrf) = zu10m(i)1218 v10m(i,nsrf) = zv10m(i)1219 1220 ! Les variables qui suivent sont plus utilise, donc peut-etre pas la peine a les mettre ajour1221 pblh(i,nsrf) = s_pblh(i)1222 plcl(i,nsrf) = s_plcl(i)1223 capCL(i,nsrf) = s_capCL(i)1224 oliqCL(i,nsrf) = s_oliqCL(i)1225 cteiCL(i,nsrf) = s_cteiCL(i)1226 pblT(i,nsrf) = s_pblT(i)1227 therm(i,nsrf) = s_therm(i)1228 trmb1(i,nsrf) = s_trmb1(i)1229 trmb2(i,nsrf) = s_trmb2(i)1230 trmb3(i,nsrf) = s_trmb3(i)1231 ENDIF1232 ENDDO1233 ENDDO1169 !!$! 1170 !!$! If a sub-surface does not exsist for a grid point, the mean value for all 1171 !!$! sub-surfaces is distributed. 1172 !!$! 1173 !!$ DO nsrf = 1, nbsrf 1174 !!$ DO i = 1, klon 1175 !!$ IF ((pctsrf_new(i,nsrf) .LT. epsfra) .OR. (t2m(i,nsrf).EQ.0.)) THEN 1176 !!$ ts(i,nsrf) = zxtsol(i) 1177 !!$ t2m(i,nsrf) = zt2m(i) 1178 !!$ q2m(i,nsrf) = zq2m(i) 1179 !!$ u10m(i,nsrf) = zu10m(i) 1180 !!$ v10m(i,nsrf) = zv10m(i) 1181 !!$ 1182 !!$! Les variables qui suivent sont plus utilise, donc peut-etre pas la peine a les mettre ajour 1183 !!$ pblh(i,nsrf) = s_pblh(i) 1184 !!$ plcl(i,nsrf) = s_plcl(i) 1185 !!$ capCL(i,nsrf) = s_capCL(i) 1186 !!$ oliqCL(i,nsrf) = s_oliqCL(i) 1187 !!$ cteiCL(i,nsrf) = s_cteiCL(i) 1188 !!$ pblT(i,nsrf) = s_pblT(i) 1189 !!$ therm(i,nsrf) = s_therm(i) 1190 !!$ trmb1(i,nsrf) = s_trmb1(i) 1191 !!$ trmb2(i,nsrf) = s_trmb2(i) 1192 !!$ trmb3(i,nsrf) = s_trmb3(i) 1193 !!$ ENDIF 1194 !!$ ENDDO 1195 !!$ ENDDO 1234 1196 1235 1197 … … 1242 1204 DO nsrf = 1, nbsrf 1243 1205 DO i = 1, klon 1244 zxqsurf(i) = zxqsurf(i) + qsurf(i,nsrf) * pctsrf _new(i,nsrf)1245 zxsnow(i) = zxsnow(i) + snow(i,nsrf) * pctsrf _new(i,nsrf)1206 zxqsurf(i) = zxqsurf(i) + qsurf(i,nsrf) * pctsrf(i,nsrf) 1207 zxsnow(i) = zxsnow(i) + snow(i,nsrf) * pctsrf(i,nsrf) 1246 1208 END DO 1247 1209 END DO 1248 1210 1249 !1250 !IM Calculer l'humidite relative a 2m (rh2m) pour diagnostique1251 !IM ajout dependance type surface1252 rh2m(:) = 0.01253 qsat2m(:) = 0.01254 1255 DO i = 1, klon1256 DO nsrf=1, nbsrf1257 zx_t1(i,nsrf) = t2m(i,nsrf)1258 IF (thermcep) THEN1259 zdelta1(i,nsrf) = MAX(0.,SIGN(1.,rtt-zx_t1(i,nsrf)))1260 zx_qs1(i,nsrf) = r2es * &1261 FOEEW(zx_t1(i,nsrf),zdelta1(i,nsrf))/paprs(i,1)1262 zx_qs1(i,nsrf) = MIN(0.5,zx_qs1(i,nsrf))1263 zcor1(i,nsrf) = 1./(1.-retv*zx_qs1(i,nsrf))1264 zx_qs1(i,nsrf) = zx_qs1(i,nsrf)*zcor1(i,nsrf)1265 END IF1266 zx_rh2m(i,nsrf) = q2m(i,nsrf)/zx_qs1(i,nsrf)1267 zx_qsat2m(i,nsrf)=zx_qs1(i,nsrf)1268 rh2m(i) = rh2m(i)+zx_rh2m(i,nsrf)*pctsrf_new(i,nsrf)1269 qsat2m(i)=qsat2m(i)+zx_qsat2m(i,nsrf)*pctsrf_new(i,nsrf)1270 END DO1271 END DO1272 1211 1273 1212 ! Some of the module declared variables are returned for printing in physiq.F … … 1322 1261 ! 1323 1262 !**************************************************************************************** 1263 ! 1264 SUBROUTINE pbl_surface_newfrac(itime, pctsrf_new, pctsrf_old, tsurf, alb1, alb2, u10m, v10m, tke) 1265 1266 ! Give default values where new fraction has appread 1267 1268 INCLUDE "indicesol.h" 1269 INCLUDE "dimsoil.h" 1270 INCLUDE "clesphys.h" 1271 1272 ! Input variables 1273 !**************************************************************************************** 1274 INTEGER, INTENT(IN) :: itime 1275 REAL, DIMENSION(klon,nbsrf), INTENT(IN) :: pctsrf_new, pctsrf_old 1276 1277 ! InOutput variables 1278 !**************************************************************************************** 1279 REAL, DIMENSION(klon,nbsrf), INTENT(INOUT) :: tsurf 1280 REAL, DIMENSION(klon,nbsrf), INTENT(INOUT) :: alb1, alb2 1281 REAL, DIMENSION(klon,nbsrf), INTENT(INOUT) :: u10m, v10m 1282 REAL, DIMENSION(klon,klev+1,nbsrf), INTENT(INOUT) :: tke 1283 1284 ! Local variables 1285 !**************************************************************************************** 1286 INTEGER :: nsrf, nsrf_comp1, nsrf_comp2, nsrf_comp3, i 1287 CHARACTER(len=80) :: abort_message 1288 CHARACTER(len=20) :: modname = 'pbl_surface_newfrac' 1289 INTEGER, DIMENSION(nbsrf) :: nfois=0, mfois=0, pfois=0 1290 LOGICAL :: debug=.FALSE. 1291 ! 1292 ! All at once !! 1293 !**************************************************************************************** 1294 1295 DO nsrf = 1, nbsrf 1296 ! First decide complement sub-surfaces 1297 SELECT CASE (nsrf) 1298 CASE(is_oce) 1299 nsrf_comp1=is_sic 1300 nsrf_comp2=is_ter 1301 nsrf_comp3=is_lic 1302 CASE(is_sic) 1303 nsrf_comp1=is_oce 1304 nsrf_comp2=is_ter 1305 nsrf_comp3=is_lic 1306 CASE(is_ter) 1307 nsrf_comp1=is_lic 1308 nsrf_comp2=is_oce 1309 nsrf_comp3=is_sic 1310 CASE(is_lic) 1311 nsrf_comp1=is_ter 1312 nsrf_comp2=is_oce 1313 nsrf_comp3=is_sic 1314 END SELECT 1315 1316 ! Initialize all new fractions 1317 DO i=1, klon 1318 IF (pctsrf_new(i,nsrf) > 0. .AND. pctsrf_old(i,nsrf) == 0.) THEN 1319 1320 IF (pctsrf_old(i,nsrf_comp1) > 0.) THEN 1321 ! Use the complement sub-surface, keeping the continents unchanged 1322 qsurf(i,nsrf) = qsurf(i,nsrf_comp1) 1323 evap(i,nsrf) = evap(i,nsrf_comp1) 1324 rugos(i,nsrf) = rugos(i,nsrf_comp1) 1325 tsurf(i,nsrf) = tsurf(i,nsrf_comp1) 1326 alb1(i,nsrf) = alb1(i,nsrf_comp1) 1327 alb2(i,nsrf) = alb2(i,nsrf_comp1) 1328 u10m(i,nsrf) = u10m(i,nsrf_comp1) 1329 v10m(i,nsrf) = v10m(i,nsrf_comp1) 1330 tke(i,:,nsrf) = tke(i,:,nsrf_comp1) 1331 mfois(nsrf) = mfois(nsrf) + 1 1332 ELSE 1333 ! The continents have changed. The new fraction receives the mean sum of the existent fractions 1334 qsurf(i,nsrf) = qsurf(i,nsrf_comp2)*pctsrf_old(i,nsrf_comp2) + qsurf(i,nsrf_comp3)*pctsrf_old(i,nsrf_comp3) 1335 evap(i,nsrf) = evap(i,nsrf_comp2) *pctsrf_old(i,nsrf_comp2) + evap(i,nsrf_comp3) *pctsrf_old(i,nsrf_comp3) 1336 rugos(i,nsrf) = rugos(i,nsrf_comp2)*pctsrf_old(i,nsrf_comp2) + rugos(i,nsrf_comp3)*pctsrf_old(i,nsrf_comp3) 1337 tsurf(i,nsrf) = tsurf(i,nsrf_comp2)*pctsrf_old(i,nsrf_comp2) + tsurf(i,nsrf_comp3)*pctsrf_old(i,nsrf_comp3) 1338 alb1(i,nsrf) = alb1(i,nsrf_comp2) *pctsrf_old(i,nsrf_comp2) + alb1(i,nsrf_comp3) *pctsrf_old(i,nsrf_comp3) 1339 alb2(i,nsrf) = alb2(i,nsrf_comp2) *pctsrf_old(i,nsrf_comp2) + alb2(i,nsrf_comp3) *pctsrf_old(i,nsrf_comp3) 1340 u10m(i,nsrf) = u10m(i,nsrf_comp2) *pctsrf_old(i,nsrf_comp2) + u10m(i,nsrf_comp3) *pctsrf_old(i,nsrf_comp3) 1341 v10m(i,nsrf) = v10m(i,nsrf_comp2) *pctsrf_old(i,nsrf_comp2) + v10m(i,nsrf_comp3) *pctsrf_old(i,nsrf_comp3) 1342 tke(i,:,nsrf) = tke(i,:,nsrf_comp2)*pctsrf_old(i,nsrf_comp2) + tke(i,:,nsrf_comp3)*pctsrf_old(i,nsrf_comp3) 1343 1344 ! Security abort. This option has never been tested. To test, comment the following line. 1345 ! abort_message='The fraction of the continents have changed!' 1346 ! CALL abort_gcm(modname,abort_message,1) 1347 nfois(nsrf) = nfois(nsrf) + 1 1348 END IF 1349 snow(i,nsrf) = 0. 1350 agesno(i,nsrf) = 0. 1351 ftsoil(i,:,nsrf) = tsurf(i,nsrf) 1352 ELSE 1353 pfois(nsrf) = pfois(nsrf)+ 1 1354 END IF 1355 END DO 1356 1357 END DO 1358 1359 IF (debug) THEN 1360 print*,'itime=,',itime, 'Pas de nouveau fraction',pfois,'fois' 1361 print*,'itime=,',itime, 'The fraction of the continents have changed',nfois,'fois' 1362 print*,'itime=,',itime, 'The fraction ocean-seaice has changed',mfois,'fois' 1363 END IF 1364 1365 END SUBROUTINE pbl_surface_newfrac 1366 1324 1367 ! 1368 !**************************************************************************************** 1369 ! 1325 1370 1326 1371 END MODULE pbl_surface_mod
Note: See TracChangeset
for help on using the changeset viewer.
