Changeset 1975 for branches/DEV_r1837_MLF/NEMO/OPA_SRC/TRA
- Timestamp:
- 2010-06-28T19:22:14+02:00 (14 years ago)
- Location:
- branches/DEV_r1837_MLF/NEMO/OPA_SRC/TRA
- Files:
-
- 3 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/DEV_r1837_MLF/NEMO/OPA_SRC/TRA/tranxt.F90
r1870 r1975 25 25 USE oce ! ocean dynamics and tracers variables 26 26 USE dom_oce ! ocean space and time domain variables 27 USE sbc_oce ! surface boundary condition: ocean 27 28 USE zdf_oce ! ??? 28 29 USE domvvl ! variable volume … … 37 38 USE lbclnk ! ocean lateral boundary conditions (or mpp link) 38 39 USE prtctl ! Print control 40 USE traqsr ! penetrative solar radiation (needed for nksr) 39 41 USE agrif_opa_update 40 42 USE agrif_opa_interp … … 45 47 PUBLIC tra_nxt ! routine called by step.F90 46 48 47 REAL(wp) :: rbcp , r1_2bcp! Brown & Campana parameters for semi-implicit hpg49 REAL(wp) :: rbcp ! Brown & Campana parameters for semi-implicit hpg 48 50 REAL(wp), DIMENSION(jpk) :: r2dt_t ! vertical profile time step, =2*rdttra (leapfrog) or =rdttra (Euler) 49 51 … … 98 100 IF(lwp) WRITE(numout,*) '~~~~~~~' 99 101 ! 100 rbcp = 0.25 * (1 + atfp) * (1 + atfp * atfp) ! Brown & Campana parameter for semi-implicit hpg 101 r1_2bcp = 1.e0 - 2.e0 * rbcp 102 rbcp = 0.25 * (1. + atfp) * (1. + atfp) * ( 1. - atfp) ! Brown & Campana parameter for semi-implicit hpg 102 103 ENDIF 103 104 ! ! set time step size (Euler/Leapfrog) … … 164 165 !! - swap tracer fields to prepare the next time_step. 165 166 !! This can be summurized for temperature as: 166 !! ztm = (rbcp*ta+(2-rbcp)*tn+rbcp*tb)ln_dynhpg_imp = T167 !! ztm = tn + rbcp * [ta -2 tn + tb ] ln_dynhpg_imp = T 167 168 !! ztm = 0 otherwise 168 !! with rbcp= (1+atfp)*(1+atfp*atfp)/4169 !! with rbcp=1/4 * (1-atfp^4) / (1-atfp) 169 170 !! tb = tn + atfp*[ tb - 2 tn + ta ] 170 171 !! tn = ta … … 177 178 !! 178 179 INTEGER :: ji, jj, jk ! dummy loop indices 179 REAL(wp) :: zt m, ztf, zac! temporary scalars180 REAL(wp) :: z sm, zsf! - -180 REAL(wp) :: zt_m, zs_m ! temporary scalars 181 REAL(wp) :: ztn, zsn ! - - 181 182 !!---------------------------------------------------------------------- 182 183 … … 187 188 ENDIF 188 189 ! 189 ! ! ----------------------- ! 190 IF( ln_dynhpg_imp ) THEN ! semi-implicite hpg case ! 191 ! ! ----------------------- ! 192 ! 193 IF( neuler == 0 .AND. kt == nit000 ) THEN ! Euler time-stepping at first time-step 194 DO jk = 1, jpkm1 ! (only swap) 195 tn(:,:,jk) = ta(:,:,jk) ! tb <-- tn 196 sn(:,:,jk) = sa(:,:,jk) 197 END DO 198 ELSE ! general case (Leapfrog + Asselin filter) 199 DO jk = 1, jpkm1 200 DO jj = 1, jpj 201 DO ji = 1, jpi 202 ztm = rbcp * ( ta(ji,jj,jk) + tb(ji,jj,jk) ) + r1_2bcp * tn(ji,jj,jk) ! mean t 203 zsm = rbcp * ( sa(ji,jj,jk) + sb(ji,jj,jk) ) + r1_2bcp * sn(ji,jj,jk) 204 ztf = atfp * ( ta(ji,jj,jk) + tb(ji,jj,jk) - 2. * tn(ji,jj,jk) ) ! Asselin filter on t 205 zsf = atfp * ( sa(ji,jj,jk) + sb(ji,jj,jk) - 2. * sn(ji,jj,jk) ) 206 tb(ji,jj,jk) = tn(ji,jj,jk) + ztf ! tb <-- filtered tn 207 sb(ji,jj,jk) = sn(ji,jj,jk) + zsf 208 tn(ji,jj,jk) = ta(ji,jj,jk) ! tn <-- ta 209 sn(ji,jj,jk) = sa(ji,jj,jk) 210 ta(ji,jj,jk) = ztm ! ta <-- mean t 211 sa(ji,jj,jk) = zsm 212 END DO 190 IF( neuler == 0 .AND. kt == nit000 ) THEN ! Euler time-stepping at first time-step 191 DO jk = 1, jpkm1 ! (only swap) 192 tn(:,:,jk) = ta(:,:,jk) ! tn <-- ta 193 sn(:,:,jk) = sa(:,:,jk) ! sn <-- sa 194 END DO 195 ELSE ! General case (Leapfrog + Asselin filter) 196 DO jk = 1, jpkm1 197 DO jj = 1, jpj 198 DO ji = 1, jpi 199 IF( ln_dynhpg_imp ) THEN ! implicit hpg: keep tn, sn in memory 200 ztn = tn(ji,jj,jk) 201 zsn = sn(ji,jj,jk) 202 ENDIF 203 ! ! time laplacian on tracers 204 ! ! used for both Asselin and Brown & Campana filters 205 zt_m = ta(ji,jj,jk) - 2. * tn(ji,jj,jk) + tb(ji,jj,jk) 206 zs_m = sa(ji,jj,jk) - 2. * sn(ji,jj,jk) + sb(ji,jj,jk) 207 ! 208 ! ! swap of arrays 209 tb(ji,jj,jk) = tn(ji,jj,jk) + atfp * zt_m ! tb <-- tn filtered 210 sb(ji,jj,jk) = sn(ji,jj,jk) + atfp * zs_m ! sb <-- sn filtered 211 tn(ji,jj,jk) = ta(ji,jj,jk) ! tn <-- ta 212 sn(ji,jj,jk) = sa(ji,jj,jk) ! sn <-- sa 213 ! ! semi imlicit hpg computation (Brown & Campana) 214 IF( ln_dynhpg_imp ) THEN 215 ta(ji,jj,jk) = ztn + rbcp * zt_m ! ta <-- Brown & Campana average 216 sa(ji,jj,jk) = zsn + rbcp * zs_m ! sa <-- Brown & Campana average 217 ENDIF 213 218 END DO 214 219 END DO 215 ENDIF 216 ! ! ----------------------- ! 217 ELSE ! explicit hpg case ! 218 ! ! ----------------------- ! 219 ! 220 IF( neuler == 0 .AND. kt == nit000 ) THEN ! Euler time-stepping at first time-step 221 DO jk = 1, jpkm1 222 tn(:,:,jk) = ta(:,:,jk) ! tn <-- ta 223 sn(:,:,jk) = sa(:,:,jk) 224 END DO 225 ELSE ! general case (Leapfrog + Asselin filter) 226 DO jk = 1, jpkm1 227 DO jj = 1, jpj 228 DO ji = 1, jpi 229 ztf = atfp * ( ta(ji,jj,jk) - 2.* tn(ji,jj,jk) + tb(ji,jj,jk) ) ! Asselin filter on t 230 zsf = atfp * ( sa(ji,jj,jk) - 2.* sn(ji,jj,jk) + sb(ji,jj,jk) ) 231 tb(ji,jj,jk) = tn(ji,jj,jk) + ztf ! tb <-- filtered tn 232 sb(ji,jj,jk) = sn(ji,jj,jk) + zsf 233 tn(ji,jj,jk) = ta(ji,jj,jk) ! tn <-- ta 234 sn(ji,jj,jk) = sa(ji,jj,jk) 235 END DO 236 END DO 237 END DO 238 ENDIF 239 ENDIF 240 ! 241 !!gm from Matthieu : unchecked 242 IF( neuler /= 0 .OR. kt /= nit000 ) THEN ! remove the forcing from the Asselin filter 243 zac = atfp * rdttra(1) 244 tb(:,:,1) = tb(:,:,1) - zac * ( qns (:,:) - qns_b (:,:) ) ! non solar surface flux 245 sb(:,:,1) = sn(:,:,1) - zac * ( emps(:,:) - emps_b(:,:) ) ! surface salt flux 246 ! 247 IF( ln_traqsr ) ! solar penetrating flux 248 DO jk = 1, nksr 249 DO jj = 1, jpj 250 DO ji = 1, jpi 251 IF( ln_sco ) THEN 252 z_cor_qsr = etot3(ji,jj,jk) * ( qsr(ji,jj) - qsrb(ji,jj) ) 253 ELSEIF( ln_zps .OR. ln_zco ) THEN 254 z_cor_qsr = ( qsr(ji,jj) - qsrb(ji,jj) ) * & 255 & ( gdsr(jk) * tmask(ji,jj,jk) - gdsr(jk+1) * tmask(ji,jj,jk+1) ) 256 ENDIF 257 zt = zt - zfact1 * z_cor_qsr 258 END DO 259 END DO 220 END DO 260 221 ENDIF 261 222 ! … … 276 237 !! - swap tracer fields to prepare the next time_step. 277 238 !! This can be summurized for temperature as: 278 !! ztm = ( rbcp*e3t_a*ta + (2-rbtp)*e3t_n*tn + rbcp*e3t_b*tb) ln_dynhpg_imp = T279 !! /( rbcp*e3t_a + (2-rbcp)*e3t_n + rbcp*e3t_b)239 !! ztm = ( e3t_n*tn + rbcp*[ e3t_b*tb - 2 e3t_n*tn + e3t_a*ta ] ) ln_dynhpg_imp = T 240 !! /( e3t_n + rbcp*[ e3t_b - 2 e3t_n + e3t_a ] ) 280 241 !! ztm = 0 otherwise 281 242 !! tb = ( e3t_n*tn + atfp*[ e3t_b*tb - 2 e3t_n*tn + e3t_a*ta ] ) … … 287 248 !! - (ta,sa) time averaged (t,s) (ln_dynhpg_imp = T) 288 249 !!---------------------------------------------------------------------- 289 INTEGER, INTENT(in) :: kt ! ocean time-step index250 INTEGER, INTENT(in) :: kt ! ocean time-step index 290 251 !! 291 INTEGER :: ji, jj, jk ! dummy loop indices 292 REAL(wp) :: ztm , ztc_f , ztf , ztca, ztcn, ztcb ! temporary scalar 293 REAL(wp) :: zsm , zsc_f , zsf , zsca, zscn, zscb ! - - 294 REAL(wp) :: ze3mr, ze3fr ! - - 295 REAL(wp) :: ze3t_b, ze3t_n, ze3t_a, ze3t_f ! - - 296 !!---------------------------------------------------------------------- 297 252 INTEGER :: ji, jj, jk ! dummy loop indices 253 REAL :: ze3t_a, ze3t_n, ze3t_b ! temporary scalars 254 REAL :: ztc_a, ztc_n, ztc_b ! - - 255 REAL :: zsc_a, zsc_n, zsc_b ! - - 256 REAL :: ztc_f, zsc_f, ztc_m, zsc_m ! - - 257 REAL :: ze3t_f, ze3t_m ! - - 258 REAL :: zfact1, zfact2 ! - - 259 !!---------------------------------------------------------------------- 260 298 261 IF( kt == nit000 ) THEN 299 262 IF(lwp) WRITE(numout,*) … … 302 265 ENDIF 303 266 304 ! ! ----------------------- ! 305 IF( ln_dynhpg_imp ) THEN ! semi-implicite hpg case ! (mean tracer computation) 306 ! ! ----------------------- ! 307 ! 308 IF( neuler == 0 .AND. kt == nit000 ) THEN ! Euler time-stepping at first time-step 309 DO jk = 1, jpkm1 ! (only swap) 310 tn(:,:,jk) = ta(:,:,jk) ! tn <-- ta 311 sn(:,:,jk) = sa(:,:,jk) 312 END DO 313 ! ! general case (Leapfrog + Asselin filter) 314 ELSE ! apply filter on thickness weighted tracer and swap 315 DO jk = 1, jpkm1 316 DO jj = 1, jpj 317 DO ji = 1, jpi 318 ze3t_b = fse3t_b(ji,jj,jk) 319 ze3t_n = fse3t_n(ji,jj,jk) 320 ze3t_a = fse3t_a(ji,jj,jk) 321 ! ! tracer content at Before, now and after 322 ztcb = tb(ji,jj,jk) * ze3t_b ; zscb = sb(ji,jj,jk) * ze3t_b 323 ztcn = tn(ji,jj,jk) * ze3t_n ; zscn = sn(ji,jj,jk) * ze3t_n 324 ztca = ta(ji,jj,jk) * ze3t_a ; zsca = sa(ji,jj,jk) * ze3t_a 325 ! 326 ! ! Asselin filter on thickness and tracer content 327 ze3t_f = atfp * ( ze3t_a + ze3t_b - 2.* ze3t_n ) 328 ztc_f = atfp * ( ztca + ztcb - 2.* ztcn ) 329 zsc_f = atfp * ( zsca + zscb - 2.* zscn ) 330 ! 331 ! ! filtered tracer including the correction 332 ze3fr = 1.e0 / ( ze3t_n + ze3t_f ) 333 ztf = ze3fr * ( ztcn + ztc_f ) 334 zsf = ze3fr * ( zscn + zsc_f ) 335 ! ! mean thickness and tracer 336 ze3mr = 1.e0 / ( rbcp * ( ze3t_a + ze3t_b ) + r1_2bcp * ze3t_n ) 337 ztm = ze3mr * ( rbcp * ( ztca + ztcb ) + r1_2bcp * ztcn ) 338 zsm = ze3mr * ( rbcp * ( zsca + zscb ) + r1_2bcp * zscn ) 339 !!gm mean e3t have to be saved and used in dynhpg or it can be recomputed in dynhpg !! 340 !!gm e3t_m(ji,jj,jk) = 1.e0 / ze3mr 341 ! ! swap of arrays 342 tb(ji,jj,jk) = ztf ! tb <-- tn + filter 343 sb(ji,jj,jk) = zsf 344 tn(ji,jj,jk) = ta(ji,jj,jk) ! tn <-- ta 345 sn(ji,jj,jk) = sa(ji,jj,jk) 346 ta(ji,jj,jk) = ztm ! ta <-- mean t 347 sa(ji,jj,jk) = zsm 348 END DO 267 IF( neuler == 0 .AND. kt == nit000 ) THEN ! Euler time-stepping at first time-step 268 ! ! (only swap) 269 DO jk = 1, jpkm1 ! tn <-- ta 270 tn(:,:,jk) = ta(:,:,jk) ! sn <-- sa 271 sn(:,:,jk) = sa(:,:,jk) 272 END DO 273 ! ! General case (Leapfrog + Asselin filter) 274 ELSE ! apply filter on thickness weighted tracer and swap 275 DO jk = 1, jpkm1 276 zfact1 = atfp * r2dt_t(jk) 277 zfact2 = zfact1 / rau0 278 DO jj = 1, jpj 279 DO ji = 1, jpi 280 ! ! scale factors at Before, now and after 281 ze3t_b = fse3t_b(ji,jj,jk) 282 ze3t_n = fse3t_n(ji,jj,jk) 283 ze3t_a = fse3t_a(ji,jj,jk) 284 ! ! tracer content at Before, now and after 285 ztc_b = tb(ji,jj,jk) * ze3t_b ; zsc_b = sb(ji,jj,jk) * ze3t_b 286 ztc_n = tn(ji,jj,jk) * ze3t_n ; zsc_n = sn(ji,jj,jk) * ze3t_n 287 ztc_a = ta(ji,jj,jk) * ze3t_a ; zsc_a = sa(ji,jj,jk) * ze3t_a 288 ! 289 ! ! Time laplacian on thickness and tracer content 290 ! ! used for both Asselin and Brown & Campana filters 291 ze3t_m = ze3t_a - 2. * ze3t_n + ze3t_b 292 ztc_m = ztc_a - 2. * ztc_n + ztc_b 293 zsc_m = zsc_a - 2. * zsc_n + zsc_b 294 ! ! Asselin Filter + correction 295 ze3t_f = ze3t_n + atfp * ze3t_m 296 ztc_f = ztc_n + atfp * ztc_m 297 zsc_f = zsc_n + atfp * zsc_m 298 299 IF( jk == 1 ) THEN 300 ze3t_f = ze3t_f - zfact2 * ( emp_b (ji,jj) - emp (ji,jj) ) 301 ztc_f = ztc_f - zfact1 * ( sbc_trd_hc_n(ji,jj) - sbc_trd_hc_b(ji,jj) ) 302 ENDIF 303 IF( ln_traqsr .AND. ( jk .LE. nksr ) ) THEN 304 ztc_f = ztc_f - zfact1 * ( qsr_trd_hc_n(ji,jj,jk) - qsr_trd_hc_b(ji,jj,jk) ) 305 ENDIF 306 ! ! swap of arrays 307 ze3t_f = 1.e0 / ze3t_f 308 tb(ji,jj,jk) = ztc_f * ze3t_f ! tb <-- tn filtered 309 sb(ji,jj,jk) = zsc_f * ze3t_f ! sb <-- sn filtered 310 tn(ji,jj,jk) = ta(ji,jj,jk) ! tn <-- ta 311 sn(ji,jj,jk) = sa(ji,jj,jk) ! sn <-- sa 312 ! ! semi imlicit hpg computation (Brown & Campana) 313 IF( ln_dynhpg_imp ) THEN 314 ze3t_m = 1.e0 / ( ze3t_n + rbcp * ze3t_m ) 315 ta(ji,jj,jk) = ( ztc_n + rbcp * ztc_m ) * ze3t_m ! ta <-- Brown & Campana average 316 sa(ji,jj,jk) = ( zsc_n + rbcp * zsc_m ) * ze3t_m ! sa <-- Brown & Campana average 317 ENDIF 349 318 END DO 350 319 END DO 351 ENDIF 352 ! ! ----------------------- ! 353 ELSE ! explicit hpg case ! (NO mean tracer computation) 354 ! ! ----------------------- ! 355 ! 356 IF( neuler == 0 .AND. kt == nit000 ) THEN ! case of Euler time-stepping at first time-step 357 DO jk = 1, jpkm1 ! (only swap) 358 tn(:,:,jk) = ta(:,:,jk) ! tn <-- ta 359 sn(:,:,jk) = sa(:,:,jk) 360 END DO 361 ! ! general case (Leapfrog + Asselin filter) 362 ELSE ! apply filter on thickness weighted tracer and swap 363 DO jk = 1, jpkm1 364 DO jj = 1, jpj 365 DO ji = 1, jpi 366 ze3t_b = fse3t_b(ji,jj,jk) 367 ze3t_n = fse3t_n(ji,jj,jk) 368 ze3t_a = fse3t_a(ji,jj,jk) 369 ! ! tracer content at Before, now and after 370 ztcb = tb(ji,jj,jk) * ze3t_b ; zscb = sb(ji,jj,jk) * ze3t_b 371 ztcn = tn(ji,jj,jk) * ze3t_n ; zscn = sn(ji,jj,jk) * ze3t_n 372 ztca = ta(ji,jj,jk) * ze3t_a ; zsca = sa(ji,jj,jk) * ze3t_a 373 ! 374 ! ! Asselin filter on thickness and tracer content 375 ze3t_f = atfp * ( ze3t_a - 2.* ze3t_n + ze3t_b ) 376 ztc_f = atfp * ( ztca - 2.* ztcn + ztcb ) 377 zsc_f = atfp * ( zsca - 2.* zscn + zscb ) 378 ! 379 ! ! filtered tracer including the correction 380 ze3fr = 1.e0 / ( ze3t_n + ze3t_f ) 381 ztf = ( ztcn + ztc_f ) * ze3fr 382 zsf = ( zscn + zsc_f ) * ze3fr 383 ! ! swap of arrays 384 tb(ji,jj,jk) = ztf ! tb <-- tn filtered 385 sb(ji,jj,jk) = zsf 386 tn(ji,jj,jk) = ta(ji,jj,jk) ! tn <-- ta 387 sn(ji,jj,jk) = sa(ji,jj,jk) 388 END DO 389 END DO 390 END DO 391 ENDIF 392 ENDIF 393 !!gm from Matthieu : unchecked 394 ! 395 IF( neuler /= 0 .OR. kt /= nit000 ) THEN ! remove the forcing from the Asselin filter 396 IF( lk_vvl ) THEN ! Varying levels 397 DO jj = 1, jpj 398 DO ji = 1, jpi 399 ! - ML - modification for varaiance diagnostics 400 zssh1 = tmask(ji,jj,jk) / ( fse3t(ji,jj,jk) + atfp * ( zfse3ta(ji,jj,jk) - 2*fse3t(ji,jj,jk) & 401 & + fse3tb(ji,jj,jk) ) ) 402 zt = zssh1 * ( tn(ji,jj,jk) + atfp * ( ta(ji,jj,jk) - 2 * tn(ji,jj,jk) + tb(ji,jj,jk) ) ) 403 zs = zssh1 * ( sn(ji,jj,jk) + atfp * ( sa(ji,jj,jk) - 2 * sn(ji,jj,jk) + sb(ji,jj,jk) ) ) 404 ! filter correction for global conservation 405 !------------------------------------------ 406 zfact1 = atfp * rdttra(1) * zssh1 407 IF (jk == 1) THEN ! remove sbc trend from time filter 408 zt = zt - zfact1 * ( qn(ji,jj) - qb(ji,jj) ) 409 !!??? zs = zs - zfact1 * ( - emp( ji,jj) + empb(ji,jj) ) * zsrau * 35.5e0 410 ENDIF 411 ! remove qsr trend from time filter 412 IF (jk <= nksr) THEN 413 IF( ln_sco ) THEN 414 z_cor_qsr = etot3(ji,jj,jk) * ( qsr(ji,jj) - qsrb(ji,jj) ) 415 ELSEIF( ln_zps .OR. ln_zco ) THEN 416 z_cor_qsr = ( qsr(ji,jj) - qsrb(ji,jj) ) * & 417 & ( gdsr(jk) * tmask(ji,jj,jk) - gdsr(jk+1) * tmask(ji,jj,jk+1) ) 418 ENDIF 419 zt = zt - zfact1 * z_cor_qsr 420 IF (jk == nksr) qsrb(ji,jj) = qsr(ji,jj) 421 ENDIF 422 ! - ML - end of modification 423 zssh1 = tmask(ji,jj,1) / zfse3ta(ji,jj,jk) 424 tn(ji,jj,jk) = ta(ji,jj,jk) * zssh1 425 sn(ji,jj,jk) = sa(ji,jj,jk) * zssh1 426 END DO 427 END DO 428 ENDIF 429 !!gm end 320 END DO 321 ENDIF 430 322 ! 431 323 END SUBROUTINE tra_nxt_vvl -
branches/DEV_r1837_MLF/NEMO/OPA_SRC/TRA/traqsr.F90
r1870 r1975 47 47 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_chl ! structure of input Chl (file informations, fields read) 48 48 49 INTEGER 49 INTEGER , PUBLIC :: nksr ! levels below which the light cannot penetrate (depth larger than 391 m) 50 50 REAL(wp), DIMENSION(3,61) :: rkrgb ! tabulated attenuation coefficients for RGB absorption 51 51 … … 97 97 REAL(wp) :: zchl, zcoef, zsi0r ! temporary scalars 98 98 REAL(wp) :: zc0, zc1, zc2, zc3 ! - - 99 REAL(wp) :: zqsr ! - -100 99 REAL(wp), DIMENSION(jpi,jpj) :: zekb, zekg, zekr ! 2D workspace 101 100 REAL(wp), DIMENSION(jpi,jpj,jpk) :: ze0, ze1 , ze2, ze3, zea ! 3D workspace … … 123 122 DO ji = fs_2, fs_jpim1 ! vector opt. 124 123 !!gm how to stecify the mean of time step here : TOP versus OPA time stepping strategy not obvious 125 ta(ji,jj,jk) = ta(ji,jj,jk) + ro0cpr * ( etot3(ji,jj,jk) - etot3(ji,jj,jk+1) ) / fse3t(ji,jj,jk)124 qsr_trd_hc_n(ji,jj,jk) = ro0cpr * ( etot3(ji,jj,jk) - etot3(ji,jj,jk+1) ) 126 125 END DO 127 126 END DO … … 155 154 zsi0r = 1.e0 / rn_si0 156 155 zcoef = ( 1. - rn_abs ) / 3.e0 ! equi-partition in R-G-B 157 !!gm bug !!! zqsr only a constant not an array 158 zqsr = 0.5 * ( qsr_b(:,:) + qsr(:,:) ) ! mean over 2 time steps 159 ze0(:,:,1) = rn_abs * zqsr 160 ze1(:,:,1) = zcoef * zqsr 161 ze2(:,:,1) = zcoef * zqsr 162 ze3(:,:,1) = zcoef * zqsr 163 zea(:,:,1) = zqsr 156 157 ze0(:,:,1) = rn_abs * qsr(:,:) 158 ze1(:,:,1) = zcoef * qsr(:,:) 159 ze2(:,:,1) = zcoef * qsr(:,:) 160 ze3(:,:,1) = zcoef * qsr(:,:) 161 zea(:,:,1) = qsr(:,:) 164 162 ! 165 163 DO jk = 2, nksr+1 ! deeper values … … 183 181 ! 184 182 DO jk = 1, nksr ! compute and add qsr trend to ta 185 ta(:,:,jk) = ta(:,:,jk) + ro0cpr * ( zea(:,:,jk) - zea(:,:,jk+1) ) / fse3t(:,:,jk)183 qsr_trd_hc_n(:,:,jk) = ro0cpr * ( zea(:,:,jk) - zea(:,:,jk+1) ) 186 184 END DO 187 185 zea(:,:,nksr+1:jpk) = 0.e0 ! below 400m set to zero … … 190 188 ELSE !* Constant Chlorophyll 191 189 DO jk = 1, nksr 192 ta(:,:,jk) = ta(:,:,jk) + etot3(:,:,jk) * 0.5 * ( qsr_b(:,:) + qsr(:,:))190 qsr_trd_hc_n(:,:,jk) = etot3(:,:,jk) * qsr(:,:) 193 191 END DO 194 192 ENDIF … … 201 199 DO jj = 2, jpjm1 202 200 DO ji = fs_2, fs_jpim1 ! vector opt. 203 ta(ji,jj,jk) = ta(ji,jj,jk) + etot3(ji,jj,jk) * 0.5 * ( qsr_b(:,:) + qsr(:,:))201 qsr_trd_hc_n(ji,jj,jk) = etot3(ji,jj,jk) * qsr(ji,jj) 204 202 END DO 205 203 END DO … … 208 206 ENDIF 209 207 ! 208 ENDIF 209 210 ! Add qsr trend to ta in all cases 211 IF( neuler == 0 .AND. kt == nit000 ) THEN 212 DO jk = 1, nksr 213 DO jj = 2, jpjm1 214 DO ji = fs_2, fs_jpim1 ! vector opt. 215 ta(ji,jj,jk) = ta(ji,jj,jk) + qsr_trd_hc_n(ji,jj,jk) / fse3t(ji,jj,jk) 216 END DO 217 END DO 218 END DO 219 ELSE 220 DO jk = 1, nksr 221 DO jj = 2, jpjm1 222 DO ji = fs_2, fs_jpim1 ! vector opt. 223 ta(ji,jj,jk) = ta(ji,jj,jk) + 0.5 * ( qsr_trd_hc_b(ji,jj,jk) + qsr_trd_hc_n(ji,jj,jk) ) / fse3t(ji,jj,jk) 224 END DO 225 END DO 226 END DO 210 227 ENDIF 211 228 … … 382 399 ! 383 400 DO jk = 1, nksr 384 etot3(:,:,jk) = ro0cpr * ( zea(:,:,jk) - zea(:,:,jk+1) ) / fse3t(:,:,jk)401 etot3(:,:,jk) = ro0cpr * ( zea(:,:,jk) - zea(:,:,jk+1) ) 385 402 END DO 386 403 etot3(:,:,nksr+1:jpk) = 0.e0 ! below 400m set to zero … … 404 421 zc0 = rn_abs * EXP( -fsdepw(ji,jj,jk )*zsi0r ) + (1.-rn_abs) * EXP( -fsdepw(ji,jj,jk )*zsi1r ) 405 422 zc1 = rn_abs * EXP( -fsdepw(ji,jj,jk+1)*zsi0r ) + (1.-rn_abs) * EXP( -fsdepw(ji,jj,jk+1)*zsi1r ) 406 etot3(ji,jj,jk) = ro0cpr * ( zc0 * tmask(ji,jj,jk) - zc1 * tmask(ji,jj,jk+1) ) / fse3t(ji,jj,jk)423 etot3(ji,jj,jk) = ro0cpr * ( zc0 * tmask(ji,jj,jk) - zc1 * tmask(ji,jj,jk+1) ) 407 424 END DO 408 425 END DO -
branches/DEV_r1837_MLF/NEMO/OPA_SRC/TRA/trasbc.F90
r1870 r1975 104 104 INTEGER, INTENT(in) :: kt ! ocean time-step index 105 105 !! 106 INTEGER :: ji, jj 107 REAL(wp) :: z ta, zsa, zsrau, zse3t! temporary scalars106 INTEGER :: ji, jj ! dummy loop indices 107 REAL(wp) :: zsrau, zse3t ! temporary scalars 108 108 !!---------------------------------------------------------------------- 109 109 … … 129 129 qsr(:,:) = 0.e0 ! qsr set to zero 130 130 IF( kt == nit000 ) THEN ! idem on before field at nit000 131 qns_b(:,:) = qns_b(:,:) + qsr_b(:,:)132 131 qsr_b(:,:) = 0.e0 133 132 ENDIF 134 133 ENDIF 135 136 ! Concentration dillution effect on (t,s)137 !! DO jj = 2, jpj138 !! DO ji = fs_2, fs_jpim1 ! vector opt.139 !!#if ! defined key_zco140 !! zse3t = 1. / fse3t(ji,jj,1)141 !!#endif142 !! IF( lk_vvl) THEN143 !! zta = ro0cpr * qns(ji,jj) * zse3t & ! temperature : heat flux144 !! & - emp(ji,jj) * zsrau * tn(ji,jj,1) * zse3t ! & cooling/heating effet of EMP flux145 !! zsa = 0.e0 ! No salinity concent./dilut. effect146 !! ELSE147 !! zta = ro0cpr * qns(ji,jj) * zse3t ! temperature : heat flux148 !! zsa = emps(ji,jj) * zsrau * sn(ji,jj,1) * zse3t ! salinity : concent./dilut. effect149 !! ENDIF150 !! ta(ji,jj,1) = ta(ji,jj,1) + zta ! add the trend to the general tracer trend151 !! sa(ji,jj,1) = sa(ji,jj,1) + zsa152 !! END DO153 !! END DO154 134 155 135 ! ! ---------------------- ! 156 136 IF( lk_vvl ) THEN ! Variable Volume case ! 157 137 ! ! ---------------------- ! 158 DO jj = 2, jpj 159 DO ji = fs_2, fs_jpim1 ! vector opt. 160 zse3t = 0.5 / fse3t(ji,jj,1) 161 ! ! temperature: heat flux + cooling/heating effet of EMP flux 162 ta(ji,jj,1) = ta(ji,jj,1) + ( ro0cpr * ( qns(ji,jj) + qns_b(ji,jj) ) & 163 & - zsrau * ( emp(ji,jj) + emp_b(ji,jj) ) * tn(ji,jj,1) ) * zse3t 164 ! ! salinity: salt flux 165 sa(ji,jj,1) = sa(ji,jj,1) + ( emps(ji,jj) + emps_b(ji,jj) ) * zse3t 166 167 !!gm BUG : in key_vvl emps must be modified to only include the salt flux due to with sea-ice freezing/melting 168 !!gm otherwise this flux will be missing ==> modification required in limsbc, limsbc_2 and CICE interface. 169 170 END DO 171 END DO 138 !!gm BUG : in key_vvl emps must be modified to only include the salt flux due to sea-ice freezing/melting 139 !!gm otherwise this flux will be missing ==> modification required in limsbc, limsbc_2 and CICE interface.s 140 IF ( neuler == 0 .AND. kt == nit000 ) THEN 141 DO jj = 2, jpj 142 DO ji = fs_2, fs_jpim1 ! vector opt. 143 ! temperature : heat flux + cooling/heating effet of EMP flux 144 sbc_trd_hc_n(ji,jj) = ro0cpr * qns(ji,jj) - zsrau * emp(ji,jj) * tn(ji,jj,1) 145 #if ! defined key_zco 146 zse3t = 1. / fse3t(ji,jj,1) 147 #endif 148 ta(ji,jj,1) = ta(ji,jj,1) + zse3t * sbc_trd_hc_n(ji,jj) 149 END DO 150 END DO 151 ELSE 152 DO jj = 2, jpj 153 DO ji = fs_2, fs_jpim1 ! vector opt. 154 ! temperature : heat flux + cooling/heating effet of EMP flux 155 sbc_trd_hc_n(ji,jj) = ro0cpr * qns(ji,jj) - zsrau * emp(ji,jj) * tn(ji,jj,1) 156 #if ! defined key_zco 157 zse3t = 1. / fse3t(ji,jj,1) 158 #endif 159 ta(ji,jj,1) = ta(ji,jj,1) + 0.5 * ( sbc_trd_hc_b(ji,jj) + sbc_trd_hc_n(ji,jj) ) * zse3t 160 END DO 161 END DO 162 ENDIF 172 163 ! ! ---------------------- ! 173 164 ELSE ! Constant Volume case ! 174 165 ! ! ---------------------- ! 175 DO jj = 2, jpj 176 DO ji = fs_2, fs_jpim1 ! vector opt. 177 zse3t = 0.5 / fse3t(ji,jj,1) 178 ! ! temperature: heat flux 179 ta(ji,jj,1) = ta(ji,jj,1) + ro0cpr * ( qns (ji,jj) + qns_b (ji,jj) ) * zse3t 180 ! ! salinity: salt flux + concent./dilut. effect (both in emps) 181 sa(ji,jj,1) = sa(ji,jj,1) + zsrau * ( emps(ji,jj) + emps_b(ji,jj) ) * sn(ji,jj,1) * zse3t 182 END DO 183 END DO 166 IF ( neuler == 0 .AND. kt == nit000 ) THEN 167 DO jj = 2, jpj 168 DO ji = fs_2, fs_jpim1 ! vector opt. 169 ! temperature : heat flux 170 sbc_trd_hc_n(ji,jj) = ro0cpr * qns(ji,jj) 171 ! salinity : salt flux + concent./dilut. effect (both in emps) 172 sbc_trd_sc_n(ji,jj) = zsrau * emps(ji,jj) * sn(ji,jj,1) 173 #if ! defined key_zco 174 zse3t = 1. / fse3t(ji,jj,1) 175 #endif 176 ta(ji,jj,1) = ta(ji,jj,1) + sbc_trd_hc_n(ji,jj) * zse3t 177 sa(ji,jj,1) = sa(ji,jj,1) + sbc_trd_sc_n(ji,jj) * zse3t 178 END DO 179 END DO 180 ELSE 181 DO jj = 2, jpj 182 DO ji = fs_2, fs_jpim1 ! vector opt. 183 ! temperature : heat flux 184 sbc_trd_hc_n(ji,jj) = ro0cpr * qns(ji,jj) 185 ! salinity : salt flux + concent./dilut. effect (both in emps) 186 sbc_trd_sc_n(ji,jj) = zsrau * emps(ji,jj) * sn(ji,jj,1) 187 #if ! defined key_zco 188 zse3t = 1. / fse3t(ji,jj,1) 189 #endif 190 ! temperature : heat flux 191 ta(ji,jj,1) = ta(ji,jj,1) + 0.5 * ( sbc_trd_hc_b(ji,jj) + sbc_trd_hc_n(ji,jj) ) * zse3t 192 sa(ji,jj,1) = sa(ji,jj,1) + 0.5 * ( sbc_trd_sc_b(ji,jj) + sbc_trd_sc_n(ji,jj) ) * zse3t 193 END DO 194 END DO 195 ENDIF 184 196 ! 185 197 ENDIF 186 198 187 IF( l_trdtra ) THEN ! save the sbc trends for diagnostic199 IF( l_trdtra ) THEN ! save the sbc trends for diagnostic 188 200 ztrdt(:,:,:) = ta(:,:,:) - ztrdt(:,:,:) 189 201 ztrds(:,:,:) = sa(:,:,:) - ztrds(:,:,:) 190 202 CALL trd_mod(ztrdt, ztrds, jptra_trd_nsr, 'TRA', kt) 191 203 ENDIF 192 ! ! control print204 ! ! control print 193 205 IF(ln_ctl) CALL prt_ctl( tab3d_1=ta, clinfo1=' sbc - Ta: ', mask1=tmask, & 194 206 & tab3d_2=sa, clinfo2= ' Sa: ', mask2=tmask, clinfo3='tra' )
Note: See TracChangeset
for help on using the changeset viewer.