Changeset 14062 for NEMO/branches/2020/dev_r13333_KERNEL-08_techene_gm_HPG_SPG/src/OCE/DYN/dynldf_lap_blp.F90
- Timestamp:
- 2020-12-03T17:39:30+01:00 (3 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2020/dev_r13333_KERNEL-08_techene_gm_HPG_SPG/src/OCE/DYN/dynldf_lap_blp.F90
r14037 r14062 5 5 !!====================================================================== 6 6 !! History : 3.7 ! 2014-01 (G. Madec, S. Masson) Original code, re-entrant laplacian 7 !! 4.0 ! 2020-04 (A. Nasser, G. Madec) Add symmetric mixing tensor 7 8 !!---------------------------------------------------------------------- 8 9 … … 19 20 USE in_out_manager ! I/O manager 20 21 USE lbclnk ! ocean lateral boundary conditions (or mpp link) 21 22 USE lib_mpp 23 22 24 IMPLICIT NONE 23 25 PRIVATE … … 47 49 !! 48 50 !! ** Action : - pu_rhs, pv_rhs increased by the harmonic operator applied on pu, pv. 51 !! 52 !! Reference : S.Griffies, R.Hallberg 2000 Mon.Wea.Rev., DOI:/ 49 53 !!---------------------------------------------------------------------- 50 54 INTEGER , INTENT(in ) :: kt ! ocean time-step index … … 57 61 REAL(wp) :: zsign ! local scalars 58 62 REAL(wp) :: zua, zva ! local scalars 59 REAL(wp), DIMENSION(jpi,jpj) :: zcur, zdiv 63 REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: zcur, zdiv 64 REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: zten, zshe ! tension (diagonal) and shearing (anti-diagonal) terms 60 65 !!---------------------------------------------------------------------- 61 66 ! … … 70 75 ENDIF 71 76 ! 72 ! ! =============== 73 DO jk = 1, jpkm1 ! Horizontal slab 74 ! ! =============== 75 DO_2D( 0, 1, 0, 1 ) 76 ! ! ahm * e3 * curl (computed from 1 to jpim1/jpjm1) 77 zcur(ji-1,jj-1) = ahmf(ji-1,jj-1,jk) * e3f(ji-1,jj-1,jk) * r1_e1e2f(ji-1,jj-1) & ! ahmf already * by fmask 78 & * ( e2v(ji ,jj-1) * pv(ji ,jj-1,jk) - e2v(ji-1,jj-1) * pv(ji-1,jj-1,jk) & 79 & - e1u(ji-1,jj ) * pu(ji-1,jj ,jk) + e1u(ji-1,jj-1) * pu(ji-1,jj-1,jk) ) 80 ! ! ahm * div (computed from 2 to jpi/jpj) 81 zdiv(ji,jj) = ahmt(ji,jj,jk) * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kbb) & ! ahmt already * by tmask 82 & * ( e2u(ji,jj)*e3u(ji,jj,jk,Kbb) * pu(ji,jj,jk) - e2u(ji-1,jj)*e3u(ji-1,jj,jk,Kbb) * pu(ji-1,jj,jk) & 83 & + e1v(ji,jj)*e3v(ji,jj,jk,Kbb) * pv(ji,jj,jk) - e1v(ji,jj-1)*e3v(ji,jj-1,jk,Kbb) * pv(ji,jj-1,jk) ) 84 END_2D 77 SELECT CASE( nn_dynldf_typ ) 78 ! 79 CASE ( np_typ_rot ) !== Vorticity-Divergence operator ==! 85 80 ! 86 DO_2D( 0, 0, 0, 0 ) ! - curl( curl) + grad( div ) 87 pu_rhs(ji,jj,jk) = pu_rhs(ji,jj,jk) + zsign * umask(ji,jj,jk) * ( & ! * by umask is mandatory for dyn_ldf_blp use 88 & - ( zcur(ji ,jj) - zcur(ji,jj-1) ) * r1_e2u(ji,jj) / e3u(ji,jj,jk,Kmm) & 89 & + ( zdiv(ji+1,jj) - zdiv(ji,jj ) ) * r1_e1u(ji,jj) ) 81 ALLOCATE( zcur(jpi,jpj) , zdiv(jpi,jpj) ) 82 ! 83 DO jk = 1, jpkm1 ! Horizontal slab 84 ! 85 DO_2D( 0, 1, 0, 1 ) 86 ! ! ahm * e3 * curl (computed from 1 to jpim1/jpjm1) 87 zcur(ji-1,jj-1) = ahmf(ji-1,jj-1,jk) * e3f(ji-1,jj-1,jk) * r1_e1e2f(ji-1,jj-1) & ! ahmf already * by fmask 88 & * ( e2v(ji ,jj-1) * pv(ji ,jj-1,jk) - e2v(ji-1,jj-1) * pv(ji-1,jj-1,jk) & 89 & - e1u(ji-1,jj ) * pu(ji-1,jj ,jk) + e1u(ji-1,jj-1) * pu(ji-1,jj-1,jk) ) 90 ! ! ahm * div (computed from 2 to jpi/jpj) 91 zdiv(ji,jj) = ahmt(ji,jj,jk) * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kbb) & ! ahmt already * by tmask 92 & * ( e2u(ji,jj)*e3u(ji,jj,jk,Kbb) * pu(ji,jj,jk) - e2u(ji-1,jj)*e3u(ji-1,jj,jk,Kbb) * pu(ji-1,jj,jk) & 93 & + e1v(ji,jj)*e3v(ji,jj,jk,Kbb) * pv(ji,jj,jk) - e1v(ji,jj-1)*e3v(ji,jj-1,jk,Kbb) * pv(ji,jj-1,jk) ) 94 END_2D 95 ! 96 DO_2D( 0, 0, 0, 0 ) ! - curl( curl) + grad( div ) 97 pu_rhs(ji,jj,jk) = pu_rhs(ji,jj,jk) + zsign * umask(ji,jj,jk) * ( & ! * by umask is mandatory for dyn_ldf_blp use 98 & - ( zcur(ji ,jj) - zcur(ji,jj-1) ) * r1_e2u(ji,jj) / e3u(ji,jj,jk,Kmm) & 99 & + ( zdiv(ji+1,jj) - zdiv(ji,jj ) ) * r1_e1u(ji,jj) ) 90 100 ! 91 pv_rhs(ji,jj,jk) = pv_rhs(ji,jj,jk) + zsign * vmask(ji,jj,jk) * ( & ! * by vmask is mandatory for dyn_ldf_blp use 92 & ( zcur(ji,jj ) - zcur(ji-1,jj) ) * r1_e1v(ji,jj) / e3v(ji,jj,jk,Kmm) & 93 & + ( zdiv(ji,jj+1) - zdiv(ji ,jj) ) * r1_e2v(ji,jj) ) 94 END_2D 95 ! ! =============== 96 END DO ! End of slab 97 ! ! =============== 101 pv_rhs(ji,jj,jk) = pv_rhs(ji,jj,jk) + zsign * vmask(ji,jj,jk) * ( & ! * by vmask is mandatory for dyn_ldf_blp use 102 & ( zcur(ji,jj ) - zcur(ji-1,jj) ) * r1_e1v(ji,jj) / e3v(ji,jj,jk,Kmm) & 103 & + ( zdiv(ji,jj+1) - zdiv(ji ,jj) ) * r1_e2v(ji,jj) ) 104 END_2D 105 ! 106 END DO ! End of slab 107 ! 108 DEALLOCATE( zcur , zdiv ) 109 ! 110 CASE ( np_typ_sym ) !== Symmetric operator ==! 111 ! 112 ALLOCATE( zten(jpi,jpj) , zshe(jpi,jpj) ) 113 ! 114 DO jk = 1, jpkm1 ! Horizontal slab 115 ! 116 DO_2D( 0, 1, 0, 1 ) 117 ! ! shearing stress component (F-point) NB : ahmf has already been multiplied by fmask 118 zshe(ji-1,jj-1) = ahmf(ji-1,jj-1,jk) & 119 & * ( e1f(ji-1,jj-1) * r1_e2f(ji-1,jj-1) & 120 & * ( pu(ji-1,jj ,jk) * r1_e1u(ji-1,jj ) - pu(ji-1,jj-1,jk) * r1_e1u(ji-1,jj-1) ) & 121 & + e2f(ji-1,jj-1) * r1_e1f(ji-1,jj-1) & 122 & * ( pv(ji ,jj-1,jk) * r1_e2v(ji ,jj-1) - pv(ji-1,jj-1,jk) * r1_e2v(ji-1,jj-1) ) ) 123 ! ! tension stress component (T-point) NB : ahmt has already been multiplied by tmask 124 zten(ji,jj) = ahmt(ji,jj,jk) & 125 & * ( e2t(ji,jj) * r1_e1t(ji,jj) & 126 & * ( pu(ji,jj,jk) * r1_e2u(ji,jj) - pu(ji-1,jj,jk) * r1_e2u(ji-1,jj) ) & 127 & - e1t(ji,jj) * r1_e2t(ji,jj) & 128 & * ( pv(ji,jj,jk) * r1_e1v(ji,jj) - pv(ji,jj-1,jk) * r1_e1v(ji,jj-1) ) ) 129 END_2D 130 ! 131 DO_2D( 0, 0, 0, 0 ) 132 pu_rhs(ji,jj,jk) = pu_rhs(ji,jj,jk) + zsign * r1_e1e2u(ji,jj) / e3u(ji,jj,jk,Kmm) & 133 & * ( ( zten(ji+1,jj ) * e2t(ji+1,jj )*e2t(ji+1,jj ) * e3t(ji+1,jj ,jk,Kmm) & 134 & - zten(ji ,jj ) * e2t(ji ,jj )*e2t(ji ,jj ) * e3t(ji ,jj ,jk,Kmm) ) * r1_e2u(ji,jj) & 135 & + ( zshe(ji ,jj ) * e1f(ji ,jj )*e1f(ji ,jj ) * e3f(ji ,jj ,jk) & 136 & - zshe(ji ,jj-1) * e1f(ji ,jj-1)*e1f(ji ,jj-1) * e3f(ji ,jj-1,jk) ) * r1_e1u(ji,jj) ) 137 ! 138 pv_rhs(ji,jj,jk) = pv_rhs(ji,jj,jk) + zsign * r1_e1e2v(ji,jj) / e3v(ji,jj,jk,Kmm) & 139 & * ( ( zshe(ji ,jj ) * e2f(ji ,jj )*e2f(ji ,jj ) * e3f(ji ,jj ,jk) & 140 & - zshe(ji-1,jj ) * e2f(ji-1,jj )*e2f(ji-1,jj ) * e3f(ji-1,jj ,jk) ) * r1_e2v(ji,jj) & 141 & - ( zten(ji ,jj+1) * e1t(ji ,jj+1)*e1t(ji ,jj+1) * e3t(ji ,jj+1,jk,Kmm) & 142 & - zten(ji ,jj ) * e1t(ji ,jj )*e1t(ji ,jj ) * e3t(ji ,jj ,jk,Kmm) ) * r1_e1v(ji,jj) ) 143 ! 144 END_2D 145 ! 146 END DO 147 ! 148 DEALLOCATE( zten , zshe ) 149 ! 150 END SELECT 98 151 ! 99 152 END SUBROUTINE dyn_ldf_lap
Note: See TracChangeset
for help on using the changeset viewer.
