-
Notifications
You must be signed in to change notification settings - Fork 2
/
press_grad.f90
106 lines (78 loc) · 1.81 KB
/
press_grad.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_grad
! pressure term of the momentum equations
use vars
use params, only: dowallx, dowally
implicit none
real *8 rdx,rdy,rdz
integer i,j,k,kb,jb,ib
real du(nx,ny,nz,3)
rdx=1./dx
rdy=1./dy
if(dostatis) then
do k=1,nzm
do j=1,ny
do i=1,nx
du(i,j,k,1)=dudt(i,j,k,na)
du(i,j,k,2)=dvdt(i,j,k,na)
du(i,j,k,3)=dwdt(i,j,k,na)
end do
end do
end do
endif
do k=1,nzm
kb=max(1,k-1)
rdz = 1./(dz*adzw(k))
do j=1,ny
jb=j-YES3D
do i=1,nx
ib=i-1
dudt(i,j,k,na)=dudt(i,j,k,na)-(p(i,j,k)-p(ib,j,k))*rdx
dvdt(i,j,k,na)=dvdt(i,j,k,na)-(p(i,j,k)-p(i,jb,k))*rdy
dwdt(i,j,k,na)=dwdt(i,j,k,na)-(p(i,j,k)-p(i,j,kb))*rdz
end do ! i
end do ! j
end do ! k
do k=1,nzm
do j=1-YES3D,ny !bloss: 0,n* fixes computation of dp/d* in stats.
do i=0,nx
p(i,j,k)=p(i,j,k)*rho(k) ! convert p'/rho to p'
end do
end do
end do
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
if(dostatis) then
do k=1,nzm
do j=1,ny
do i=1,nx
du(i,j,k,1)=dudt(i,j,k,na)-du(i,j,k,1)
du(i,j,k,2)=dvdt(i,j,k,na)-du(i,j,k,2)
du(i,j,k,3)=dwdt(i,j,k,na)-du(i,j,k,3)
end do
end do
end do
call stat_tke(du,tkelepress)
call stat_mom(du,momlepress)
call setvalue(twlepres,nzm,0.)
call setvalue(qwlepres,nzm,0.)
call stat_sw1(du,twlepres,qwlepres)
endif
call task_barrier()
end subroutine press_grad