Neko 1.99.1
A portable framework for high-order spectral element flow simulations
Loading...
Searching...
No Matches
rhs_maker_cpu.f90
Go to the documentation of this file.
5 use field, only : field_t
6 use num_types, only : rp, c_rp
8 implicit none
9 private
10
11 type, public, extends(rhs_maker_sumab_t) :: rhs_maker_sumab_cpu_t
12 contains
13 procedure, nopass :: compute_fluid => rhs_maker_sumab_cpu
15
16 type, public, extends(rhs_maker_ext_t) :: rhs_maker_ext_cpu_t
17 contains
18 procedure, nopass :: compute_fluid => rhs_maker_ext_cpu
19 procedure, nopass :: compute_scalar => scalar_rhs_maker_ext_cpu
20 end type rhs_maker_ext_cpu_t
21
22 type, public, extends(rhs_maker_bdf_t) :: rhs_maker_bdf_cpu_t
23 contains
24 procedure, nopass :: compute_fluid => rhs_maker_bdf_cpu
25 procedure, nopass :: compute_scalar => scalar_rhs_maker_bdf_cpu
26 end type rhs_maker_bdf_cpu_t
27
28 type, public, extends(rhs_maker_oifs_t) :: rhs_maker_oifs_cpu_t
29 contains
30 procedure, nopass :: compute_fluid => rhs_maker_oifs_cpu
31 procedure, nopass :: compute_scalar => scalar_rhs_maker_oifs_cpu
33
34contains
35
36 subroutine rhs_maker_sumab_cpu(u, v, w, uu, vv, ww, uulag, vvlag, wwlag, &
37 ab, nab)
38 type(field_t), intent(inout) :: u,v, w
39 type(field_t), intent(inout) :: uu, vv, ww
40 type(field_series_t), intent(inout) :: uulag, vvlag, wwlag
41 real(kind=rp), dimension(3), intent(in) :: ab
42 integer, intent(in) :: nab
43 integer :: i, n
44
45 n = uu%dof%size()
46
47 if (nab .eq. 3) then
48 do concurrent(i = 1:n)
49 u%x(i,1,1,1) = ab(1) * uu%x(i,1,1,1) + &
50 ab(2) * uulag%lf(1)%x(i,1,1,1) + ab(3) * uulag%lf(2)%x(i,1,1,1)
51 v%x(i,1,1,1) = ab(1) * vv%x(i,1,1,1) + &
52 ab(2) * vvlag%lf(1)%x(i,1,1,1) + ab(3) * vvlag%lf(2)%x(i,1,1,1)
53 w%x(i,1,1,1) = ab(1) * ww%x(i,1,1,1) + &
54 ab(2) * wwlag%lf(1)%x(i,1,1,1) + ab(3) * wwlag%lf(2)%x(i,1,1,1)
55 end do
56 else
57 do concurrent(i = 1:n)
58 u%x(i,1,1,1) = ab(1) * uu%x(i,1,1,1) + ab(2) * uulag%lf(1)%x(i,1,1,1)
59 v%x(i,1,1,1) = ab(1) * vv%x(i,1,1,1) + ab(2) * vvlag%lf(1)%x(i,1,1,1)
60 w%x(i,1,1,1) = ab(1) * ww%x(i,1,1,1) + ab(2) * wwlag%lf(1)%x(i,1,1,1)
61 end do
62 end if
63
64 end subroutine rhs_maker_sumab_cpu
65
66 subroutine rhs_maker_ext_cpu(fx_lag, fy_lag, fz_lag, &
67 fx_laglag, fy_laglag, fz_laglag, fx, fy, fz, &
68 rho, ext_coeffs, n)
69 type(field_t), intent(inout) :: fx_lag, fy_lag, fz_lag
70 type(field_t), intent(inout) :: fx_laglag, fy_laglag, fz_laglag
71 real(kind=rp), intent(in) :: rho, ext_coeffs(4)
72 integer, intent(in) :: n
73 real(kind=rp), intent(inout) :: fx(n), fy(n), fz(n)
74 integer :: i
75 type(field_t), pointer :: temp1, temp2, temp3
76 integer :: temp_indices(3)
77
78 call neko_scratch_registry%request_field(temp1, temp_indices(1))
79 call neko_scratch_registry%request_field(temp2, temp_indices(2))
80 call neko_scratch_registry%request_field(temp3, temp_indices(3))
81
82 do concurrent(i = 1:n)
83 temp1%x(i,1,1,1) = ext_coeffs(2) * fx_lag%x(i,1,1,1) + &
84 ext_coeffs(3) * fx_laglag%x(i,1,1,1)
85 temp2%x(i,1,1,1) = ext_coeffs(2) * fy_lag%x(i,1,1,1) + &
86 ext_coeffs(3) * fy_laglag%x(i,1,1,1)
87 temp3%x(i,1,1,1) = ext_coeffs(2) * fz_lag%x(i,1,1,1) + &
88 ext_coeffs(3) * fz_laglag%x(i,1,1,1)
89 end do
90
91 do concurrent(i = 1:n)
92 fx_laglag%x(i,1,1,1) = fx_lag%x(i,1,1,1)
93 fy_laglag%x(i,1,1,1) = fy_lag%x(i,1,1,1)
94 fz_laglag%x(i,1,1,1) = fz_lag%x(i,1,1,1)
95 fx_lag%x(i,1,1,1) = fx(i)
96 fy_lag%x(i,1,1,1) = fy(i)
97 fz_lag%x(i,1,1,1) = fz(i)
98 end do
99
100 do concurrent(i = 1:n)
101 fx(i) = (ext_coeffs(1) * fx(i) + temp1%x(i,1,1,1)) * rho
102 fy(i) = (ext_coeffs(1) * fy(i) + temp2%x(i,1,1,1)) * rho
103 fz(i) = (ext_coeffs(1) * fz(i) + temp3%x(i,1,1,1)) * rho
104 end do
105
106 call neko_scratch_registry%relinquish_field(temp_indices)
107
108 end subroutine rhs_maker_ext_cpu
109
110 subroutine scalar_rhs_maker_ext_cpu(fs_lag, fs_laglag, fs, rho, &
111 ext_coeffs, n)
112 type(field_t), intent(inout) :: fs_lag
113 type(field_t), intent(inout) :: fs_laglag
114 real(kind=rp), intent(in) :: rho, ext_coeffs(4)
115 integer, intent(in) :: n
116 real(kind=rp), intent(inout) :: fs(n)
117 integer :: i
118 type(field_t), pointer :: temp1
119 integer :: temp_index
120
121 call neko_scratch_registry%request_field(temp1, temp_index)
122
123 do concurrent(i = 1:n)
124 temp1%x(i,1,1,1) = ext_coeffs(2) * fs_lag%x(i,1,1,1) + &
125 ext_coeffs(3) * fs_laglag%x(i,1,1,1)
126 end do
127
128 do concurrent(i = 1:n)
129 fs_laglag%x(i,1,1,1) = fs_lag%x(i,1,1,1)
130 fs_lag%x(i,1,1,1) = fs(i)
131 end do
132
133 do concurrent(i = 1:n)
134 fs(i) = (ext_coeffs(1) * fs(i) + temp1%x(i,1,1,1)) * rho
135 end do
136
137 call neko_scratch_registry%relinquish_field(temp_index)
138 end subroutine scalar_rhs_maker_ext_cpu
139
140 subroutine rhs_maker_bdf_cpu(ulag, vlag, wlag, bfx, bfy, bfz, &
141 u, v, w, B, rho, dt, bd, nbd, n)
142 integer, intent(in) :: n, nbd
143 type(field_t), intent(in) :: u, v, w
144 type(field_series_t), intent(in) :: ulag, vlag, wlag
145 real(kind=rp), intent(inout) :: bfx(n), bfy(n), bfz(n)
146 real(kind=rp), intent(in) :: b(n)
147 real(kind=rp), intent(in) :: dt, rho, bd(4)
148 type(field_t), pointer :: tb1, tb2, tb3
149 integer :: temp_indices(3)
150 integer :: i, ilag
151
152 call neko_scratch_registry%request_field(tb1, temp_indices(1))
153 call neko_scratch_registry%request_field(tb2, temp_indices(2))
154 call neko_scratch_registry%request_field(tb3, temp_indices(3))
155
156 do concurrent(i = 1:n)
157 tb1%x(i,1,1,1) = u%x(i,1,1,1) * b(i) * bd(2)
158 tb2%x(i,1,1,1) = v%x(i,1,1,1) * b(i) * bd(2)
159 tb3%x(i,1,1,1) = w%x(i,1,1,1) * b(i) * bd(2)
160 end do
161
162 do ilag = 2, nbd
163 do concurrent(i = 1:n)
164 tb1%x(i,1,1,1) = tb1%x(i,1,1,1) + &
165 (ulag%lf(ilag-1)%x(i,1,1,1) * b(i) * bd(ilag+1))
166 tb2%x(i,1,1,1) = tb2%x(i,1,1,1) + &
167 (vlag%lf(ilag-1)%x(i,1,1,1) * b(i) * bd(ilag+1))
168 tb3%x(i,1,1,1) = tb3%x(i,1,1,1) + &
169 (wlag%lf(ilag-1)%x(i,1,1,1) * b(i) * bd(ilag+1))
170 end do
171 end do
172
173 do concurrent(i = 1:n)
174 bfx(i) = bfx(i) + tb1%x(i,1,1,1) * (rho / dt)
175 bfy(i) = bfy(i) + tb2%x(i,1,1,1) * (rho / dt)
176 bfz(i) = bfz(i) + tb3%x(i,1,1,1) * (rho / dt)
177 end do
178
179 call neko_scratch_registry%relinquish_field(temp_indices)
180
181 end subroutine rhs_maker_bdf_cpu
182
183 subroutine scalar_rhs_maker_bdf_cpu(s_lag, fs, s, B, rho, dt, bd, nbd, n)
184 integer, intent(in) :: n, nbd
185 type(field_t), intent(in) :: s
186 type(field_series_t), intent(in) :: s_lag
187 real(kind=rp), intent(inout) :: fs(n)
188 real(kind=rp), intent(in) :: b(n)
189 real(kind=rp), intent(in) :: dt, rho, bd(4)
190 integer :: i, ilag
191 type(field_t), pointer :: temp1
192 integer :: temp_indices
193
194 call neko_scratch_registry%request_field(temp1, temp_indices)
195
196 do concurrent(i = 1:n)
197 temp1%x(i,1,1,1) = s%x(i,1,1,1) * b(i) * bd(2)
198 end do
199
200 do ilag = 2, nbd
201 do concurrent(i = 1:n)
202 temp1%x(i,1,1,1) = temp1%x(i,1,1,1) + &
203 (s_lag%lf(ilag-1)%x(i,1,1,1) * b(i) * bd(ilag+1))
204 end do
205 end do
206
207 do concurrent(i = 1:n)
208 fs(i) = fs(i) + temp1%x(i,1,1,1) * (rho / dt)
209 end do
210
211 call neko_scratch_registry%relinquish_field(temp_indices)
212 end subroutine scalar_rhs_maker_bdf_cpu
213
214 subroutine rhs_maker_oifs_cpu(phi_x, phi_y, phi_z, bf_x, bf_y, bf_z, &
215 rho, dt, n)
216 real(kind=rp), intent(in) :: rho, dt
217 integer, intent(in) :: n
218 real(kind=rp), intent(inout) :: bf_x(n), bf_y(n), bf_z(n)
219 real(kind=rp), intent(inout) :: phi_x(n), phi_y(n), phi_z(n)
220 integer :: i
221
222 do concurrent(i = 1:n)
223 bf_x(i) = bf_x(i) + phi_x(i) * (rho / dt)
224 bf_y(i) = bf_y(i) + phi_y(i) * (rho / dt)
225 bf_z(i) = bf_z(i) + phi_z(i) * (rho / dt)
226 end do
227
228 end subroutine rhs_maker_oifs_cpu
229
230 subroutine scalar_rhs_maker_oifs_cpu(phi_s, bf_s, rho, dt, n)
231 real(kind=rp), intent(in) :: rho, dt
232 integer, intent(in) :: n
233 real(kind=rp), intent(inout) :: bf_s(n)
234 real(kind=rp), intent(inout) :: phi_s(n)
235 integer :: i
236
237 do concurrent(i = 1:n)
238 bf_s(i) = bf_s(i) + phi_s(i) * (rho / dt)
239 end do
240
241 end subroutine scalar_rhs_maker_oifs_cpu
242
243end module rhs_maker_cpu
244
Contains the field_serties_t type.
Defines a field.
Definition field.f90:34
integer, parameter, public c_rp
Definition num_types.f90:13
integer, parameter, public rp
Global precision used in computations.
Definition num_types.f90:12
subroutine scalar_rhs_maker_ext_cpu(fs_lag, fs_laglag, fs, rho, ext_coeffs, n)
subroutine rhs_maker_sumab_cpu(u, v, w, uu, vv, ww, uulag, vvlag, wwlag, ab, nab)
subroutine scalar_rhs_maker_oifs_cpu(phi_s, bf_s, rho, dt, n)
subroutine rhs_maker_bdf_cpu(ulag, vlag, wlag, bfx, bfy, bfz, u, v, w, b, rho, dt, bd, nbd, n)
subroutine rhs_maker_oifs_cpu(phi_x, phi_y, phi_z, bf_x, bf_y, bf_z, rho, dt, n)
subroutine scalar_rhs_maker_bdf_cpu(s_lag, fs, s, b, rho, dt, bd, nbd, n)
subroutine rhs_maker_ext_cpu(fx_lag, fy_lag, fz_lag, fx_laglag, fy_laglag, fz_laglag, fx, fy, fz, rho, ext_coeffs, n)
Routines to generate the right-hand sides for the convection-diffusion equation. Employs the EXT/BDF ...
Definition rhs_maker.f90:38
Defines a registry for storing and requesting temporary fields This can be used when you have a funct...
type(scratch_registry_t), target, public neko_scratch_registry
Global scratch registry.
Stores a series (sequence) of fields, logically connected to a base field, and arranged according to ...
Abstract type to add contributions to F from lagged BD terms.
Definition rhs_maker.f90:59
Abstract type to sum up contributions to kth order extrapolation scheme.
Definition rhs_maker.f90:52
Abstract type to add contributions of kth order OIFS scheme.
Definition rhs_maker.f90:66
Abstract type to compute extrapolated velocity field for the pressure equation.
Definition rhs_maker.f90:46