38 type(
field_t),
intent(inout) :: u,v, w
39 type(
field_t),
intent(inout) :: uu, vv, ww
41 real(kind=
rp),
dimension(3),
intent(in) :: ab
42 integer,
intent(in) :: nab
48 u%x(i,1,1,1) = ab(1) * uu%x(i,1,1,1) + ab(2) * uulag%lf(1)%x(i,1,1,1)
49 v%x(i,1,1,1) = ab(1) * vv%x(i,1,1,1) + ab(2) * vvlag%lf(1)%x(i,1,1,1)
50 w%x(i,1,1,1) = ab(1) * ww%x(i,1,1,1) + ab(2) * wwlag%lf(1)%x(i,1,1,1)
55 u%x(i,1,1,1) = u%x(i,1,1,1) + ab(3) * uulag%lf(2)%x(i,1,1,1)
56 v%x(i,1,1,1) = v%x(i,1,1,1) + ab(3) * vvlag%lf(2)%x(i,1,1,1)
57 w%x(i,1,1,1) = w%x(i,1,1,1) + ab(3) * wwlag%lf(2)%x(i,1,1,1)
64 fx_laglag, fy_laglag, fz_laglag, fx, fy, fz, &
66 type(
field_t),
intent(inout) :: fx_lag, fy_lag, fz_lag
67 type(
field_t),
intent(inout) :: fx_laglag, fy_laglag, fz_laglag
68 real(kind=
rp),
intent(inout) :: rho, ext_coeffs(4)
69 integer,
intent(in) :: n
70 real(kind=
rp),
intent(inout) :: fx(n), fy(n), fz(n)
72 type(
field_t),
pointer :: temp1, temp2, temp3
73 integer :: temp_indices(3)
80 temp1%x(i,1,1,1) = ext_coeffs(2) * fx_lag%x(i,1,1,1) + &
81 ext_coeffs(3) * fx_laglag%x(i,1,1,1)
82 temp2%x(i,1,1,1) = ext_coeffs(2) * fy_lag%x(i,1,1,1) + &
83 ext_coeffs(3) * fy_laglag%x(i,1,1,1)
84 temp3%x(i,1,1,1) = ext_coeffs(2) * fz_lag%x(i,1,1,1) + &
85 ext_coeffs(3) * fz_laglag%x(i,1,1,1)
89 fx_laglag%x(i,1,1,1) = fx_lag%x(i,1,1,1)
90 fy_laglag%x(i,1,1,1) = fy_lag%x(i,1,1,1)
91 fz_laglag%x(i,1,1,1) = fz_lag%x(i,1,1,1)
92 fx_lag%x(i,1,1,1) = fx(i)
93 fy_lag%x(i,1,1,1) = fy(i)
94 fz_lag%x(i,1,1,1) = fz(i)
98 fx(i) = (ext_coeffs(1) * fx(i) + temp1%x(i,1,1,1)) * rho
99 fy(i) = (ext_coeffs(1) * fy(i) + temp2%x(i,1,1,1)) * rho
100 fz(i) = (ext_coeffs(1) * fz(i) + temp3%x(i,1,1,1)) * rho
109 type(
field_t),
intent(inout) :: fs_lag
110 type(
field_t),
intent(inout) :: fs_laglag
111 real(kind=
rp),
intent(inout) :: rho, ext_coeffs(4)
112 integer,
intent(in) :: n
113 real(kind=
rp),
intent(inout) :: fs(n)
115 type(
field_t),
pointer :: temp1
116 integer :: temp_index
121 temp1%x(i,1,1,1) = ext_coeffs(2) * fs_lag%x(i,1,1,1) + &
122 ext_coeffs(3) * fs_laglag%x(i,1,1,1)
126 fs_laglag%x(i,1,1,1) = fs_lag%x(i,1,1,1)
127 fs_lag%x(i,1,1,1) = fs(i)
131 fs(i) = (ext_coeffs(1) * fs(i) + temp1%x(i,1,1,1)) * rho
138 u, v, w, B, rho, dt, bd, nbd, n)
139 integer,
intent(in) :: n, nbd
140 type(
field_t),
intent(in) :: u, v, w
142 real(kind=
rp),
intent(inout) :: bfx(n), bfy(n), bfz(n)
143 real(kind=
rp),
intent(in) :: b(n)
144 real(kind=
rp),
intent(in) :: dt, rho, bd(4)
145 type(
field_t),
pointer :: tb1, tb2, tb3
146 type(
field_t),
pointer :: ta1, ta2, ta3
147 integer :: temp_indices(6)
158 tb1%x(i,1,1,1) = u%x(i,1,1,1) * b(i) * bd(2)
159 tb2%x(i,1,1,1) = v%x(i,1,1,1) * b(i) * bd(2)
160 tb3%x(i,1,1,1) = w%x(i,1,1,1) * b(i) * bd(2)
165 ta1%x(i,1,1,1) = ulag%lf(ilag-1)%x(i,1,1,1) * b(i) * bd(ilag+1)
166 ta2%x(i,1,1,1) = vlag%lf(ilag-1)%x(i,1,1,1) * b(i) * bd(ilag+1)
167 ta3%x(i,1,1,1) = wlag%lf(ilag-1)%x(i,1,1,1) * b(i) * bd(ilag+1)
171 tb1%x(i,1,1,1) = tb1%x(i,1,1,1) + ta1%x(i,1,1,1)
172 tb2%x(i,1,1,1) = tb2%x(i,1,1,1) + ta2%x(i,1,1,1)
173 tb3%x(i,1,1,1) = tb3%x(i,1,1,1) + ta3%x(i,1,1,1)
178 bfx(i) = bfx(i) + tb1%x(i,1,1,1) * (rho / dt)
179 bfy(i) = bfy(i) + tb2%x(i,1,1,1) * (rho / dt)
180 bfz(i) = bfz(i) + tb3%x(i,1,1,1) * (rho / dt)
188 integer,
intent(in) :: n, nbd
191 real(kind=
rp),
intent(inout) :: fs(n)
192 real(kind=
rp),
intent(in) :: b(n)
193 real(kind=
rp),
intent(in) :: dt, rho, bd(4)
195 type(
field_t),
pointer :: temp1, temp2
196 integer :: temp_indices(2)
202 temp2%x(i,1,1,1) = s%x(i,1,1,1) * b(i) * bd(2)
207 temp1%x(i,1,1,1) = s_lag%lf(ilag-1)%x(i,1,1,1) * b(i) * bd(ilag+1)
211 temp2%x(i,1,1,1) = temp2%x(i,1,1,1) + temp1%x(i,1,1,1)
216 fs(i) = fs(i) + temp2%x(i,1,1,1) * (rho / dt)
224 real(kind=
rp),
intent(in) :: rho, dt
225 integer,
intent(in) :: n
226 real(kind=
rp),
intent(inout) :: bf_x(n), bf_y(n), bf_z(n)
227 real(kind=
rp),
intent(inout) :: phi_x(n), phi_y(n), phi_z(n)
231 bf_x(i) = bf_x(i) + phi_x(i) * (rho / dt)
232 bf_y(i) = bf_y(i) + phi_y(i) * (rho / dt)
233 bf_z(i) = bf_z(i) + phi_z(i) * (rho / dt)
239 real(kind=
rp),
intent(in) :: rho, dt
240 integer,
intent(in) :: n
241 real(kind=
rp),
intent(inout) :: bf_s(n)
242 real(kind=
rp),
intent(inout) :: phi_s(n)
246 bf_s(i) = bf_s(i) + phi_s(i) * (rho / dt)
integer, parameter, public c_rp
integer, parameter, public rp
Global precision used in computations.
subroutine rhs_maker_bdf_sx(ulag, vlag, wlag, bfx, bfy, bfz, u, v, w, B, rho, dt, bd, nbd, n)
subroutine rhs_maker_oifs_sx(phi_x, phi_y, phi_z, bf_x, bf_y, bf_z, rho, dt, n)
subroutine scalar_rhs_maker_oifs_sx(phi_s, bf_s, rho, dt, n)
subroutine scalar_rhs_maker_bdf_sx(s_lag, fs, s, B, rho, dt, bd, nbd, n)
subroutine rhs_maker_ext_sx(fx_lag, fy_lag, fz_lag, fx_laglag, fy_laglag, fz_laglag, fx, fy, fz, rho, ext_coeffs, n)
subroutine rhs_maker_sumab_sx(u, v, w, uu, vv, ww, uulag, vvlag, wwlag, ab, nab)
subroutine scalar_rhs_maker_ext_sx(fs_lag, fs_laglag, fs, 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 add contributions of kth order OIFS scheme.
Abstract type to compute extrapolated velocity field for the pressure equation.