Changeset 474 for trunk/NEMO/OPA_SRC/DYN/dynldf_bilap.F90
- Timestamp:
- 2006-05-11T17:24:19+02:00 (18 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMO/OPA_SRC/DYN/dynldf_bilap.F90
r455 r474 86 86 REAL(wp) :: zua, zva, zbt, ze2u, ze2v ! temporary scalar 87 87 REAL(wp), DIMENSION(jpi,jpj) :: & 88 zuf, zut, zlu, zlv, zcu, zcv ! temporary workspace 88 zcu, zcv ! temporary workspace 89 REAL(wp), DIMENSION(jpi,jpj,jpk) :: & 90 zuf, zut, zlu, zlv ! temporary workspace 89 91 !!---------------------------------------------------------------------- 90 92 !! OPA 8.5, LODYC-IPSL (2002) … … 96 98 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~' 97 99 ENDIF 98 zuf(:,:) = 0.e0 99 zut(:,:) = 0.e0 100 zlu(:,:) = 0.e0 101 zlv(:,:) = 0.e0 100 101 !!bug gm this should be enough 102 !!$ zuf(:,:,jpk) = 0.e0 103 !!$ zut(:,:,jpk) = 0.e0 104 !!$ zlu(:,:,jpk) = 0.e0 105 !!$ zlv(:,:,jpk) = 0.e0 106 zuf(:,:,:) = 0.e0 107 zut(:,:,:) = 0.e0 108 zlu(:,:,:) = 0.e0 109 zlv(:,:,:) = 0.e0 102 110 103 111 ! ! =============== … … 108 116 109 117 IF( ln_sco .OR. ln_zps ) THEN ! s-coordinate or z-coordinate with partial steps 110 zuf(:,: ) = rotb(:,:,jk) * fse3f(:,:,jk)118 zuf(:,:,jk) = rotb(:,:,jk) * fse3f(:,:,jk) 111 119 DO jj = 2, jpjm1 112 120 DO ji = fs_2, fs_jpim1 ! vector opt. 113 zlu(ji,jj ) = - ( zuf(ji,jj) - zuf(ji,jj-1) ) / ( e2u(ji,jj) * fse3u(ji,jj,jk) ) &121 zlu(ji,jj,jk) = - ( zuf(ji,jj,jk) - zuf(ji,jj-1,jk) ) / ( e2u(ji,jj) * fse3u(ji,jj,jk) ) & 114 122 & + ( hdivb(ji+1,jj,jk) - hdivb(ji,jj,jk) ) / e1u(ji,jj) 115 123 116 zlv(ji,jj ) = + ( zuf(ji,jj) - zuf(ji-1,jj) ) / ( e1v(ji,jj) * fse3v(ji,jj,jk) ) &124 zlv(ji,jj,jk) = + ( zuf(ji,jj,jk) - zuf(ji-1,jj,jk) ) / ( e1v(ji,jj) * fse3v(ji,jj,jk) ) & 117 125 & + ( hdivb(ji,jj+1,jk) - hdivb(ji,jj,jk) ) / e2v(ji,jj) 118 126 END DO … … 121 129 DO jj = 2, jpjm1 122 130 DO ji = fs_2, fs_jpim1 ! vector opt. 123 zlu(ji,jj ) = - ( rotb (ji ,jj,jk) - rotb (ji,jj-1,jk) ) / e2u(ji,jj) &131 zlu(ji,jj,jk) = - ( rotb (ji ,jj,jk) - rotb (ji,jj-1,jk) ) / e2u(ji,jj) & 124 132 & + ( hdivb(ji+1,jj,jk) - hdivb(ji,jj ,jk) ) / e1u(ji,jj) 125 133 126 zlv(ji,jj ) = + ( rotb (ji,jj ,jk) - rotb (ji-1,jj,jk) ) / e1v(ji,jj) &134 zlv(ji,jj,jk) = + ( rotb (ji,jj ,jk) - rotb (ji-1,jj,jk) ) / e1v(ji,jj) & 127 135 & + ( hdivb(ji,jj+1,jk) - hdivb(ji ,jj,jk) ) / e2v(ji,jj) 128 136 END DO 129 137 END DO 130 138 ENDIF 131 132 ! Boundary conditions on the laplacian (zlu,zlv) 133 CALL lbc_lnk( zlu, 'U', -1. ) 134 CALL lbc_lnk( zlv, 'V', -1. ) 135 136 139 ENDDO 140 141 ! Boundary conditions on the laplacian (zlu,zlv) 142 CALL lbc_lnk( zlu, 'U', -1. ) 143 CALL lbc_lnk( zlv, 'V', -1. ) 144 145 146 DO jk = 1, jpkm1 147 137 148 ! Third derivative 138 149 ! ---------------- 139 150 140 151 ! Multiply by the eddy viscosity coef. (at u- and v-points) 141 zlu(:,: ) = zlu(:,:) * fsahmu(:,:,jk)142 zlv(:,: ) = zlv(:,:) * fsahmv(:,:,jk)152 zlu(:,:,jk) = zlu(:,:,jk) * fsahmu(:,:,jk) 153 zlv(:,:,jk) = zlv(:,:,jk) * fsahmv(:,:,jk) 143 154 144 155 ! Contravariant "laplacian" 145 zcu(:,:) = e1u(:,:) * zlu(:,: )146 zcv(:,:) = e2v(:,:) * zlv(:,: )156 zcu(:,:) = e1u(:,:) * zlu(:,:,jk) 157 zcv(:,:) = e2v(:,:) * zlv(:,:,jk) 147 158 148 159 ! Laplacian curl ( * e3f if s-coordinates or z-coordinate with partial steps) 149 160 DO jj = 1, jpjm1 150 161 DO ji = 1, fs_jpim1 ! vector opt. 151 zuf(ji,jj ) = fmask(ji,jj,jk) * ( zcv(ji+1,jj ) - zcv(ji,jj) &162 zuf(ji,jj,jk) = fmask(ji,jj,jk) * ( zcv(ji+1,jj ) - zcv(ji,jj) & 152 163 & - zcu(ji ,jj+1) + zcu(ji,jj) ) & 153 164 #if defined key_zco … … 163 174 DO ji = 1, fs_jpim1 ! vector opt. 164 175 #if defined key_zco 165 zlu(ji,jj ) = e2u(ji,jj) * zlu(ji,jj)166 zlv(ji,jj ) = e1v(ji,jj) * zlv(ji,jj)167 #else 168 zlu(ji,jj ) = e2u(ji,jj) * fse3u(ji,jj,jk) * zlu(ji,jj)169 zlv(ji,jj ) = e1v(ji,jj) * fse3v(ji,jj,jk) * zlv(ji,jj)176 zlu(ji,jj,jk) = e2u(ji,jj) * zlu(ji,jj,jk) 177 zlv(ji,jj,jk) = e1v(ji,jj) * zlv(ji,jj,jk) 178 #else 179 zlu(ji,jj,jk) = e2u(ji,jj) * fse3u(ji,jj,jk) * zlu(ji,jj,jk) 180 zlv(ji,jj,jk) = e1v(ji,jj) * fse3v(ji,jj,jk) * zlv(ji,jj,jk) 170 181 #endif 171 182 END DO … … 180 191 zbt = e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) 181 192 #endif 182 zut(ji,jj ) = ( zlu(ji,jj) - zlu(ji-1,jj) &183 & + zlv(ji,jj) - zlv(ji ,jj-1)) / zbt193 zut(ji,jj,jk) = ( zlu(ji,jj,jk) - zlu(ji-1,jj ,jk) & 194 & + zlv(ji,jj,jk) - zlv(ji ,jj-1,jk) ) / zbt 184 195 END DO 185 196 END DO 197 END DO 186 198 187 199 188 200 ! boundary conditions on the laplacian curl and div (zuf,zut) 201 !!bug gm no need to do this 2 following lbc... 189 202 CALL lbc_lnk( zuf, 'F', 1. ) 190 203 CALL lbc_lnk( zut, 'T', 1. ) 191 204 192 205 DO jk = 1, jpkm1 206 193 207 ! Bilaplacian 194 208 ! ----------- … … 204 218 #endif 205 219 ! horizontal biharmonic diffusive trends 206 zua = - ( zuf(ji ,jj ) - zuf(ji,jj-1) ) / ze2u &207 & + ( zut(ji+1,jj ) - zut(ji,jj) ) / e1u(ji,jj)208 209 zva = + ( zuf(ji,jj ) - zuf(ji-1,jj) ) / ze2v &210 & + ( zut(ji,jj+1 ) - zut(ji ,jj) ) / e2v(ji,jj)220 zua = - ( zuf(ji ,jj,jk) - zuf(ji,jj-1,jk) ) / ze2u & 221 & + ( zut(ji+1,jj,jk) - zut(ji,jj ,jk) ) / e1u(ji,jj) 222 223 zva = + ( zuf(ji,jj ,jk) - zuf(ji-1,jj,jk) ) / ze2v & 224 & + ( zut(ji,jj+1,jk) - zut(ji ,jj,jk) ) / e2v(ji,jj) 211 225 ! add it to the general momentum trends 212 226 ua(ji,jj,jk) = ua(ji,jj,jk) + zua
Note: See TracChangeset
for help on using the changeset viewer.