-
Notifications
You must be signed in to change notification settings - Fork 2
/
press_rhs.f90
103 lines (78 loc) · 2.05 KB
/
press_rhs.f90
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
subroutine press_rhs
! right-hand-side of the Poisson equation for pressure
use vars
use params, only: dowallx, dowally
implicit none
real *8 dta,rdx,rdy,rdz,btat,ctat,rup,rdn
integer i,j,k,ic,jc,kc
if(dowallx.and.mod(rank,nsubdomains_x).eq.0) then
do k=1,nzm
do j=1,ny
dudt(1,j,k,na) = 0.
end do
end do
end if
if(dowally.and.RUN3D.and.rank.lt.nsubdomains_x) then
do k=1,nzm
do i=1,nx
dvdt(i,1,k,na) = 0.
end do
end do
end if
if(dompi) then
call task_bound_duvdt()
else
call bound_duvdt()
endif
dta=1./dt3(na)/at
rdx=1./dx
rdy=1./dy
btat=bt/at
ctat=ct/at
if(RUN3D) then
do k=1,nzm
kc=k+1
rdz=1./(adz(k)*dz)
rup = rhow(kc)/rho(k)*rdz
rdn = rhow(k)/rho(k)*rdz
do j=1,ny
jc=j+1
do i=1,nx
ic=i+1
p(i,j,k)=(rdx*(u(ic,j,k)-u(i,j,k))+ &
rdy*(v(i,jc,k)-v(i,j,k))+ &
(w(i,j,kc)*rup-w(i,j,k)*rdn) )*dta + &
(rdx*(dudt(ic,j,k,na)-dudt(i,j,k,na))+ &
rdy*(dvdt(i,jc,k,na)-dvdt(i,j,k,na))+ &
(dwdt(i,j,kc,na)*rup-dwdt(i,j,k,na)*rdn) ) + &
btat*(rdx*(dudt(ic,j,k,nb)-dudt(i,j,k,nb))+ &
rdy*(dvdt(i,jc,k,nb)-dvdt(i,j,k,nb))+ &
(dwdt(i,j,kc,nb)*rup-dwdt(i,j,k,nb)*rdn) ) + &
ctat*(rdx*(dudt(ic,j,k,nc)-dudt(i,j,k,nc))+ &
rdy*(dvdt(i,jc,k,nc)-dvdt(i,j,k,nc))+ &
(dwdt(i,j,kc,nc)*rup-dwdt(i,j,k,nc)*rdn) )
end do
end do
end do
else
j=1
do k=1,nzm
kc=k+1
rdz=1./(adz(k)*dz)
rup = rhow(kc)/rho(k)*rdz
rdn = rhow(k)/rho(k)*rdz
do i=1,nx
ic=i+1
p(i,j,k)=(rdx*(u(ic,j,k)-u(i,j,k))+ &
(w(i,j,kc)*rup-w(i,j,k)*rdn) )*dta + &
(rdx*(dudt(ic,j,k,na)-dudt(i,j,k,na))+ &
(dwdt(i,j,kc,na)*rup-dwdt(i,j,k,na)*rdn) ) + &
btat*(rdx*(dudt(ic,j,k,nb)-dudt(i,j,k,nb))+ &
(dwdt(i,j,kc,nb)*rup-dwdt(i,j,k,nb)*rdn) ) + &
ctat*(rdx*(dudt(ic,j,k,nc)-dudt(i,j,k,nc))+ &
(dwdt(i,j,kc,nc)*rup-dwdt(i,j,k,nc)*rdn) )
end do
end do
endif
call task_barrier()
end subroutine press_rhs