46 use,
intrinsic :: iso_c_binding, only : c_ptr, c_int
64 wa1_d, wa2_d, wa3_d, f_u_d, f_v_d, f_w_d, &
65 B_d, h1_d, mu, rho, n) &
66 bind(c, name=
'pnpn_prs_res_part1_hip')
67 use,
intrinsic :: iso_c_binding
70 type(c_ptr),
value :: ta1_d, ta2_d, ta3_d
71 type(c_ptr),
value :: wa1_d, wa2_d, wa3_d
72 type(c_ptr),
value :: f_u_d, f_v_d, f_w_d
73 type(c_ptr),
value :: B_d, h1_d
81 bind(c, name=
'pnpn_prs_res_part2_hip')
82 use,
intrinsic :: iso_c_binding
84 type(c_ptr),
value :: p_res_d, wa1_d, wa2_d, wa3_d
91 bind(c, name=
'pnpn_prs_res_part3_hip')
92 use,
intrinsic :: iso_c_binding
95 type(c_ptr),
value :: p_res_d, ta1_d, ta2_d, ta3_d
103 ta1_d, ta2_d, ta3_d, f_u_d, f_v_d, f_w_d, n) &
104 bind(c, name=
'pnpn_vel_res_update_hip')
105 use,
intrinsic :: iso_c_binding
107 type(c_ptr),
value :: u_res_d, v_res_d, w_res_d
108 type(c_ptr),
value :: ta1_d, ta2_d, ta3_d
109 type(c_ptr),
value :: f_u_d, f_v_d, f_w_d
116 wa1_d, wa2_d, wa3_d, f_u_d, f_v_d, f_w_d, &
117 B_d, h1_d, mu, rho, n) &
118 bind(c, name=
'pnpn_prs_res_part1_cuda')
119 use,
intrinsic :: iso_c_binding
122 type(c_ptr),
value :: ta1_d, ta2_d, ta3_d
123 type(c_ptr),
value :: wa1_d, wa2_d, wa3_d
124 type(c_ptr),
value :: f_u_d, f_v_d, f_w_d
125 type(c_ptr),
value :: B_d, h1_d
126 real(c_rp) :: mu, rho
133 bind(c, name=
'pnpn_prs_res_part2_cuda')
134 use,
intrinsic :: iso_c_binding
136 type(c_ptr),
value :: p_res_d, wa1_d, wa2_d, wa3_d
143 bind(c, name=
'pnpn_prs_res_part3_cuda')
144 use,
intrinsic :: iso_c_binding
147 type(c_ptr),
value :: p_res_d, ta1_d, ta2_d, ta3_d
155 ta1_d, ta2_d, ta3_d, f_u_d, f_v_d, f_w_d, n) &
156 bind(c, name=
'pnpn_vel_res_update_cuda')
157 use,
intrinsic :: iso_c_binding
159 type(c_ptr),
value :: u_res_d, v_res_d, w_res_d
160 type(c_ptr),
value :: ta1_d, ta2_d, ta3_d
161 type(c_ptr),
value :: f_u_d, f_v_d, f_w_d
168 wa1_d, wa2_d, wa3_d, f_u_d, f_v_d, f_w_d, &
169 B_d, h1_d, mu, rho, n) &
170 bind(c, name=
'pnpn_prs_res_part1_opencl')
171 use,
intrinsic :: iso_c_binding
174 type(c_ptr),
value :: ta1_d, ta2_d, ta3_d
175 type(c_ptr),
value :: wa1_d, wa2_d, wa3_d
176 type(c_ptr),
value :: f_u_d, f_v_d, f_w_d
177 type(c_ptr),
value :: B_d, h1_d
178 real(c_rp) :: mu, rho
185 bind(c, name=
'pnpn_prs_res_part2_opencl')
186 use,
intrinsic :: iso_c_binding
188 type(c_ptr),
value :: p_res_d, wa1_d, wa2_d, wa3_d
195 bind(c, name=
'pnpn_prs_res_part3_opencl')
196 use,
intrinsic :: iso_c_binding
199 type(c_ptr),
value :: p_res_d, ta1_d, ta2_d, ta3_d
207 ta1_d, ta2_d, ta3_d, f_u_d, f_v_d, f_w_d, n) &
208 bind(c, name=
'pnpn_vel_res_update_opencl')
209 use,
intrinsic :: iso_c_binding
211 type(c_ptr),
value :: u_res_d, v_res_d, w_res_d
212 type(c_ptr),
value :: ta1_d, ta2_d, ta3_d
213 type(c_ptr),
value :: f_u_d, f_v_d, f_w_d
223 f_x, f_y, f_z, c_Xh, gs_Xh, bc_prs_surface, bc_sym_surface, Ax, bd, dt,&
225 type(
field_t),
intent(inout) :: p, u, v, w
226 type(
field_t),
intent(inout) :: u_e, v_e, w_e
227 type(
field_t),
intent(inout) :: p_res
228 type(
field_t),
intent(inout) :: f_x, f_y, f_z
229 type(
coef_t),
intent(inout) :: c_Xh
230 type(
gs_t),
intent(inout) :: gs_Xh
233 class(
ax_t),
intent(inout) :: Ax
234 real(kind=
rp),
intent(inout) :: bd
235 real(kind=
rp),
intent(in) :: dt
236 real(kind=
rp),
intent(in) :: mu
237 real(kind=
rp),
intent(in) :: rho
238 real(kind=
rp) :: dtbd
240 type(
field_t),
pointer :: ta1, ta2, ta3, wa1, wa2, wa3, work1, work2
241 integer :: temp_indices(8)
255 call curl(ta1, ta2, ta3, u_e, v_e, w_e, work1, work2, c_xh)
256 call curl(wa1, wa2, wa3, ta1, ta2, ta3, work1, work2, c_xh)
261 wa1%x_d, wa2%x_d, wa3%x_d, f_x%x_d, f_y%x_d, f_z%x_d, &
262 c_xh%B_d, c_xh%h1_d, mu, rho, n)
266 wa1%x_d, wa2%x_d, wa3%x_d, f_x%x_d, f_y%x_d, f_z%x_d, &
267 c_xh%B_d, c_xh%h1_d, mu, rho, n)
270 wa1%x_d, wa2%x_d, wa3%x_d, f_x%x_d, f_z%x_d, f_z%x_d, &
271 c_xh%B_d, c_xh%h1_d, mu, rho, n)
276 call gs_xh%op(ta1, gs_op_add)
277 call gs_xh%op(ta2, gs_op_add)
278 call gs_xh%op(ta3, gs_op_add)
280 call device_opcolv(ta1%x_d, ta2%x_d, ta3%x_d, c_xh%Binv_d, gdim, n)
282 call cdtp(wa1%x, ta1%x, c_xh%drdx, c_xh%dsdx, c_xh%dtdx, c_xh)
283 call cdtp(wa2%x, ta2%x, c_xh%drdy, c_xh%dsdy, c_xh%dtdy, c_xh)
284 call cdtp(wa3%x, ta3%x, c_xh%drdz, c_xh%dsdz, c_xh%dtdz, c_xh)
286 call ax%compute(p_res%x,p%x,c_xh,p%msh,p%Xh)
303 call bc_sym_surface%apply_surfvec_dev(wa1%x_d, wa2%x_d, wa3%x_d, ta1%x_d, ta2%x_d, ta3%x_d)
319 call bc_prs_surface%apply_surfvec_dev(ta1%x_d, ta2%x_d, ta3%x_d, &
335 p, f_x, f_y, f_z, c_Xh, msh, Xh, mu, rho, bd, dt, n)
336 class(
ax_t),
intent(in) :: Ax
337 type(
mesh_t),
intent(inout) :: msh
338 type(
space_t),
intent(inout) :: Xh
339 type(
field_t),
intent(inout) :: p, u, v, w
340 type(
field_t),
intent(inout) :: u_res, v_res, w_res
341 type(
field_t),
intent(inout) :: f_x, f_y, f_z
342 type(
coef_t),
intent(inout) :: c_Xh
343 real(kind=
rp),
intent(in) :: mu
344 real(kind=
rp),
intent(in) :: rho
345 real(kind=
rp),
intent(in) :: bd
346 real(kind=
rp),
intent(in) :: dt
347 integer,
intent(in) :: n
348 integer :: temp_indices(3)
349 type(
field_t),
pointer :: ta1, ta2, ta3
355 call ax%compute(u_res%x, u%x, c_xh, msh, xh)
356 call ax%compute(v_res%x, v%x, c_xh, msh, xh)
357 call ax%compute(w_res%x, w%x, c_xh, msh, xh)
363 call opgrad(ta1%x, ta2%x, ta3%x, p%x, c_xh)
367 ta1%x_d, ta2%x_d, ta3%x_d, f_x%x_d, f_y%x_d, f_z%x_d, n)
370 ta1%x_d, ta2%x_d, ta3%x_d, f_x%x_d, f_y%x_d, f_z%x_d, n)
373 ta1%x_d, ta2%x_d, ta3%x_d, f_x%x_d, f_y%x_d, f_z%x_d, n)
Defines a Matrix-vector product.
subroutine, public device_rzero(a_d, n)
subroutine, public device_cfill(a_d, c, n)
subroutine, public device_opcolv(a1_d, a2_d, a3_d, c_d, gdim, n)
Dirichlet condition applied in the facet normal direction.
integer, parameter, public c_rp
integer, parameter, public rp
Global precision used in computations.
subroutine, public cdtp(dtx, x, dr, ds, dt, coef)
Apply D^T to a scalar field, where D is the derivative matrix.
subroutine, public opgrad(ux, uy, uz, u, coef, es, ee)
Compute the gradient of a scalar field.
subroutine, public curl(w1, w2, w3, u1, u2, u3, work1, work2, coef)
subroutine pnpn_vel_res_device_compute(Ax, u, v, w, u_res, v_res, w_res, p, f_x, f_y, f_z, c_Xh, msh, Xh, mu, rho, bd, dt, n)
subroutine pnpn_prs_res_device_compute(p, p_res, u, v, w, u_e, v_e, w_e, f_x, f_y, f_z, c_Xh, gs_Xh, bc_prs_surface, bc_sym_surface, Ax, bd, dt, mu, rho)
Defines Pressure and velocity residuals in the Pn-Pn formulation.
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.
Defines a function space.
void pnpn_prs_res_part3_opencl(void *p_res, void *ta1, void *ta2, void *ta3, real *dtbd, int *n)
void pnpn_prs_res_part1_opencl(void *ta1, void *ta2, void *ta3, void *wa1, void *wa2, void *wa3, void *f_u, void *f_v, void *f_w, void *B, void *h1, real *mu, real *rho, int *n)
void pnpn_prs_res_part2_opencl(void *p_res, void *wa1, void *wa2, void *wa3, int *n)
void pnpn_vel_res_update_opencl(void *u_res, void *v_res, void *w_res, void *ta1, void *ta2, void *ta3, void *f_u, void *f_v, void *f_w, int *n)
void pnpn_prs_res_part2_cuda(void *p_res, void *wa1, void *wa2, void *wa3, int *n)
void pnpn_prs_res_part3_cuda(void *p_res, void *ta1, void *ta2, void *ta3, real *dtbd, int *n)
void pnpn_vel_res_update_cuda(void *u_res, void *v_res, void *w_res, void *ta1, void *ta2, void *ta3, void *f_u, void *f_v, void *f_w, int *n)
void pnpn_prs_res_part1_cuda(void *ta1, void *ta2, void *ta3, void *wa1, void *wa2, void *wa3, void *f_u, void *f_v, void *f_w, void *B, void *h1, real *mu, real *rho, int *n)
Base type for a matrix-vector product providing .
Coefficients defined on a given (mesh, ) tuple. Arrays use indices (i,j,k,e): element e,...
Dirichlet condition in facet normal direction.
Abstract type to compute pressure residual.
Abstract type to compute velocity residual.
The function space for the SEM solution fields.