copyright (C) 2001 MSC-RPN COMM %%%MC2%%%
***s/r adveul
*
subroutine adveul ( ua,va,wa,ba,pa,hua,tra, 1,1
$ lminx,lmaxx,lminy,lmaxy,lnk )
implicit none
*
integer lminx, lmaxx, lminy, lmaxy, lnk
real
$ ua(lminx:lmaxx,lminy:lmaxy,*), va(lminx:lmaxx,lminy:lmaxy,*),
$ wa(lminx:lmaxx,lminy:lmaxy,*), ba(lminx:lmaxx,lminy:lmaxy,*),
$ pa(lminx:lmaxx,lminy:lmaxy,0:lnk),hua(lminx:lmaxx,lminy:lmaxy,*),
$ tra(lminx:lmaxx,lminy:lmaxy,lnk,*)
*
*
#include "grd.cdk"
#include "consdyn_8.cdk"
#include "yomdyn.cdk"
#include "yomdyn1.cdk"
#include "dynpar.cdk"
#include "dynmem.cdk"
#include "topo.cdk"
#include "nbcpu.cdk"
*
integer i, j, k, n, kp1, kkk, km1, id, jd, iff, jf
real*8 ken(minx:maxx,miny:maxy),
$ vor(minx:maxx,miny:maxy),
$ udw(minx:maxx,miny:maxy),
$ udb(minx:maxx,miny:maxy),
$ udh(minx:maxx,miny:maxy),
$ vdw(minx:maxx,miny:maxy),
$ vdb(minx:maxx,miny:maxy),
$ vdh(minx:maxx,miny:maxy),
$ wdu(minx:maxx,miny:maxy,gnk),
$ wdv(minx:maxx,miny:maxy,gnk),
$ wdp(minx:maxx,miny:maxy,gnk)
real w0(minx:maxx,miny:maxy,gnk)
real*8 p25, pt5, frtsvn, uu0b, vv0b, ubarz, vbarz, coe
parameter (p25=0.25d0, pt5=0.5d0, frtsvn=4.d0/7.d0 )
*
*----------------------------------------------------------------------
*
if (flextop) print*,'not a valid option yet'
if (flextop) stop
*
id = 1
jd = 1
iff = ldni-east
jf = ldnj-north
*
call sw2w3(w0,uu0,vv0,ww0,sbxy,gg1,gg2,g0wr,dhdt,
$ minx,maxx,miny,maxy,gnk)
*
do k = gnk, 0, -1
kp1 = min( k+1, gnk)
kkk = min( max(k,1), gnk-1)
km1 = max( k-1, 1 )
coe = pt5
ccccc encore en doute cccc
ccc if (k.eq.gnk-1 .or. k.eq.2 ) coe = frtsvn
ccccc encore en doute cccc
*
if (k.eq.gnk) then
do j = jd,jf
do i = id,iff
wdu(i,j,gnk) = 0.
wdv(i,j,gnk) = 0.
end do
end do
endif
*
* horizontal wind components
*
if (k.gt.0 .and. k.lt.gnk) then
*
do j = jd+south-1,jf
do i = id+west-1,iff
ken(i,j) = p25 * ( uu0(i+1,j,k)**2 + uu0(i,j,k)**2
$ + vv0(i,j+1,k)**2 + vv0(i,j,k)**2 )
end do
end do
*
do j = jd,jf+1
do i = id,iff+1
vor(i,j) = (vv0(i,j,k)-vv0(i-1,j,k))*odxu(1)
$ - (uu0(i,j,k)-uu0(i,j-1,k))*odyv(j)
end do
end do
*
do j = jd, jf
do i = id+west, iff
vv0b = p25 * ( vv0(i,j+1,k) + vv0(i-1,j+1,k)
$ + vv0(i,j ,k) + vv0(i-1,j ,k) )
wdu(i,j,k) = pt5 * ( w0(i,j,k) + w0(i-1,j,k) )
$ * ( uu0(i,j,k) - uu0(i,j,km1) )
ua(i,j,k) = - sby(i,j) *
$ ( ( ken(i,j) - ken(i-1,j) ) * odxu(1)
$ - ( vor(i,j) + vor(i,j+1) ) * pt5 * vv0b )
$ - pt5 * ( wdu(i,j,k+1) + wdu(i,j,k) )
end do
end do
*
do j = jd+south, jf
do i = id, iff
uu0b = p25 * ( uu0(i+1,j,k) + uu0(i+1,j-1,k)
$ + uu0(i ,j,k) + uu0(i ,j-1,k) )
wdv(i,j,k) = pt5 * ( w0(i,j,k) + w0(i,j-1,k) )
$ * ( vv0(i,j,k) - vv0(i,j,km1) )
va(i,j,k) = - sbx(i,j) *
$ ( ( ken(i,j) - ken(i,j-1) ) * odyv(j)
$ + ( vor(i,j) + vor(i+1,j) ) * pt5 * uu0b )
$ - pt5 * ( wdv(i,j,k+1) + wdv(i,j,k) )
end do
end do
*
endif
*
* vertical motion, buoyancy, humidity and tracers
*
if(k.gt.0) then
*
do j = jd,jf+1
do i = id,iff+1
ubarz = pt5 * ( uu0(i,j,kkk) + uu0(i,j,km1) )
udw(i,j) = ubarz *( ww0(i,j,k)- ww0(i-1,j,k) )*odxu(1)
udb(i,j) = ubarz *( bb0(i,j,k)- bb0(i-1,j,k) )*odxu(1)
udh(i,j) = ubarz *( hu0(i,j,k)- hu0(i-1,j,k) )*odxu(1)
vbarz = pt5 * ( vv0(i,j,kkk) + vv0(i,j,km1) )
vdw(i,j) = vbarz *( ww0(i,j,k)- ww0(i,j-1,k) )*odyv(j)
vdb(i,j) = vbarz *( bb0(i,j,k)- bb0(i,j-1,k) )*odyv(j)
vdh(i,j) = vbarz *( hu0(i,j,k)- hu0(i,j-1,k) )*odyv(j)
end do
end do
*
do j = jd,jf
do i = id,iff
wdp(i,j,k) = w0(i,j,k) * ( pp0(i,j,k) - pp0(i,j,k-1) )
wa(i,j,k) = - sbxy(i,j) * pt5 * ( udw(i+1,j) + udw(i,j)
$ + vdw(i,j+1) + vdw(i,j) )
$ - pt5 * w0(i,j,k) * ( ww0(i,j,kp1) - ww0(i,j,km1) )
ba(i,j,k) = - sbxy(i,j) * pt5 * ( udb(i+1,j) + udb(i,j)
$ + vdb(i,j+1) + vdb(i,j) )
$ - coe * w0(i,j,k) * ( bb0(i,j,kp1) - bb0(i,j,km1) )
hua(i,j,k) = - sbxy(i,j) * pt5 * ( udh(i+1,j) + udh(i,j)
$ + vdh(i,j+1) + vdh(i,j) )
$ - coe * w0(i,j,k) * ( hu0(i,j,kp1) - hu0(i,j,km1) )
end do
end do
*
do n=1,ntr
do j = jd,jf+1
do i = id,iff+1
ubarz = pt5 * ( uu0(i,j,kkk) + uu0(i,j,km1) )
udh(i,j) = ubarz *( tr0(i,j,k,n)- tr0(i-1,j,k,n) )*odxu(1)
vbarz = pt5 * ( vv0(i,j,kkk) + vv0(i,j,km1) )
vdh(i,j) = vbarz *( tr0(i,j,k,n)- tr0(i,j-1,k,n) )*odyv(j)
end do
end do
do j=jd,jf
do i=id,iff
tra(i,j,k,n) = - sbxy(i,j) * pt5 * ( udh(i+1,j) + udh(i,j)
$ + vdh(i,j+1) + vdh(i,j) )
$ - coe * w0(i,j,k) * ( tr0(i,j,kp1,n) - tr0(i,j,km1,n) )
end do
end do
end do
*
endif
*
* pressure (horizontal advection)
*
do j = jd,jf+1
do i = id,iff+1
udh(i,j) = uu0(i,j,kkk)*( pp0(i,j,k) - pp0(i-1,j,k))*odxu(1)
vdh(i,j) = vv0(i,j,kkk)*( pp0(i,j,k) - pp0(i,j-1,k))*odyv(j)
end do
end do
*
do j = jd,jf
do i = id,iff
pa(i,j,k) = - sbxy(i,j) * pt5 * ( udh(i+1,j) + udh(i,j)
$ + vdh(i,j+1) + vdh(i,j) )
end do
end do
*
end do
*
* thermodynamic variable and generalized pressure
*
do k = 1, gnk
if (k.eq.1) then
do j = jd,jf
do i = id,iff
ba(i,j,k) = ba(i,j,k) - gama_star * ( - wdp(i,j,k)
$ + pt5 * ( pa(i,j,k) + pa(i,j,k-1) ) )
end do
end do
else
do j = jd,jf
do i = id,iff
ba(i,j,k) = ba(i,j,k) - gama_star * ( - wdp(i,j,k)
$ + pt5 * ( pa(i,j,k) + pa(i,j,k-1) ) )
pa(i,j,k-1) = c2r_star * ( pa(i,j,k-1)
$ - pt5 * ( wdp(i,j,k) + wdp(i,j,k-1) ) )
end do
end do
endif
end do
*
*----------------------------------------------------------------------
return
end