- Timestamp:
- 2016-11-25T16:41:40+01:00 (7 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2015/dev_r5936_INGV1_WAVE/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfqiao.F90
r7171 r7340 6 6 !! History : 3.6 ! 2014-10 (E. Clementi) Original code 7 7 !!---------------------------------------------------------------------- 8 !!----------------------------------------------------------------------9 !! qiao_init10 8 !! zdf_qiao : compute Qiao parameters 11 9 !!---------------------------------------------------------------------- 12 10 13 USE iom ! I/O manager library14 11 USE in_out_manager ! I/O manager 15 12 USE lib_mpp ! distribued memory computing library … … 18 15 USE sbcwave ! wave module 19 16 USE dom_oce 17 USE lbclnk ! ocean lateral boundary conditions (or mpp link) 20 18 21 !!----------------------------------------------------------------------22 !! qiao_init : compute QBv: Qiao terms to be added to vertical eddy23 !! diffusivity and viscosity coefficients24 !!----------------------------------------------------------------------25 26 19 IMPLICIT NONE 27 20 PRIVATE 28 21 29 PUBLIC zdf_qiao ! routine called in zdf_ric22 PUBLIC zdf_qiao ! routine called in step 30 23 31 REAL(wp), PUBLIC,ALLOCATABLE,DIMENSION (:,:,:) :: QBv, QBvu, QBvv24 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:,:) :: qbv, qbvu, qbvv 32 25 33 26 !! * Substitutions 34 27 # include "domzgr_substitute.h90" 28 # include "vectopt_loop_substitute.h90" 35 29 !!---------------------------------------------------------------------- 36 30 !! NEMO/OPA 4.0 , NEMO Consortium (2011) … … 45 39 !! *** ROUTINE zdf_qiao *** 46 40 !! 47 !! ** Purpose :Compute the Qiao term ( QBv) to be added to41 !! ** Purpose :Compute the Qiao term (qbv) to be added to 48 42 !! vertical viscosity and diffusivity coeffs. 49 43 !! 50 !! ** Method : QBv = alpha * A * Us(0) * exp (3 * k * z)44 !! ** Method :qbv = alpha * A * Us(0) * exp (3 * k * z) 51 45 !! 52 46 !! ** action :Compute the Qiao wave dependent term … … 56 50 INTEGER, INTENT( in ) :: kt ! ocean time step 57 51 ! 58 INTEGER :: jj, ji, jk52 INTEGER :: jj, ji, jk ! dummy loop indices 59 53 !!--------------------------------------------------------------------- 60 54 ! 61 !62 ! ! -------------------- !63 55 IF( kt == nit000 ) THEN ! First call kt=nit000 ! 64 ALLOCATE(QBv(jpi,jpj,jpk)) ! -------------------- ! 65 ALLOCATE(QBvu(jpi,jpj,jpk)) 66 ALLOCATE(QBvv(jpi,jpj,jpk)) 56 IF( .NOT. ( ln_wave .AND. ln_sdw ) ) & 57 & CALL ctl_stop ( 'Ask for wave Qiao enhanced turbulence but ln_wave & 58 & and ln_sdw have to be activated') 59 IF( zdf_qiao_alloc() /= 0 ) & 60 & CALL ctl_stop( 'STOP', 'zdf_qiao : unable to allocate arrays' ) 67 61 ENDIF 68 62 69 QBv (:,:,:) = 0.070 QBvu(:,:,:) = 0.071 QBvv(:,:,:) = 0.072 73 63 ! 74 ! Compute the Qiao term Bv ( QBv) to be added to64 ! Compute the Qiao term Bv (qbv) to be added to 75 65 ! vertical viscosity and diffusivity 76 ! QBv = alpha * A * Us(0) * exp (3 * k * z)66 ! qbv = alpha * A * Us(0) * exp (3 * k * z) 77 67 ! alpha here is set to 1 78 68 !--------------------------------------------------------------------------------- 79 69 ! 80 IF ( ln_wave ) THEN 81 DO jk = 1, jpk 82 DO jj = 1, jpjm1 83 DO ji = 1, jpim1 84 QBv(ji,jj,jk) = 1.0 * 0.353553 * swh(ji,jj) * tsd2d(ji,jj) * & 85 & exp(3.0 * wnum(ji,jj) * & 86 & (-MIN( fsdept(ji ,jj ,jk) , fsdept(ji+1,jj ,jk), & 87 & fsdept(ji ,jj+1,jk) , fsdept(ji+1,jj+1,jk)))) 88 END DO 70 DO jk = 1, jpk 71 DO jj = 1, jpjm1 72 DO ji = 1, fs_jpim1 73 qbv(ji,jj,jk) = 1.0 * 0.353553 * swh(ji,jj) * tsd2d(ji,jj) * & 74 & EXP(3.0 * wnum(ji,jj) * & 75 & (-MIN( fsdepw(ji ,jj ,jk), fsdepw(ji+1,jj ,jk), & 76 & fsdepw(ji ,jj+1,jk), fsdepw(ji+1,jj+1,jk)))) & 77 & * wmask(ji,jj,jk) 89 78 END DO 90 79 END DO 91 92 QBv(jpi,:,:)=QBv(jpim1,:,:) 93 QBv(:,jpj,:)=QBv(:,jpjm1,:) 94 95 ! 96 ! Interpolate Qiao parameter QBv into the grid_U and grid_V 97 !------------------------------------------------- 98 ! 99 DO jk = 1, jpk 100 DO jj = 1, jpjm1 101 DO ji = 1, jpim1 102 QBvu(ji,jj,jk) = 0.5 * umask(ji,jj,jk) * & 103 & ( QBv(ji ,jj,jk) * tmask(ji ,jj,jk) & 104 & + QBv(ji+1,jj,jk) * tmask(ji+1,jj,jk) ) 105 QBvv(ji,jj,jk) = 0.5 * vmask(ji,jj,jk) * & 106 & ( QBv(ji,jj ,jk) * tmask(ji,jj ,jk) & 107 & + QBv(ji,jj+1,jk) * tmask(ji,jj+1,jk) ) 108 END DO 80 END DO 81 ! 82 CALL lbc_lnk( qbv, 'W', 1. ) ! Lateral boundary conditions 83 84 ! 85 ! Interpolate Qiao parameter qbv into the grid_U and grid_V 86 !---------------------------------------------------------- 87 ! 88 DO jk = 1, jpk 89 DO jj = 1, jpjm1 90 DO ji = 1, fs_jpim1 91 qbvu(ji,jj,jk) = 0.5 * wumask(ji,jj,jk) * & 92 & ( qbv(ji,jj,jk) + qbv(ji+1,jj ,jk) ) 93 qbvv(ji,jj,jk) = 0.5 * wvmask(ji,jj,jk) * & 94 & ( qbv(ji,jj,jk) + qbv(ji ,jj+1,jk) ) 109 95 END DO 110 96 END DO 111 ! 112 QBvu(jpi,:,:)=QBvu(jpim1,:,:) 113 QBvu(:,jpj,:)=QBvu(:,jpjm1,:) 114 QBvv(jpi,:,:)=QBvv(jpim1,:,:) 115 QBvv(:,jpj,:)=QBvv(:,jpjm1,:) 97 END DO 98 ! 99 CALL lbc_lnk( qbvu, 'U', 1. ) ; CALL lbc_lnk( qbvv, 'V', 1. ) ! Lateral boundary conditions 116 100 117 ELSE 118 CALL ctl_stop( 'STOP', 'To use Qiao formulation you have to set: ln_wave=.true.') 119 ENDIF 120 ! 101 ! Enhance vertical mixing coeff. 102 !------------------------------- 103 ! 104 DO jk = 1, jpkm1 105 DO jj = 1, jpj 106 DO ji = 1, jpi 107 avmu(ji,jj,jk) = ( avmu(ji,jj,jk) + qbvu(ji,jj,jk) ) * umask(ji,jj,jk) 108 avmv(ji,jj,jk) = ( avmv(ji,jj,jk) + qbvv(ji,jj,jk) ) * vmask(ji,jj,jk) 109 avt (ji,jj,jk) = ( avt (ji,jj,jk) + qbv (ji,jj,jk) ) * tmask(ji,jj,jk) 110 END DO 111 END DO 112 END DO 113 ! 121 114 END SUBROUTINE zdf_qiao 115 116 INTEGER FUNCTION zdf_qiao_alloc() 117 !!---------------------------------------------------------------------- 118 !! *** FUNCTION zdf_qiao_alloc *** 119 !!---------------------------------------------------------------------- 120 ALLOCATE( qbv(jpi,jpj,jpk), qbvu(jpi,jpj,jpk), qbvv(jpi,jpj,jpk), & 121 & STAT = zdf_qiao_alloc ) 122 ! 123 IF( lk_mpp ) CALL mpp_sum ( zdf_qiao_alloc ) 124 IF( zdf_qiao_alloc > 0 ) CALL ctl_warn('zdf_qiao_alloc: allocation of arrays failed.') 125 ! 126 END FUNCTION zdf_qiao_alloc 122 127 123 128 !!======================================================================
Note: See TracChangeset
for help on using the changeset viewer.
