29 subroutine rhs_maker_sumab_cpu(u, v, w, uu, vv, ww, uulag, vvlag, wwlag, ab, nab)
30 type(
field_t),
intent(inout) :: u,v, w
31 type(
field_t),
intent(inout) :: uu, vv, ww
33 real(kind=
rp),
dimension(3),
intent(in) :: ab
34 integer,
intent(in) :: nab
39 do concurrent(i = 1:n)
40 u%x(i,1,1,1) = ab(1) * uu%x(i,1,1,1) + ab(2) * uulag%lf(1)%x(i,1,1,1)
41 v%x(i,1,1,1) = ab(1) * vv%x(i,1,1,1) + ab(2) * vvlag%lf(1)%x(i,1,1,1)
42 w%x(i,1,1,1) = ab(1) * ww%x(i,1,1,1) + ab(2) * wwlag%lf(1)%x(i,1,1,1)
46 do concurrent(i = 1:n)
47 u%x(i,1,1,1) = u%x(i,1,1,1) + ab(3) * uulag%lf(2)%x(i,1,1,1)
48 v%x(i,1,1,1) = v%x(i,1,1,1) + ab(3) * vvlag%lf(2)%x(i,1,1,1)
49 w%x(i,1,1,1) = w%x(i,1,1,1) + ab(3) * wwlag%lf(2)%x(i,1,1,1)
56 fx_laglag, fy_laglag, fz_laglag, fx, fy, fz, &
58 type(
field_t),
intent(inout) :: fx_lag, fy_lag, fz_lag
59 type(
field_t),
intent(inout) :: fx_laglag, fy_laglag, fz_laglag
60 real(kind=
rp),
intent(inout) :: rho, ext_coeffs(4)
61 integer,
intent(in) :: n
62 real(kind=
rp),
intent(inout) :: fx(n), fy(n), fz(n)
64 type(
field_t),
pointer :: temp1, temp2, temp3
65 integer :: temp_indices(3)
71 do concurrent(i = 1:n)
72 temp1%x(i,1,1,1) = ext_coeffs(2) * fx_lag%x(i,1,1,1) + &
73 ext_coeffs(3) * fx_laglag%x(i,1,1,1)
74 temp2%x(i,1,1,1) = ext_coeffs(2) * fy_lag%x(i,1,1,1) + &
75 ext_coeffs(3) * fy_laglag%x(i,1,1,1)
76 temp3%x(i,1,1,1) = ext_coeffs(2) * fz_lag%x(i,1,1,1) + &
77 ext_coeffs(3) * fz_laglag%x(i,1,1,1)
80 do concurrent(i = 1:n)
81 fx_laglag%x(i,1,1,1) = fx_lag%x(i,1,1,1)
82 fy_laglag%x(i,1,1,1) = fy_lag%x(i,1,1,1)
83 fz_laglag%x(i,1,1,1) = fz_lag%x(i,1,1,1)
84 fx_lag%x(i,1,1,1) = fx(i)
85 fy_lag%x(i,1,1,1) = fy(i)
86 fz_lag%x(i,1,1,1) = fz(i)
89 do concurrent(i = 1:n)
90 fx(i) = (ext_coeffs(1) * fx(i) + temp1%x(i,1,1,1)) * rho
91 fy(i) = (ext_coeffs(1) * fy(i) + temp2%x(i,1,1,1)) * rho
92 fz(i) = (ext_coeffs(1) * fz(i) + temp3%x(i,1,1,1)) * rho
100 type(
field_t),
intent(inout) :: fs_lag
101 type(
field_t),
intent(inout) :: fs_laglag
102 real(kind=
rp),
intent(inout) :: rho, ext_coeffs(4)
103 integer,
intent(in) :: n
104 real(kind=
rp),
intent(inout) :: fs(n)
106 type(
field_t),
pointer :: temp1
107 integer :: temp_index
111 do concurrent(i = 1:n)
112 temp1%x(i,1,1,1) = ext_coeffs(2) * fs_lag%x(i,1,1,1) + &
113 ext_coeffs(3) * fs_laglag%x(i,1,1,1)
116 do concurrent(i = 1:n)
117 fs_laglag%x(i,1,1,1) = fs_lag%x(i,1,1,1)
118 fs_lag%x(i,1,1,1) = fs(i)
121 do concurrent(i = 1:n)
122 fs(i) = (ext_coeffs(1) * fs(i) + temp1%x(i,1,1,1)) * rho
129 u, v, w, B, rho, dt, bd, nbd, n)
130 integer,
intent(in) :: n, nbd
131 type(
field_t),
intent(in) :: u, v, w
133 real(kind=
rp),
intent(inout) :: bfx(n), bfy(n), bfz(n)
134 real(kind=
rp),
intent(in) :: b(n)
135 real(kind=
rp),
intent(in) :: dt, rho, bd(4)
136 type(
field_t),
pointer :: tb1, tb2, tb3
137 type(
field_t),
pointer :: ta1, ta2, ta3
138 integer :: temp_indices(6)
148 do concurrent(i = 1:n)
149 tb1%x(i,1,1,1) = u%x(i,1,1,1) * b(i) * bd(2)
150 tb2%x(i,1,1,1) = v%x(i,1,1,1) * b(i) * bd(2)
151 tb3%x(i,1,1,1) = w%x(i,1,1,1) * b(i) * bd(2)
155 do concurrent(i = 1:n)
156 ta1%x(i,1,1,1) = ulag%lf(ilag-1)%x(i,1,1,1) * b(i) * bd(ilag+1)
157 ta2%x(i,1,1,1) = vlag%lf(ilag-1)%x(i,1,1,1) * b(i) * bd(ilag+1)
158 ta3%x(i,1,1,1) = wlag%lf(ilag-1)%x(i,1,1,1) * b(i) * bd(ilag+1)
161 do concurrent(i = 1:n)
162 tb1%x(i,1,1,1) = tb1%x(i,1,1,1) + ta1%x(i,1,1,1)
163 tb2%x(i,1,1,1) = tb2%x(i,1,1,1) + ta2%x(i,1,1,1)
164 tb3%x(i,1,1,1) = tb3%x(i,1,1,1) + ta3%x(i,1,1,1)
168 do concurrent(i = 1:n)
169 bfx(i) = bfx(i) + tb1%x(i,1,1,1) * (rho / dt)
170 bfy(i) = bfy(i) + tb2%x(i,1,1,1) * (rho / dt)
171 bfz(i) = bfz(i) + tb3%x(i,1,1,1) * (rho / dt)
179 integer,
intent(in) :: n, nbd
182 real(kind=
rp),
intent(inout) :: fs(n)
183 real(kind=
rp),
intent(in) :: b(n)
184 real(kind=
rp),
intent(in) :: dt, rho, bd(4)
186 type(
field_t),
pointer :: temp1, temp2
187 integer :: temp_indices(2)
192 do concurrent(i = 1:n)
193 temp2%x(i,1,1,1) = s%x(i,1,1,1) * b(i) * bd(2)
197 do concurrent(i = 1:n)
198 temp1%x(i,1,1,1) = s_lag%lf(ilag-1)%x(i,1,1,1) * b(i) * bd(ilag+1)
201 do concurrent(i = 1:n)
202 temp2%x(i,1,1,1) = temp2%x(i,1,1,1) + temp1%x(i,1,1,1)
206 do concurrent(i = 1:n)
207 fs(i) = fs(i) + temp2%x(i,1,1,1) * (rho / dt)
integer, parameter, public c_rp
integer, parameter, public rp
Global precision used in computations.
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 rhs_maker_bdf_cpu(ulag, vlag, wlag, bfx, bfy, bfz, u, v, w, B, rho, dt, bd, nbd, 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 ...
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.
Abstract type to add contributions to F from lagged BD terms.
Abstract type to sum up contributions to kth order extrapolation scheme.
Abstract type to compute extrapolated velocity field for the pressure equation.