-
Notifications
You must be signed in to change notification settings - Fork 69
/
Copy pathshifted_inflow.f90
103 lines (88 loc) · 3.64 KB
/
shifted_inflow.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
!!
!! Copyright (C) 2020 Johns Hopkins University
!!
!! This file is part of lesgo.
!!
!! lesgo is free software: you can redistribute it and/or modify
!! it under the terms of the GNU General Public License as published by
!! the Free Software Foundation, either version 3 of the License, or
!! (at your option) any later version.
!!
!! lesgo is distributed in the hope that it will be useful,
!! but WITHOUT ANY WARRANTY; without even the implied warranty of
!! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
!! GNU General Public License for more details.
!!
!! You should have received a copy of the GNU General Public License
!! along with lesgo. If not, see <http://www.gnu.org/licenses/>.
!!
!*******************************************************************************
module shifted_inflow
!*******************************************************************************
use types, only : rprec
use fringe
implicit none
private
public shifted_inflow_init, inflow_shifted
type(fringe_t) :: sample_fringe
type(fringe_t) :: apply_fringe
real(rprec), allocatable, dimension(:,:,:) :: u_s, v_s, w_s
#ifdef PPSCALARS
real(rprec), allocatable, dimension(:,:,:) :: theta_s
#endif
contains
!*******************************************************************************
subroutine shifted_inflow_init
!*******************************************************************************
use param, only : fringe_region_end, fringe_region_len, sampling_region_end
use param, only : shift_n, ny, nz
sample_fringe = fringe_t(sampling_region_end, fringe_region_len)
apply_fringe = fringe_t(fringe_region_end, fringe_region_len)
! Allocate the sample block
allocate(u_s(sample_fringe%nx, ny, nz ))
allocate(v_s(sample_fringe%nx, ny, nz ))
allocate(w_s(sample_fringe%nx, ny, nz ))
#ifdef PPSCALARS
allocate(theta_s(sample_fringe%nx, ny, nz))
#endif
! Only allow positive shifts than are less than ny/2
shift_n = modulo(abs(shift_n), ny)
if (shift_n > ny/2) shift_n = ny-shift_n
if (shift_n == 0) shift_n = 1
end subroutine shifted_inflow_init
!*******************************************************************************
subroutine inflow_shifted
!*******************************************************************************
use param, only : shift_n, ny, nz, coord
use sim_param, only : u, v, w
#ifdef PPSCALARS
use scalars, only : theta
#endif
integer :: i, i_w
! Sample and shift velocity
u_s(:,shift_n+1:ny,:) = u(sample_fringe%iwrap(:),1:ny-shift_n,1:nz)
u_s(:,1:shift_n,:) = u(sample_fringe%iwrap(:),ny-shift_n+1:ny,1:nz)
v_s(:,shift_n+1:ny,:) = v(sample_fringe%iwrap(:),1:ny-shift_n,1:nz)
v_s(:,1:shift_n,:) = v(sample_fringe%iwrap(:),ny-shift_n+1:ny,1:nz)
w_s(:,shift_n+1:ny,:) = w(sample_fringe%iwrap(:),1:ny-shift_n,1:nz)
w_s(:,1:shift_n,:) = w(sample_fringe%iwrap(:),ny-shift_n+1:ny,1:nz)
#ifdef PPSCALARS
theta_s(:,shift_n+1:ny,:) = theta(sample_fringe%iwrap(:),1:ny-shift_n,1:nz)
theta_s(:,1:shift_n,:) = theta(sample_fringe%iwrap(:),ny-shift_n+1:ny,1:nz)
#endif
! Apply inflow conditions
do i = 1, apply_fringe%nx
i_w = apply_fringe%iwrap(i)
u(i_w,1:ny,1:nz) = apply_fringe%alpha(i) * u(i_w,1:ny,1:nz) &
+ apply_fringe%beta(i) * u_s(i,1:ny,1:nz)
v(i_w,1:ny,1:nz) = apply_fringe%alpha(i) * v(i_w,1:ny,1:nz) &
+ apply_fringe%beta(i) * v_s(i,1:ny,1:nz)
w(i_w,1:ny,1:nz) = apply_fringe%alpha(i) * w(i_w,1:ny,1:nz) &
+ apply_fringe%beta(i) * w_s(i,1:ny,1:nz)
#ifdef PPSCALARS
theta(i_w,1:ny,1:nz) = apply_fringe%alpha(i) * theta(i_w,1:ny,1:nz) &
+ apply_fringe%beta(i) * theta_s(i,1:ny,1:nz)
#endif
end do
end subroutine inflow_shifted
end module shifted_inflow