Calculates the vertical background diffusivities/viscosities.
310 type(ocean_grid_type),
intent(in) :: G
311 type(verticalGrid_type),
intent(in) :: GV
312 real,
dimension(SZI_(G),SZJ_(G),SZK_(GV)),
intent(in) :: h
313 type(thermo_var_ptrs),
intent(in) :: tv
314 real,
dimension(SZI_(G),SZK_(GV)),
intent(in) :: N2_lay
316 real,
dimension(SZI_(G),SZK_(GV)),
intent(out) :: Kd_lay
318 real,
dimension(SZI_(G),SZK_(GV)+1),
intent(out) :: Kd_int
320 real,
dimension(SZI_(G),SZK_(GV)+1),
intent(out) :: Kv_bkgnd
322 integer,
intent(in) :: j
323 type(unit_scale_type),
intent(in) :: US
324 type(bkgnd_mixing_cs),
pointer :: CS
328 real,
dimension(SZK_(GV)+1) :: depth_int
329 real,
dimension(SZK_(GV)+1) :: Kd_col
330 real,
dimension(SZK_(GV)+1) :: Kv_col
331 real,
dimension(SZI_(G)) :: Kd_sfc
332 real,
dimension(SZI_(G)) :: depth
342 real :: bckgrnd_vdc_psin
343 real :: bckgrnd_vdc_psis
344 integer :: i, k, is, ie, js, je, nz
346 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
349 deg_to_rad = atan(1.0)/45.0
357 if (cs%Bryan_Lewis_diffusivity)
then
362 depth_int(k) = depth_int(k-1) + gv%H_to_m*h(i,j,k-1)
365 call cvmix_init_bkgnd(max_nlev=nz, &
367 bl1 = cs%Bryan_Lewis_c1, &
368 bl2 = cs%Bryan_Lewis_c2, &
369 bl3 = cs%Bryan_Lewis_c3, &
370 bl4 = cs%Bryan_Lewis_c4, &
371 prandtl = cs%prandtl_bkgnd)
373 kd_col(:) = 0.0 ; kv_col(:) = 0.0
374 call cvmix_coeffs_bkgnd(mdiff_out=kv_col, tdiff_out=kd_col, nlev=nz, max_nlev=nz)
378 kv_bkgnd(i,k) = us%m2_s_to_Z2_T*kv_col(k)
379 kd_int(i,k) = us%m2_s_to_Z2_T*kd_col(k)
382 kd_lay(i,k) = kd_lay(i,k) + 0.5 * us%m2_s_to_Z2_T * (kd_col(k) + kd_col(k+1))
386 elseif (cs%horiz_varying_background)
then
389 bckgrnd_vdc_psis = cs%bckgrnd_vdc_psim * exp(-(0.4*(g%geoLatT(i,j)+28.9))**2)
390 bckgrnd_vdc_psin = cs%bckgrnd_vdc_psim * exp(-(0.4*(g%geoLatT(i,j)-28.9))**2)
391 kd_int(i,1) = (cs%bckgrnd_vdc_eq + bckgrnd_vdc_psin) + bckgrnd_vdc_psis
393 if (g%geoLatT(i,j) < -10.0)
then
394 kd_int(i,1) = kd_int(i,1) + cs%bckgrnd_vdc1
395 elseif (g%geoLatT(i,j) <= 10.0)
then
396 kd_int(i,1) = kd_int(i,1) + cs%bckgrnd_vdc1 * (g%geoLatT(i,j)/10.0)**2
398 kd_int(i,1) = kd_int(i,1) + cs%bckgrnd_vdc1
402 if ( (g%geoLatT(i,j) < -1.0) .and. (g%geoLatT(i,j) > -4.0) .and. &
403 ( mod(g%geoLonT(i,j)+360.0,360.0) > 103.0) .and. &
404 ( mod(g%geoLonT(i,j)+360.0,360.0) < 134.0) )
then
405 kd_int(i,1) = cs%bckgrnd_vdc_Banda
409 if ( (g%geoLatT(i,j) <= -4.0) .and. (g%geoLatT(i,j) > -7.0) .and. &
410 ( mod(g%geoLonT(i,j)+360.0,360.0) > 106.0) .and. &
411 ( mod(g%geoLonT(i,j)+360.0,360.0) < 140.0) )
then
412 kd_int(i,1) = cs%bckgrnd_vdc_Banda
416 if ( (g%geoLatT(i,j) <= -7.0) .and. (g%geoLatT(i,j) > -8.3) .and. &
417 ( mod(g%geoLonT(i,j)+360.0,360.0) > 111.0) .and. &
418 ( mod(g%geoLonT(i,j)+360.0,360.0) < 142.0) )
then
419 kd_int(i,1) = cs%bckgrnd_vdc_Banda
424 do k=1,nz+1 ;
do i=is,ie
425 kd_int(i,k) = kd_int(i,1)
426 kv_bkgnd(i,k) = kd_int(i,1) * cs%prandtl_bkgnd
428 do k=1,nz ;
do i=is,ie
429 kd_lay(i,k) = kd_int(i,1)
432 elseif (cs%Henyey_IGW_background_new)
then
433 i_x30 = 2.0 / invcosh(cs%N0_2Omega*2.0)
434 i_2omega = 0.5 / cs%omega
435 do k=1,nz ;
do i=is,ie
436 abs_sinlat = max(min_sinlat, abs(sin(g%geoLatT(i,j)*deg_to_rad)))
437 n_2omega = max(abs_sinlat, sqrt(n2_lay(i,k))*i_2omega)
438 n02_n2 = (cs%N0_2Omega/n_2omega)**2
439 kd_lay(i,k) = max(cs%Kd_min, cs%Kd * &
440 ((abs_sinlat * invcosh(n_2omega/abs_sinlat)) * i_x30)*n02_n2)
444 kd_int(i,1) = 0.0; kv_bkgnd(i,1) = 0.0
445 kd_int(i,nz+1) = 0.0; kv_bkgnd(i,nz+1) = 0.0
447 do k=2,nz ;
do i=is,ie
448 kd_int(i,k) = 0.5*(kd_lay(i,k-1) + kd_lay(i,k))
449 kv_bkgnd(i,k) = kd_int(i,k) * cs%prandtl_bkgnd
453 if (cs%Henyey_IGW_background)
then
454 i_x30 = 2.0 / invcosh(cs%N0_2Omega*2.0)
456 abs_sinlat = abs(sin(g%geoLatT(i,j)*deg_to_rad))
457 kd_sfc(i) = max(cs%Kd_min, cs%Kd * &
458 ((abs_sinlat * invcosh(cs%N0_2Omega / max(min_sinlat, abs_sinlat))) * i_x30) )
460 elseif (cs%Kd_tanh_lat_fn)
then
465 kd_sfc(i) = max(cs%Kd_min, cs%Kd * (1.0 + &
466 cs%Kd_tanh_lat_scale * 0.5*tanh((abs(g%geoLatT(i,j)) - 35.0)/5.0) ))
475 if ((.not.cs%bulkmixedlayer) .and. (cs%Kd /= cs%Kdml))
then
478 i_hmix = 1.0 / (cs%Hmix + gv%H_subroundoff*gv%H_to_Z)
479 do i=is,ie ; depth(i) = 0.0 ;
enddo
480 do k=1,nz ;
do i=is,ie
481 depth_c = depth(i) + 0.5*gv%H_to_Z*h(i,j,k)
482 if (cs%Kd_via_Kdml_bug)
then
485 if (depth_c <= cs%Hmix)
then ; kd_int(i,k) = cs%Kdml
486 elseif (depth_c >= 2.0*cs%Hmix)
then ; kd_int(i,k) = kd_sfc(i)
488 kd_lay(i,k) = ((kd_sfc(i) - cs%Kdml) * i_hmix) * depth_c + (2.0*cs%Kdml - kd_sfc(i))
491 if (depth_c <= cs%Hmix)
then ; kd_lay(i,k) = cs%Kdml
492 elseif (depth_c >= 2.0*cs%Hmix)
then ; kd_lay(i,k) = kd_sfc(i)
494 kd_lay(i,k) = ((kd_sfc(i) - cs%Kdml) * i_hmix) * depth_c + (2.0*cs%Kdml - kd_sfc(i))
498 depth(i) = depth(i) + gv%H_to_Z*h(i,j,k)
502 do k=1,nz ;
do i=is,ie
503 kd_lay(i,k) = kd_sfc(i)
509 kd_int(i,1) = 0.0; kv_bkgnd(i,1) = 0.0
510 kd_int(i,nz+1) = 0.0; kv_bkgnd(i,nz+1) = 0.0
512 do k=2,nz ;
do i=is,ie
513 kd_int(i,k) = 0.5*(kd_lay(i,k-1) + kd_lay(i,k))
514 kv_bkgnd(i,k) = kd_int(i,k) * cs%prandtl_bkgnd