- Timestamp:
- 2011-07-12T19:12:16+02:00 (13 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2011/dev_r2784_CMCC1_topbfm/NEMOGCM/NEMO/OPA_SRC/SBC/sbcblk_core.F90
r2777 r2799 44 44 PUBLIC blk_ice_core ! routine called in sbc_ice_lim module 45 45 46 INTEGER , PARAMETER :: jpfld = 9 ! maximum number of files to read47 46 INTEGER , PARAMETER :: jp_wndi = 1 ! index of 10m wind velocity (i-component) (m/s) at T-point 48 47 INTEGER , PARAMETER :: jp_wndj = 2 ! index of 10m wind velocity (j-component) (m/s) at T-point … … 54 53 INTEGER , PARAMETER :: jp_snow = 8 ! index of snow (solid prcipitation) (kg/m2/s) 55 54 INTEGER , PARAMETER :: jp_tdif = 9 ! index of tau diff associated to HF tau (N/m2) at T-point 55 #if defined key_orca_r025 56 INTEGER , PARAMETER :: jp_swc = 10 ! index of GEWEX correction for SW radiation at T-point 57 INTEGER , PARAMETER :: jp_lwc = 11 ! index of GEWEX correction for LW radiation at T-point 58 INTEGER , PARAMETER :: jp_prc = 12 ! index of PMWC correction forat T-point 59 INTEGER , PARAMETER :: jpfld = 12 ! maximum number of files to read 60 #else 61 INTEGER , PARAMETER :: jpfld = 9 ! maximum number of files to read 62 #endif 56 63 57 64 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf ! structure of input fields (file informations, fields read) … … 71 78 REAL(wp) :: rn_pfac = 1. ! multiplication factor for precipitation 72 79 80 #if defined key_orca_r025 81 LOGICAL :: ln_printdia= .TRUE. ! logical flag for height of air temp. and hum 82 LOGICAL :: ln_netsw = .TRUE. ! logical flag for height of air temp. and hum 83 LOGICAL :: ln_core_graceopt=.FALSE., ln_core_spinup=.FALSE. 84 LOGICAL :: ln_gwxc = .TRUE. 85 LOGICAL :: ln_corad_antar =.FALSE., ln_corad_arc =.FALSE. , ln_cotair_arc = .FALSE. 86 LOGICAL :: ln_coprecip =.FALSE. 87 REAL(wp) :: rn_qns_bias = 0._wp ! heat flux bias 88 89 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: area 90 REAL(wp) :: araux 91 92 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: zqlw, zqsb ! long wave and sensible heat fluxes 93 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: zqla, zevap ! latent heat fluxes and evaporation 94 95 REAL(wp), PARAMETER :: zalph = 2.408724e-06_wp, & 96 & zbet = -0.006936579_wp, zgam = 449.9094_wp ! GRACE regression coefficients 97 #endif 98 73 99 !! * Substitutions 74 100 # include "domzgr_substitute.h90" … … 80 106 !!---------------------------------------------------------------------- 81 107 CONTAINS 82 108 83 109 SUBROUTINE sbc_blk_core( kt ) 84 110 !!--------------------------------------------------------------------- … … 112 138 !! - emp, emps evaporation minus precipitation 113 139 !!---------------------------------------------------------------------- 140 #if defined key_orca_r025 141 USE ice_2 142 #endif 114 143 INTEGER, INTENT(in) :: kt ! ocean time step 115 144 !! … … 117 146 INTEGER :: ifpr ! dummy loop indice 118 147 INTEGER :: jfld ! dummy loop arguments 148 INTEGER :: ji, jj 119 149 !! 120 150 CHARACTER(len=100) :: cn_dir ! Root directory for location of core files … … 123 153 TYPE(FLD_N) :: sn_qlw , sn_tair, sn_prec, sn_snow ! " " 124 154 TYPE(FLD_N) :: sn_tdif ! " " 155 TYPE(FLD_N) :: sn_swc, sn_lwc ! " " 156 TYPE(FLD_N) :: sn_prc 157 INTEGER :: iter_shapiro = 250 158 REAL :: zzlat, zzlat1, zzlat2, zfm, zfrld, ztmp 159 REAL(wp), DIMENSION(jpi,jpj):: xyt,z_qsr,z_qlw,z_qsr1,z_qlw1, z_hum, z_tair 160 REAL(wp), DIMENSION(jpi,jpj):: zqsr_lr, zqsr_hr, zqlw_lr, zqlw_hr, zprec_hr, zprec_lr 161 CHARACTER(len=20) :: c_kind='ORCA_GLOB' 162 #if defined key_orca_r025 163 NAMELIST/namsbc_core/ cn_dir , ln_2m , ln_taudif, rn_pfac, & 164 & sn_wndi, sn_wndj, sn_humi , sn_qsr , & 165 & sn_qlw , sn_tair, sn_prec , sn_snow, sn_tdif & 166 & sn_swc , sn_lwc , sn_prc , ln_gwxc , & 167 & ln_corad_antar, ln_corad_arc, ln_cotair_arc, ln_coprecip , & 168 & rn_qns_bias, ln_printdia, ln_netsw, ln_core_graceopt,ln_core_spinup 169 !!--------------------------------------------------------------------- 170 #else 125 171 NAMELIST/namsbc_core/ cn_dir , ln_2m , ln_taudif, rn_pfac, & 126 172 & sn_wndi, sn_wndj, sn_humi , sn_qsr , & 127 173 & sn_qlw , sn_tair, sn_prec , sn_snow, sn_tdif 128 174 !!--------------------------------------------------------------------- 175 #endif 129 176 130 177 ! ! ====================== ! 131 178 IF( kt == nit000 ) THEN ! First call kt=nit000 ! 132 179 ! ! ====================== ! 180 #if defined key_orca_r025 181 ! !== allocate sbc arrays 182 IF( sbc_blk_core_alloc() /= 0 ) CALL ctl_stop( 'STOP', 'sbc_blk_core_alloc : unable to allocate arrays' ) 183 #endif 184 133 185 ! set file information (default values) 134 186 cn_dir = './' ! directory in which the model is executed … … 146 198 sn_snow = FLD_N( 'snow' , -1 , 'snow' , .true. , .false. , 'yearly' , '' , '' ) 147 199 sn_tdif = FLD_N( 'taudif' , 24 , 'taudif' , .true. , .false. , 'yearly' , '' , '' ) 200 #if defined key_orca_r025 201 sn_swc = FLD_N( 'swc' , 24 , 'swc' , .true. , .false. , 'yearly' , '' , '' ) 202 sn_lwc = FLD_N( 'lwc' , 24 , 'lwc' , .true. , .false. , 'yearly' , '' , '' ) 203 sn_prc = FLD_N( 'prc' , 24 , 'prc' , .true. , .false. , 'yearly' , '' , '' ) 204 #endif 148 205 ! 149 206 REWIND( numnam ) ! read in namlist namsbc_core … … 163 220 slf_i(jp_prec) = sn_prec ; slf_i(jp_snow) = sn_snow 164 221 slf_i(jp_tdif) = sn_tdif 222 #if defined key_orca_r025 223 slf_i(jp_swc ) = sn_swc 224 slf_i(jp_lwc ) = sn_lwc 225 slf_i(jp_prc ) = sn_prc 226 #endif 165 227 ! 166 228 lhftau = ln_taudif ! do we use HF tau information? 167 229 jfld = jpfld - COUNT( (/.NOT. lhftau/) ) 230 #if defined key_orca_r025 231 IF( .NOT. ln_gwxc ) jfld = jfld - 2 232 IF( .NOT. ln_coprecip ) jfld = jfld - 1 233 #endif 168 234 ! 169 235 ALLOCATE( sf(jfld), STAT=ierror ) ! set sf structure … … 176 242 CALL fld_fill( sf, slf_i, cn_dir, 'sbc_blk_core', 'flux formulation for ocean surface boundary condition', 'namsbc_core' ) 177 243 ! 244 #if defined key_orca_r025 245 IF( ln_printdia .OR. ln_core_graceopt ) THEN 246 area = (e1t * e2t) 247 araux = sum ( area * tmask(:,:,1) ) 248 IF(lk_mpp) CALL mpp_sum ( araux ) 249 ENDIF 250 #endif 178 251 ENDIF 179 252 180 253 CALL fld_read( kt, nn_fsbc, sf ) ! input fields provided at the current time-step 181 254 182 ! ! surface ocean fluxes computed with CLIO bulk formulea 183 IF( MOD( kt - 1, nn_fsbc ) == 0 ) CALL blk_oce_core( sf, sst_m, ssu_m, ssv_m ) 184 ! 255 IF( MOD( kt - 1, nn_fsbc ) == 0 ) THEN 256 257 #if defined key_orca_r025 258 ! Introduce ERA-Interim filtering and correction 259 260 IF( ln_gwxc ) THEN 261 262 call Shapiro_1D(sf(jp_qsr)%fnow(:,:,1),iter_shapiro, c_kind, zqsr_lr) 263 zqsr_hr(:,:)=sf(jp_qsr)%fnow(:,:,1)-zqsr_lr(:,:) ! We get large scale and small scale 264 265 call Shapiro_1D(sf(jp_qlw)%fnow(:,:,1),iter_shapiro, c_kind, zqlw_lr) 266 zqlw_hr(:,:)=sf(jp_qlw)%fnow(:,:,1)-zqlw_lr(:,:) ! We get large scale and small scale 267 268 z_qsr1(:,:)=zqsr_lr(:,:)*sf(jp_swc)%fnow(:,:,1) + zqsr_hr(:,:) 269 z_qlw1(:,:)=zqlw_lr(:,:)*sf(jp_lwc)%fnow(:,:,1) + zqlw_hr(:,:) 270 271 DO jj=1,jpj 272 DO ji=1,jpi 273 z_qsr1(ji,jj)=max(z_qsr1(ji,jj),0.0) 274 z_qlw1(ji,jj)=max(z_qlw1(ji,jj),0.0) 275 END DO 276 END DO 277 278 ENDIF 279 280 IF( ln_coprecip ) THEN 281 282 call Shapiro_1D(sf(jp_prec)%fnow(:,:,1),iter_shapiro,c_kind,zprec_lr) 283 zprec_hr(:,:)=sf(jp_prec)%fnow(:,:,1)-zprec_lr(:,:) ! We get large scale and small scale 284 285 DO jj=1,jpj 286 DO ji=1,jpi 287 IF( zprec_lr(ji,jj) .GT. 0._wp ) THEN 288 ztmp = LOG( ( 1000._wp + sf(jp_prc)%fnow(ji,jj,1) ) * EXP( zprec_lr(ji,jj) ) / 1000._wp ) 289 sf(jp_prec)%fnow(ji,jj,1) = max(ztmp+zprec_hr(ji,jj),0.0) 290 ENDIF 291 END DO 292 END DO 293 294 ENDIF 295 296 IF ( ln_corad_antar ) THEN ! correction of SW and LW in the Southern Ocean 297 298 z_qsr(:,:)=0.8*z_qsr1(:,:) 299 z_qlw(:,:)=1.1*z_qlw1(:,:) 300 xyt(:,:) = 0.e0 301 zzlat1 = -65. 302 zzlat2 = -60. 303 DO jj = 1, jpj 304 DO ji = 1, jpi 305 zzlat = gphit(ji,jj) 306 IF ( zzlat >= zzlat1 .AND. zzlat <= zzlat2 ) THEN 307 xyt(ji,jj) = (zzlat2-zzlat)/(zzlat2-zzlat1) 308 ELSE IF ( zzlat < zzlat1 ) THEN 309 xyt(ji,jj) = 1 310 ENDIF 311 END DO 312 END DO 313 z_qsr1(:,:)=z_qsr(:,:)*xyt(:,:)+(1.0-xyt(:,:))*z_qsr1(:,:) 314 z_qlw1(:,:)=z_qlw(:,:)*xyt(:,:)+(1.0-xyt(:,:))*z_qlw1(:,:) 315 316 ENDIF 317 318 IF ( ln_corad_arc ) THEN ! correction of SW in the Arctic Ocean 319 320 z_qsr(:,:)=0.7*z_qsr1(:,:) 321 xyt(:,:) = 0.e0 322 zzlat1 = 78. 323 zzlat2 = 82. 324 DO jj = 1, jpj 325 DO ji = 1, jpi 326 zzlat = gphit(ji,jj) 327 IF ( zzlat >= zzlat1 .AND. zzlat <= zzlat2 ) THEN 328 xyt(ji,jj) = (zzlat-zzlat1)/(zzlat2-zzlat1) 329 ELSE IF ( zzlat > zzlat2 ) THEN 330 xyt(ji,jj) = 1 331 ENDIF 332 END DO 333 END DO 334 z_qsr1(:,:)=z_qsr(:,:)*xyt(:,:)+(1.0-xyt(:,:))*z_qsr1(:,:) 335 336 ENDIF 337 338 sf(jp_qsr)%fnow(:,:,1)=z_qsr1(:,:) 339 sf(jp_qlw)%fnow(:,:,1)=z_qlw1(:,:) 340 341 IF ( ln_cotair_arc ) THEN ! correction of Air Temperature in the Arctic Ocean 342 343 z_tair(:,:)=sf(jp_tair)%fnow(:,:,1) - 2.0 344 xyt(:,:) = 0.e0 ; zzlat1 = 78. ; zzlat2 = 82. 345 DO jj = 1, jpj 346 DO ji = 1, jpi 347 zzlat = gphit(ji,jj) ; zfrld=frld(ji,jj) 348 IF ( zzlat >= zzlat1 .AND. zzlat <= zzlat2 .AND. zfrld < 0.85 ) THEN 349 xyt(ji,jj) = (zzlat-zzlat1)/(zzlat2-zzlat1) 350 ELSE IF ( zzlat > zzlat2 .AND. zfrld < 0.85 ) THEN 351 xyt(ji,jj) = 1 352 ENDIF 353 END DO 354 END DO 355 sf(jp_tair)%fnow(:,:,1)=z_tair(:,:)*xyt(:,:)+(1.0-xyt(:,:))*sf(jp_tair)%fnow(:,:,1) 356 357 ENDIF 358 359 #endif 360 361 CALL blk_oce_core( sf, sst_m, ssu_m, ssv_m ) ! compute the surface ocean fluxes using CLIO bulk formulea 362 363 ENDIF 364 ! ! using CORE bulk formulea 185 365 END SUBROUTINE sbc_blk_core 186 366 … … 315 495 IF( lhftau ) THEN 316 496 !CDIR COLLAPSE 497 #if defined key_orca_r025 498 ! Changed!!! Multiply by QSCAT correction 499 zwnd_i(:,:) = zwnd_i(:,:) * sf(jp_tdif)%fnow(:,:,1) 500 zwnd_j(:,:) = zwnd_j(:,:) * sf(jp_tdif)%fnow(:,:,1) 501 #endif 317 502 taum(:,:) = taum(:,:) + sf(jp_tdif)%fnow(:,:,1) 318 503 ENDIF … … 938 1123 ! 939 1124 END FUNCTION psi_h 1125 1126 #if defined key_orca_r025 1127 INTEGER FUNCTION sbc_blk_core_alloc() 1128 !!---------------------------------------------------------------------- 1129 !! *** ROUTINE sbc_blk_core_alloc *** 1130 !!---------------------------------------------------------------------- 1131 ALLOCATE( area(jpi,jpj) , zqlw(jpi,jpj) , & 1132 & zqsb(jpi,jpj) , zqla(jpi,jpj) , & 1133 & zevap(jpi,jpj) , STAT=sbc_blk_core_alloc ) 1134 zqlw=0._wp 1135 zqsb=0._wp 1136 ! 1137 IF( lk_mpp ) CALL mpp_sum ( sbc_blk_core_alloc ) 1138 IF( sbc_blk_core_alloc > 0 ) CALL ctl_warn('sbc_blk_core_alloc: allocation of arrays failed') 1139 END FUNCTION sbc_blk_core_alloc 1140 #endif 1141 1142 SUBROUTINE Shapiro_1D(rla_varin,id_np, cd_overlap, rlpa_varout) !GIG 1143 !!===================================================================== 1144 !! 1145 !! Description: This function applies a 1D Shapiro filter 1146 !! (3 points filter) horizontally to a 2D field 1147 !! in regular grid 1148 !! Arguments : 1149 !! rla_varin : Input variable to filter 1150 !! zla_mask : Input mask variable 1151 !! id_np : Number of Shapiro filter iterations 1152 !! cd_overlap : Logical argument for periodical condition 1153 !! (global ocean case) 1154 !! rlpa_varout : Output filtered variable 1155 !! 1156 !! History : 08/2009 S. CAILLEAU : from 1st version of N. FERRY 1157 !! 09/2009 C. REGNIER : Corrections 1158 !! 1159 !!===================================================================== 1160 IMPLICIT NONE 1161 INTEGER, INTENT(IN) :: id_np 1162 REAL(wp), DIMENSION(jpi,jpj), INTENT(IN) :: rla_varin !GIG 1163 CHARACTER(len=20), INTENT(IN) :: cd_overlap !GIG 1164 REAL(wp), DIMENSION(jpi,jpj), INTENT(OUT) :: rlpa_varout !GIG 1165 1166 REAL(wp), DIMENSION(jpi,jpj) :: rlpa_varout_tmp 1167 REAL, PARAMETER :: rl_alpha = 1./2. ! fixed stability coefficient (isotrope case) 1168 REAL, parameter :: rap_aniso_diff_XY=2.25 ! anisotrope case 1169 REAL :: alphax,alphay, znum, zden,test 1170 INTEGER :: ji, jj, jn, nn 1171 ! 1172 !! rap_aniso_diff_XY=2.25 : valeur trouvée empiriquement pour 140 itération pour le filtre de Shapiro et 1173 !! pour un rapport d'anisotopie de 1.5 : on filtre de plus rapidement en x qu'eny. 1174 !------------------------------------------------------------------------------ 1175 ! 1176 ! Loop on several filter iterations 1177 1178 ! Global ocean case 1179 IF (( cd_overlap == 'MERCA_GLOB' ) .OR. & 1180 ( cd_overlap == 'REGULAR_GLOB' ) .OR. & 1181 ( cd_overlap == 'ORCA_GLOB' )) THEN 1182 rlpa_varout(:,:) = rla_varin(:,:) 1183 rlpa_varout_tmp(:,:) = rlpa_varout(:,:) 1184 ! 1185 1186 alphax=1./2. 1187 alphay=1./2. 1188 ! Dx/Dy=rap_aniso_diff_XY , D_ = vitesse de diffusion 1189 ! 140 passes du fitre, Lx/Ly=1.5, le rap_aniso_diff_XY correspondant est: 1190 IF ( rap_aniso_diff_XY .GE. 1. ) alphay=alphay/rap_aniso_diff_XY 1191 IF ( rap_aniso_diff_XY .LT. 1. ) alphax=alphax*rap_aniso_diff_XY 1192 1193 DO jn = 1,id_np ! number of passes of the filter 1194 DO ji = 2,jpim1 1195 DO jj = 2,jpjm1 1196 ! We crop on the coast 1197 znum = rlpa_varout_tmp(ji,jj) & 1198 + 0.25*alphax*(rlpa_varout_tmp(ji-1,jj )-rlpa_varout_tmp(ji,jj))*tmask(ji-1,jj ,1) & 1199 + 0.25*alphax*(rlpa_varout_tmp(ji+1,jj )-rlpa_varout_tmp(ji,jj))*tmask(ji+1,jj ,1) & 1200 + 0.25*alphay*(rlpa_varout_tmp(ji ,jj-1)-rlpa_varout_tmp(ji,jj))*tmask(ji ,jj-1,1) & 1201 + 0.25*alphay*(rlpa_varout_tmp(ji ,jj+1)-rlpa_varout_tmp(ji,jj))*tmask(ji ,jj+1,1) 1202 rlpa_varout(ji,jj)=znum*tmask(ji,jj,1)+rla_varin(ji,jj)*(1.-tmask(ji,jj,1)) 1203 ENDDO ! end loop ji 1204 ENDDO ! end loop jj 1205 ! 1206 ! 1207 ! Periodical condition in case of cd_overlap (global ocean) 1208 ! - on a mercator projection grid we consider that singular point at poles 1209 ! are a mean of the values at points of the previous latitude 1210 ! - on ORCA and regular grid we copy the values at points of the previous latitude 1211 IF ( cd_overlap == 'MERCAT_GLOB' ) THEN 1212 !GIG case unchecked 1213 rlpa_varout(1,1) = SUM(rlpa_varout(:,2)) / jpi 1214 rlpa_varout(jpi,jpj) = SUM(rlpa_varout(:,jpj-1)) / jpi 1215 ELSE 1216 call lbc_lnk(rlpa_varout, 'T', 1.) ! Boundary condition 1217 ENDIF 1218 rlpa_varout_tmp(:,:) = rlpa_varout(:,:) 1219 ENDDO ! end loop jn 1220 ENDIF 1221 1222 ! 1223 END SUBROUTINE Shapiro_1D 1224 940 1225 941 1226 !!======================================================================
Note: See TracChangeset
for help on using the changeset viewer.
