00001
00002
00003
00004 #include "cppdefs.h"
00005 # if defined SOLVE3D && defined SPONGE_VIS2 && !defined UV_VIS2
00006
00007 # define MIX_S_UV_SPG
00008 # if defined M3CLIMATOLOGY && !defined AGRIF
00009 # define CLIMAT_SPONGE
00010 # endif
00011
00012 # ifdef MIX_S_UV_SPG
00013
00014
00015
00016
00017
00018
00019
00020
00021 subroutine uv3dmix_spg (tile)
00022 implicit none
00023 integer tile, trd, omp_get_thread_num
00024 # include "param.h"
00025 # include "private_scratch.h"
00026 # include "compute_tile_bounds.h"
00027 trd=omp_get_thread_num()
00028 call uv3dmix_spg_tile (Istr,Iend,Jstr,Jend,
00029 & A2d(1,1,trd), A2d(1,2,trd),
00030 & A2d(1,3,trd), A2d(1,4,trd))
00031 return
00032 end
00033
00034 subroutine uv3dmix_spg_tile (Istr,Iend,Jstr,Jend, UFx,UFe,VFx,VFe)
00035
00036
00037
00038
00039
00040
00041
00042
00043
00044
00045
00046
00047
00048
00049
00050
00051
00052
00053
00054
00055
00056
00057
00058
00059
00060
00061
00062 implicit none
00063 integer Istr,Iend,Jstr,Jend, i,j,k, indx
00064 real UFe(PRIVATE_2D_SCRATCH_ARRAY),
00065 & UFx(PRIVATE_2D_SCRATCH_ARRAY), cff,
00066 & VFe(PRIVATE_2D_SCRATCH_ARRAY), cff1,
00067 & VFx(PRIVATE_2D_SCRATCH_ARRAY)
00068 # include "param.h"
00069 # include "scalars.h"
00070 # include "grid.h"
00071 # include "ocean3d.h"
00072 # include "coupling.h"
00073 # include "mixing.h"
00074 # include "climat.h"
00075 # ifdef DIAGNOSTICS_UV
00076 # include "diagnostics.h"
00077 # endif
00078
00079 # include "compute_auxiliary_bounds.h"
00080
00081 indx=3-nstp
00082
00083
00084
00085
00086
00087 do k=1,N
00088 do j=JstrV-1,Jend
00089 do i=IstrU-1,Iend
00090 cff=Hz(i,j,k)*visc2_r(i,j)*
00091 # ifdef CLIMAT_SPONGE
00092 & ( on_r(i,j)*pm(i,j)*
00093 & ( pn_u(i+1,j)*(u(i+1,j,k,nstp)-uclm(i+1,j,k))
00094 & -pn_u(i ,j)*(u(i ,j,k,nstp)-uclm(i ,j,k)) )
00095 & -om_r(i,j)*pn(i,j)*
00096 & ( pm_v(i,j+1)*(v(i,j+1,k,nstp)-vclm(i,j+1,k))
00097 & -pm_v(i,j )*(v(i,j ,k,nstp)-vclm(i,j ,k)) ) )
00098 # else
00099 & ( on_r(i,j)*pm(i,j)*( pn_u(i+1,j)*u(i+1,j,k,nstp)
00100 & -pn_u(i ,j)*u(i ,j,k,nstp) )
00101 & -om_r(i,j)*pn(i,j)*( pm_v(i,j+1)*v(i,j+1,k,nstp)
00102 & -pm_v(i,j )*v(i,j ,k,nstp) ) )
00103 # endif
00104 UFx(i,j)=on_r(i,j)*on_r(i,j)*cff
00105 VFe(i,j)=om_r(i,j)*om_r(i,j)*cff
00106 enddo
00107 enddo
00108 do j=Jstr,Jend+1
00109 do i=Istr,Iend+1
00110 cff=0.0625*visc2_p(i,j)*
00111 & (Hz(i-1,j,k)+Hz(i,j,k)+Hz(i-1,j-1,k)+Hz(i,j-1,k))*
00112 # ifdef CLIMAT_SPONGE
00113 & ( (pm(i-1,j)+pm(i,j)+pm(i-1,j-1)+pm(i,j-1))*on_p(i,j)
00114 & *( pn_v(i ,j)*(v(i ,j,k,nstp)-vclm(i ,j,k))
00115 & -pn_v(i-1,j)*(v(i-1,j,k,nstp)-vclm(i-1,j,k)) )
00116 & +(pn(i-1,j)+pn(i,j)+pn(i-1,j-1)+pn(i,j-1))*om_p(i,j)
00117 & *( pm_u(i,j )*(u(i,j ,k,nstp)-uclm(i ,j,k))
00118 & -pm_u(i,j-1)*(u(i,j-1,k,nstp)-uclm(i-1,j,k)) ) )
00119 # else
00120 & ( (pm(i-1,j)+pm(i,j)+pm(i-1,j-1)+pm(i,j-1))*on_p(i,j)
00121 & *( pn_v(i ,j)*v(i ,j,k,nstp)
00122 & -pn_v(i-1,j)*v(i-1,j,k,nstp) )
00123 & +(pn(i-1,j)+pn(i,j)+pn(i-1,j-1)+pn(i,j-1))*om_p(i,j)
00124 & *( pm_u(i,j )*u(i,j ,k,nstp)
00125 & -pm_u(i,j-1)*u(i,j-1,k,nstp) ) )
00126 # endif
00127 # ifdef MASKING
00128 & *pmask(i,j)
00129 # endif
00130 UFe(i,j)=om_p(i,j)*om_p(i,j)*cff
00131 VFx(i,j)=on_p(i,j)*on_p(i,j)*cff
00132 enddo
00133 enddo
00134
00135
00136
00137
00138
00139
00140 do j=Jstr,Jend
00141 do i=IstrU,Iend
00142 cff=pn_u(i,j)*(UFx(i,j)-UFx(i-1,j))
00143 & +pm_u(i,j)*(UFe(i,j+1)-UFe(i,j))
00144 cff1=pm_u(i,j)*pn_u(i,j)*cff
00145 rufrc(i,j)=rufrc(i,j) + cff
00146 u(i,j,k,indx)=u(i,j,k,indx) + dt*cff1
00147 # ifdef DIAGNOSTICS_UV
00148 MHmix(i,j,k,1) = MHmix(i,j,k,1)+
00149 & 2*cff1/(Hz(i-1,j,k)+Hz(i,j,k))
00150 # endif
00151 enddo
00152 enddo
00153 do j=JstrV,Jend
00154 do i=Istr,Iend
00155 cff=pn_v(i,j)*(VFx(i+1,j)-VFx(i,j))
00156 & -pm_v(i,j)*(VFe(i,j)-VFe(i,j-1))
00157 cff1=pm_v(i,j)*pn_v(i,j)*cff
00158 rvfrc(i,j)=rvfrc(i,j) + cff
00159 v(i,j,k,indx)=v(i,j,k,indx) + dt*cff1
00160 # ifdef DIAGNOSTICS_UV
00161 MHmix(i,j,k,2) = MHmix(i,j,k,2)+
00162 & 2*cff1/(Hz(i,j,k)+Hz(i,j-1,k))
00163 # endif
00164 enddo
00165 enddo
00166 enddo
00167 return
00168 end
00169
00170 # else
00171
00172
00173
00174
00175
00176
00177
00178
00179 subroutine uv3dmix_spg (tile)
00180 implicit none
00181 integer tile, trd, omp_get_thread_num
00182 #include "param.h"
00183 #include "private_scratch.h"
00184 #include "compute_tile_bounds.h"
00185 trd=omp_get_thread_num()
00186 call uv3dmix_spg_tile (Istr,Iend,Jstr,Jend,
00187 & A2d(1, 1,trd), A2d(1, 2,trd), A2d(1, 3,trd), A2d(1,4,trd),
00188 & A2d(1, 5,trd), A2d(1, 7,trd), A2d(1, 9,trd), A2d(1,11,trd),
00189 & A2d(1,13,trd), A2d(1,15,trd), A2d(1,17,trd), A2d(1,19,trd),
00190 & A2d(1,21,trd), A2d(1,23,trd), A2d(1,25,trd), A2d(1,27,trd))
00191 return
00192 end
00193
00194 subroutine uv3dmix_spg_tile (Istr,Iend,Jstr,Jend, UFx, UFe,
00195 & VFx,VFe, UFs,VFs, dnUdx, dmUde,
00196 & dUdz, dnVdx, dmVde, dVdz,
00197 & dZdx_r, dZdx_p, dZde_r, dZde_p)
00198 implicit none
00199 integer Istr,Iend,Jstr,Jend, i,j,k, k1,k2, indx
00200 real UFe(PRIVATE_2D_SCRATCH_ARRAY),
00201 & VFe(PRIVATE_2D_SCRATCH_ARRAY), cff,
00202 & UFx(PRIVATE_2D_SCRATCH_ARRAY), cff1,
00203 & VFx(PRIVATE_2D_SCRATCH_ARRAY), cff2,
00204 & UFs(PRIVATE_2D_SCRATCH_ARRAY,2), cff3,
00205 & VFs(PRIVATE_2D_SCRATCH_ARRAY,2), cff4,
00206 & dmUde(PRIVATE_2D_SCRATCH_ARRAY,2), cff5,
00207 & dmVde(PRIVATE_2D_SCRATCH_ARRAY,2), cff6,
00208 & dnUdx(PRIVATE_2D_SCRATCH_ARRAY,2), cff7,
00209 & dnVdx(PRIVATE_2D_SCRATCH_ARRAY,2), cff8,
00210 & dUdz(PRIVATE_2D_SCRATCH_ARRAY,2),
00211 & dVdz(PRIVATE_2D_SCRATCH_ARRAY,2), dmUdz,
00212 & dZde_p(PRIVATE_2D_SCRATCH_ARRAY,2), dnUdz,
00213 & dZde_r(PRIVATE_2D_SCRATCH_ARRAY,2), dmVdz,
00214 & dZdx_p(PRIVATE_2D_SCRATCH_ARRAY,2), dnVdz,
00215 & dZdx_r(PRIVATE_2D_SCRATCH_ARRAY,2)
00216
00217 # include "param.h"
00218 # include "scalars.h"
00219 # include "grid.h"
00220 # include "ocean3d.h"
00221 # include "coupling.h"
00222 # include "mixing.h"
00223 # ifdef DIAGNOSTICS_UV
00224 # include "diagnostics.h"
00225 # endif
00226
00227 # include "compute_auxiliary_bounds.h"
00228
00229 indx=3-nstp
00230
00231 k2=1
00232 do k=0,N,+1
00233 k1=k2
00234 k2=3-k1
00235 if (k.lt.N) then
00236 do j=Jstr-1,Jend+1
00237 do i=IstrU-1,Iend+1
00238 UFx(i,j)=(z_r(i,j,k+1)-z_r(i-1,j,k+1))*pm_u(i,j)
00239 # ifdef MASKING
00240 & *umask(i,j)
00241 # endif
00242 enddo
00243 enddo
00244 do j=JstrV-1,Jend+1
00245 do i=Istr-1,Iend+1
00246 VFe(i,j)=(z_r(i,j,k+1)-z_r(i,j-1,k+1))*pn_v(i,j)
00247 # ifdef MASKING
00248 & *vmask(i,j)
00249 # endif
00250 enddo
00251 enddo
00252 do j=JstrV-1,Jend
00253 do i=IstrU-1,Iend
00254 dnUdx(i,j,k2)=pm(i,j)*
00255 # ifdef CLIMAT_SPONGE
00256 & ( pn_u(i+1,j)*(u(i+1,j,k+1,nstp)-uclm(i+1,j,k+1))
00257 & -pn_u(i ,j)*(u(i ,j,k+1,nstp)-uclm(i ,j,k+1)) )
00258 # else
00259 & (pn_u(i+1,j)*u(i+1,j,k+1,nstp)
00260 & -pn_u(i ,j)*u(i ,j,k+1,nstp))
00261 # endif
00262 # ifdef MASKING
00263 & *rmask(i,j)
00264 # endif
00265 dmVde(i,j,k2)=pn(i,j)*
00266 # ifdef CLIMAT_SPONGE
00267 & ( pm_v(i,j+1)*(v(i,j+1,k+1,nstp)-vclm(i,j+1,k+1))
00268 & -pm_v(i,j )*(v(i,j ,k+1,nstp)-vclm(i,j ,k+1)) )
00269 # else
00270 & (pm_v(i,j+1)*v(i,j+1,k+1,nstp)
00271 & -pm_v(i,j )*v(i,j ,k+1,nstp))
00272 # endif
00273 # ifdef MASKING
00274 & *rmask(i,j)
00275 # endif
00276 dZdx_r(i,j,k2)=0.5*(UFx(i,j)+UFx(i+1,j))
00277 dZde_r(i,j,k2)=0.5*(VFe(i,j)+VFe(i,j+1))
00278 enddo
00279 enddo
00280 do j=Jstr,Jend+1
00281 do i=Istr,Iend+1
00282 dmUde(i,j,k2)=0.25*(pn(i,j)+pn(i-1,j)+pn(i,j-1)
00283 & +pn(i-1,j-1))
00284 # ifdef CLIMAT_SPONGE
00285 & *( pm_u(i,j )*(u(i,j ,k+1,nstp)-uclm(i+1,j,k+1))
00286 & -pm_u(i,j-1)*(u(i,j-1,k+1,nstp)-uclm(i ,j,k+1)) )
00287 # else
00288 & *(pm_u(i,j )*u(i,j ,k+1,nstp)
00289 & -pm_u(i,j-1)*u(i,j-1,k+1,nstp))
00290 # endif
00291 # ifdef MASKING
00292 & *pmask(i,j)
00293 # endif
00294 dnVdx(i,j,k2)=0.25*(pm(i,j)+pm(i-1,j)+pm(i,j-1)
00295 & +pm(i-1,j-1))
00296 # ifdef CLIMAT_SPONGE
00297 & *( pn_v(i ,j)*(v(i ,j,k+1,nstp)-vclm(i ,j,k+1))
00298 & -pn_v(i-1,j)*(v(i-1,j,k+1,nstp)-vclm(i-1,j,k+1)) )
00299 # else
00300 & *(pn_v(i ,j)*v(i ,j,k+1,nstp)
00301 & -pn_v(i-1,j)*v(i-1,j,k+1,nstp))
00302 # endif
00303 # ifdef MASKING
00304 & *pmask(i,j)
00305 # endif
00306 dZde_p(i,j,k2)=0.5*(VFe(i-1,j)+VFe(i,j))
00307 dZdx_p(i,j,k2)=0.5*(UFx(i,j-1)+UFx(i,j))
00308 enddo
00309 enddo
00310 endif
00311
00312 if (k.eq.0 .or. k.eq.N) then
00313 do j=Jstr-1,Jend+1
00314 do i=IstrU-1,Iend+1
00315 dUdz(i,j,k2)=0.
00316 UFs(i,j,k2)=0.
00317 enddo
00318 enddo
00319 do j=JstrV-1,Jend+1
00320 do i=Istr-1,Iend+1
00321 dVdz(i,j,k2)=0.
00322 VFs(i,j,k2)=0.
00323 enddo
00324 enddo
00325 else
00326 do j=Jstr-1,Jend+1
00327 do i=IstrU-1,Iend+1
00328 # ifdef CLIMAT_SPONGE
00329 dUdz(i,j,k2)=2.*(u(i,j,k+1,nstp)-u(i,j,k,nstp)
00330 & -uclm(i,j,k+1)+uclm(i,j,k))
00331 # else
00332 dUdz(i,j,k2)=2.*(u(i,j,k+1,nstp)-u(i,j,k,nstp))
00333 # endif
00334 & /( z_r(i-1,j,k+1)-z_r(i-1,j,k)
00335 & +z_r(i ,j,k+1)-z_r(i ,j,k))
00336 enddo
00337 enddo
00338 do j=JstrV-1,Jend+1
00339 do i=Istr-1,Iend+1
00340 # ifdef CLIMAT_SPONGE
00341 dVdz(i,j,k2)=2.*(v(i,j,k+1,nstp)-v(i,j,k,nstp)
00342 & -vclm(i,j,k+1)+vclm(i,j,k))
00343 # else
00344 dVdz(i,j,k2)=2.*(v(i,j,k+1,nstp)-v(i,j,k,nstp))
00345 # endif
00346 & /( z_r(i,j-1,k+1)-z_r(i,j-1,k)
00347 & +z_r(i,j ,k+1)-z_r(i,j ,k))
00348 enddo
00349 enddo
00350 endif
00351
00352
00353
00354
00355 if (k.gt.0) then
00356 do j=JstrV-1,Jend
00357 do i=IstrU-1,Iend
00358 cff=visc2_r(i,j)*Hz(i,j,k)*(
00359 & om_r(i,j)*( dnUdx(i,j,k1) - 0.5*pn(i,j)*(
00360 & min(dZdx_r(i,j,k1),0.)*(dUdz(i,j,k1)+dUdz(i+1,j,k2))
00361 & +max(dZdx_r(i,j,k1),0.)*(dUdz(i,j,k2)+dUdz(i+1,j,k1))
00362 & ))
00363 & -om_r(i,j)*( dmVde(i,j,k1) - 0.5*pm(i,j)*(
00364 & min(dZde_r(i,j,k1),0.)*(dVdz(i,j,k1)+dVdz(i,j+1,k2))
00365 & +max(dZde_r(i,j,k1),0.)*(dVdz(i,j,k2)+dVdz(i,j+1,k1))
00366 & )))
00367 # ifdef MASKING
00368 & *rmask(i,j)
00369 # endif
00370 UFx(i,j)=om_r(i,j)*om_r(i,j)*cff
00371 VFe(i,j)=om_r(i,j)*om_r(i,j)*cff
00372 enddo
00373 enddo
00374 do j=Jstr,Jend+1
00375 do i=Istr,Iend+1
00376 cff=0.25*visc2_p(i,j)
00377 & *(Hz(i,j,k)+Hz(i-1,j,k)+Hz(i,j-1,k)+Hz(i-1,j-1,k))
00378 & *( om_p(i,j)*( dnVdx(i,j,k1)-0.125*( pn(i,j)+pn(i-1,j)
00379 & +pn(i,j-1)+pn(i-1,j-1))
00380 & *( min(dZdx_p(i,j,k1),0.)*(dVdz(i-1,j,k1)+dVdz(i,j,k2))
00381 & +max(dZdx_p(i,j,k1),0.)*(dVdz(i-1,j,k2)+dVdz(i,j,k1))
00382 & ))
00383 & + om_p(i,j)*( dmUde(i,j,k1)-0.125*( pm(i,j)+pm(i-1,j)
00384 & +pm(i,j-1)+pm(i-1,j-1))
00385 & *( min(dZde_p(i,j,k1),0.)*(dUdz(i,j-1,k1)+dUdz(i,j,k2))
00386 & +max(dZde_p(i,j,k1),0.)*(dUdz(i,j-1,k2)+dUdz(i,j,k1))
00387 & )))
00388 # ifdef MASKING
00389 & *pmask(i,j)
00390 # endif
00391 UFe(i,j)=om_p(i,j)*om_p(i,j)*cff
00392 VFx(i,j)=om_p(i,j)*om_p(i,j)*cff
00393 enddo
00394 enddo
00395
00396
00397
00398
00399 if (k.lt.N) then
00400 do j=Jstr,Jend
00401 do i=IstrU,Iend
00402 cff1=pn_u(i,j)
00403 cff2=pm_u(i,j)
00404 cff=0.25*( dVdz(i,j,k2)+dVdz(i-1,j,k2)
00405 & +dVdz(i,j+1,k2)+dVdz(i-1,j+1,k2))
00406 dnUdz=cff1*dUdz(i,j,k2)
00407 dmUdz=cff2*dUdz(i,j,k2)
00408 dnVdz=cff1*cff
00409 dmVdz=cff2*cff
00410
00411 cff1=min(dZdx_r(i-1,j,k1),0.)
00412 cff2=min(dZdx_r(i ,j,k2),0.)
00413 cff3=max(dZdx_r(i-1,j,k2),0.)
00414 cff4=max(dZdx_r(i ,j,k1),0.)
00415 cff5=min(dZde_r(i-1,j,k1),0.)
00416 cff6=min(dZde_r(i ,j,k2),0.)
00417 cff7=max(dZde_r(i-1,j,k2),0.)
00418 cff8=max(dZde_r(i ,j,k1),0.)
00419
00420 cff=om_u(i,j)*( cff1*(cff1*dnUdz-dnUdx(i-1,j,k1))
00421 & +cff2*(cff2*dnUdz-dnUdx(i ,j,k2))
00422 & +cff3*(cff3*dnUdz-dnUdx(i-1,j,k2))
00423 & +cff4*(cff4*dnUdz-dnUdx(i ,j,k1))
00424 & )
00425 & -om_u(i,j)*( cff1*(cff5*dmVdz-dmVde(i-1,j,k1))
00426 & +cff2*(cff6*dmVdz-dmVde(i ,j,k2))
00427 & +cff3*(cff7*dmVdz-dmVde(i-1,j,k2))
00428 & +cff4*(cff8*dmVdz-dmVde(i ,j,k1))
00429 & )
00430 cff1=min(dZde_p(i,j ,k1),0.)
00431 cff2=min(dZde_p(i,j+1,k2),0.)
00432 cff3=max(dZde_p(i,j ,k2),0.)
00433 cff4=max(dZde_p(i,j+1,k1),0.)
00434 cff5=min(dZdx_p(i,j ,k1),0.)
00435 cff6=min(dZdx_p(i,j+1,k2),0.)
00436 cff7=max(dZdx_p(i,j ,k2),0.)
00437 cff8=max(dZdx_p(i,j+1,k1),0.)
00438
00439 cff=cff + om_u(i,j)*(
00440 & cff1*(cff1*dmUdz-dmUde(i,j ,k1))
00441 & +cff2*(cff2*dmUdz-dmUde(i,j+1,k2))
00442 & +cff3*(cff3*dmUdz-dmUde(i,j ,k2))
00443 & +cff4*(cff4*dmUdz-dmUde(i,j+1,k1))
00444 & )
00445 & +om_u(i,j)*( cff1*(cff5*dnVdz-dnVdx(i,j ,k1))
00446 & +cff2*(cff6*dnVdz-dnVdx(i,j+1,k2))
00447 & +cff3*(cff7*dnVdz-dnVdx(i,j ,k2))
00448 & +cff4*(cff8*dnVdz-dnVdx(i,j+1,k1))
00449 & )
00450 UFs(i,j,k2)=0.25*(visc2_r(i-1,j)+visc2_r(i,j))*cff
00451 enddo
00452 enddo
00453
00454 do j=JstrV,Jend
00455 do i=Istr,Iend
00456 cff1=pn_v(i,j)
00457 cff2=pm_v(i,j)
00458 cff=0.25*( dUdz(i,j,k2)+dUdz(i+1,j,k2)
00459 & +dUdz(i,j-1,k2)+dUdz(i+1,j-1,k2))
00460 dnUdz=cff1*cff
00461 dmUdz=cff2*cff
00462 dnVdz=cff1*dVdz(i,j,k2)
00463 dmVdz=cff2*dVdz(i,j,k2)
00464
00465 cff1=min(dZdx_p(i ,j,k1),0.)
00466 cff2=min(dZdx_p(i+1,j,k2),0.)
00467 cff3=max(dZdx_p(i ,j,k2),0.)
00468 cff4=max(dZdx_p(i+1,j,k1),0.)
00469 cff5=min(dZde_p(i ,j,k1),0.)
00470 cff6=min(dZde_p(i+1,j,k2),0.)
00471 cff7=max(dZde_p(i ,j,k2),0.)
00472 cff8=max(dZde_p(i+1,j,k1),0.)
00473
00474 cff=om_v(i,j)*( cff1*(cff1*dnVdz-dnVdx(i ,j,k1))
00475 & +cff2*(cff2*dnVdz-dnVdx(i+1,j,k2))
00476 & +cff3*(cff3*dnVdz-dnVdx(i ,j,k2))
00477 & +cff4*(cff4*dnVdz-dnVdx(i+1,j,k1))
00478 & )
00479 & +om_v(i,j)*( cff1*(cff5*dmUdz-dmUde(i ,j,k1))
00480 & +cff2*(cff6*dmUdz-dmUde(i+1,j,k2))
00481 & +cff3*(cff7*dmUdz-dmUde(i ,j,k2))
00482 & +cff4*(cff8*dmUdz-dmUde(i+1,j,k1))
00483 & )
00484 cff1=min(dZde_r(i,j-1,k1),0.)
00485 cff2=min(dZde_r(i,j ,k2),0.)
00486 cff3=max(dZde_r(i,j-1,k2),0.)
00487 cff4=max(dZde_r(i,j ,k1),0.)
00488 cff5=min(dZdx_r(i,j-1,k1),0.)
00489 cff6=min(dZdx_r(i,j ,k2),0.)
00490 cff7=max(dZdx_r(i,j-1,k2),0.)
00491 cff8=max(dZdx_r(i,j ,k1),0.)
00492
00493 cff=cff+om_v(i,j)*(
00494 & cff1*(cff1*dmVdz-dmVde(i,j-1,k1))
00495 & +cff2*(cff2*dmVdz-dmVde(i,j ,k2))
00496 & +cff3*(cff3*dmVdz-dmVde(i,j-1,k2))
00497 & +cff4*(cff4*dmVdz-dmVde(i,j ,k1))
00498 & )
00499 & -om_v(i,j)*( cff1*(cff5*dnUdz-dnUdx(i,j-1,k1))
00500 & +cff2*(cff6*dnUdz-dnUdx(i,j ,k2))
00501 & +cff3*(cff7*dnUdz-dnUdx(i,j-1,k2))
00502 & +cff4*(cff8*dnUdz-dnUdx(i,j ,k1))
00503 & )
00504 VFs(i,j,k2)=0.25*(visc2_r(i,j-1)+visc2_r(i,j))*cff
00505 enddo
00506 enddo
00507 endif
00508
00509
00510
00511
00512
00513
00514 do j=Jstr,Jend
00515 do i=IstrU,Iend
00516 cff=pn_u(i,j)*(UFx(i,j)-UFx(i-1,j))
00517 & +pm_u(i,j)*(UFe(i,j+1)-UFe(i,j))
00518 cff1=pm_u(i,j)*pn_u(i,j)*cff
00519 rufrc(i,j)=rufrc(i,j) + cff
00520 u(i,j,k,indx)=u(i,j,k,indx) + dt*(cff1+UFs(i,j,k2)
00521 & -UFs(i,j,k1))
00522 # ifdef DIAGNOSTICS_UV
00523 MHmix(i,j,k,1) = MHmix(i,j,k,1)+
00524 & 2*(cff1+UFs(i,j,k2)-UFs(i,j,k1))
00525 & /(Hz(i-1,j,k)+Hz(i,j,k))
00526 # endif
00527 enddo
00528 enddo
00529 do j=JstrV,Jend
00530 do i=Istr,Iend
00531 cff=pn_v(i,j)*(VFx(i+1,j)-VFx(i,j))
00532 & -pm_v(i,j)*(VFe(i,j)-VFe(i,j-1))
00533 cff1=pm_v(i,j)*pn_v(i,j)*cff
00534 rvfrc(i,j)=rvfrc(i,j) + cff
00535 v(i,j,k,indx)=v(i,j,k,indx) + dt*(cff1+VFs(i,j,k2)
00536 & -VFs(i,j,k1))
00537 # ifdef DIAGNOSTICS_UV
00538 MHmix(i,j,k,2) = MHmix(i,j,k,2)+
00539 & 2*(cff1+VFs(i,j,k2)-VFs(i,j,k1))
00540 & /(Hz(i,j,k)+Hz(i,j-1,k))
00541 # endif
00542 enddo
00543 enddo
00544 endif
00545 enddo
00546 return
00547 end
00548
00549
00550 # endif /* MIX_S_UV_SPG */
00551
00552 #else
00553 subroutine uv3dmix_spg_empty
00554 end
00555 #endif
00556