15 use,
intrinsic :: iso_c_binding, only : c_ptr, c_int
39 wa1_d, wa2_d, wa3_d, f_u_d, f_v_d, f_w_d, &
40 B_d, h1_d, mu, rho, n) &
41 bind(c, name =
'pnpn_prs_res_part1_hip')
42 use,
intrinsic :: iso_c_binding
45 type(c_ptr),
value :: ta1_d, ta2_d, ta3_d
46 type(c_ptr),
value :: wa1_d, wa2_d, wa3_d
47 type(c_ptr),
value :: f_u_d, f_v_d, f_w_d
48 type(c_ptr),
value :: B_d, h1_d
56 bind(c, name =
'pnpn_prs_res_part2_hip')
57 use,
intrinsic :: iso_c_binding
59 type(c_ptr),
value :: p_res_d, wa1_d, wa2_d, wa3_d
66 bind(c, name =
'pnpn_prs_res_part3_hip')
67 use,
intrinsic :: iso_c_binding
70 type(c_ptr),
value :: p_res_d, ta1_d, ta2_d, ta3_d
78 ta1_d, ta2_d, ta3_d, f_u_d, f_v_d, f_w_d, n) &
79 bind(c, name =
'pnpn_vel_res_update_hip')
80 use,
intrinsic :: iso_c_binding
82 type(c_ptr),
value :: u_res_d, v_res_d, w_res_d
83 type(c_ptr),
value :: ta1_d, ta2_d, ta3_d
84 type(c_ptr),
value :: f_u_d, f_v_d, f_w_d
91 wa1_d, wa2_d, wa3_d, f_u_d, f_v_d, f_w_d, &
92 B_d, h1_d, rho_d, n) &
93 bind(c, name =
'pnpn_prs_stress_res_part1_cuda')
94 use,
intrinsic :: iso_c_binding
97 type(c_ptr),
value :: ta1_d, ta2_d, ta3_d
98 type(c_ptr),
value :: wa1_d, wa2_d, wa3_d
99 type(c_ptr),
value :: f_u_d, f_v_d, f_w_d
100 type(c_ptr),
value :: B_d, h1_d, rho_d
107 bind(c, name =
'pnpn_prs_res_part2_cuda')
108 use,
intrinsic :: iso_c_binding
110 type(c_ptr),
value :: p_res_d, wa1_d, wa2_d, wa3_d
117 wa1_d, wa2_d, wa3_d, dtbd, n) &
118 bind(c, name =
'pnpn_prs_stress_res_part3_cuda')
119 use,
intrinsic :: iso_c_binding
122 type(c_ptr),
value :: p_res_d, ta1_d, ta2_d, ta3_d
123 type(c_ptr),
value :: wa1_d, wa2_d, wa3_d
131 ta1_d, ta2_d, ta3_d, f_u_d, f_v_d, f_w_d, n) &
132 bind(c, name =
'pnpn_vel_res_update_cuda')
133 use,
intrinsic :: iso_c_binding
135 type(c_ptr),
value :: u_res_d, v_res_d, w_res_d
136 type(c_ptr),
value :: ta1_d, ta2_d, ta3_d
137 type(c_ptr),
value :: f_u_d, f_v_d, f_w_d
144 wa1_d, wa2_d, wa3_d, f_u_d, f_v_d, f_w_d, &
145 B_d, h1_d, mu, rho, n) &
146 bind(c, name =
'pnpn_prs_res_part1_opencl')
147 use,
intrinsic :: iso_c_binding
150 type(c_ptr),
value :: ta1_d, ta2_d, ta3_d
151 type(c_ptr),
value :: wa1_d, wa2_d, wa3_d
152 type(c_ptr),
value :: f_u_d, f_v_d, f_w_d
153 type(c_ptr),
value :: B_d, h1_d
154 real(c_rp) :: mu, rho
161 bind(c, name =
'pnpn_prs_res_part2_opencl')
162 use,
intrinsic :: iso_c_binding
164 type(c_ptr),
value :: p_res_d, wa1_d, wa2_d, wa3_d
171 wa1_d, wa2_d, wa3_d, dtbd, n) &
172 bind(c, name =
'pnpn_prs_res_part3_opencl')
173 use,
intrinsic :: iso_c_binding
176 type(c_ptr),
value :: p_res_d, ta1_d, ta2_d, ta3_d
177 type(c_ptr),
value :: wa1_d, wa2_d, wa3_d
185 ta1_d, ta2_d, ta3_d, f_u_d, f_v_d, f_w_d, n) &
186 bind(c, name =
'pnpn_vel_res_update_opencl')
187 use,
intrinsic :: iso_c_binding
189 type(c_ptr),
value :: u_res_d, v_res_d, w_res_d
190 type(c_ptr),
value :: ta1_d, ta2_d, ta3_d
191 type(c_ptr),
value :: f_u_d, f_v_d, f_w_d
200 w_e, f_x, f_y, f_z, c_Xh, gs_Xh, bc_prs_surface, bc_sym_surface, Ax, bd,&
202 type(
field_t),
intent(inout) :: p, u, v, w
203 type(
field_t),
intent(inout) :: u_e, v_e, w_e
204 type(
field_t),
intent(inout) :: p_res
205 type(
field_t),
intent(inout) :: f_x, f_y, f_z
206 type(
coef_t),
intent(inout) :: c_Xh
207 type(
gs_t),
intent(inout) :: gs_Xh
210 class(
ax_t),
intent(inout) :: Ax
211 real(kind=
rp),
intent(inout) :: bd
212 real(kind=
rp),
intent(in) :: dt
213 type(
field_t),
intent(in) :: mu
214 type(
field_t),
intent(in) :: rho
215 real(kind=
rp) :: dtbd
216 integer :: n, nelv, lxyz, gdim
219 type(
field_t),
pointer :: ta1, ta2, ta3, wa1, wa2, wa3, work1, work2, work3
220 type(
field_t),
pointer :: s11, s22, s33, s12, s13, s23
221 integer :: temp_indices(15)
253 call curl(ta1, ta2, ta3, u_e, v_e, w_e, work1, work2, c_xh)
254 call curl(wa1, wa2, wa3, ta1, ta2, ta3, work1, work2, c_xh)
262 call strain_rate(s11%x, s22%x, s33%x, s12%x, s13%x, s23%x, &
268 call dudxyz(ta1%x, mu%x, c_xh%drdx, c_xh%dsdx, c_xh%dtdx, c_xh)
269 call dudxyz(ta2%x, mu%x, c_xh%drdy, c_xh%dsdy, c_xh%dtdy, c_xh)
270 call dudxyz(ta3%x, mu%x, c_xh%drdz, c_xh%dsdz, c_xh%dtdz, c_xh)
277 call device_vdot3 (work1%x_d, ta1%x_d, ta2%x_d, ta3%x_d, &
278 s11%x_d, s12%x_d, s13%x_d, n)
280 call device_vdot3 (work2%x_d, ta1%x_d, ta2%x_d, ta3%x_d, &
281 s12%x_d, s22%x_d, s23%x_d, lxyz)
283 call device_vdot3 (work3%x_d, ta1%x_d, ta2%x_d, ta3%x_d, &
284 s13%x_d, s23%x_d, s33%x_d, lxyz)
296 wa1%x_d, wa2%x_d, wa3%x_d, f_x%x_d, f_y%x_d, f_z%x_d, &
297 c_xh%B_d, c_xh%h1_d, rho%x_d, n)
299 call neko_error(
'No device backend configured')
302 call gs_xh%op(ta1, gs_op_add)
303 call gs_xh%op(ta2, gs_op_add)
304 call gs_xh%op(ta3, gs_op_add)
306 call device_opcolv(ta1%x_d, ta2%x_d, ta3%x_d, c_xh%Binv_d, gdim, n)
309 call cdtp(wa1%x, ta1%x, c_xh%drdx, c_xh%dsdx, c_xh%dtdx, c_xh)
310 call cdtp(wa2%x, ta2%x, c_xh%drdy, c_xh%dsdy, c_xh%dtdy, c_xh)
311 call cdtp(wa3%x, ta3%x, c_xh%drdz, c_xh%dsdz, c_xh%dtdz, c_xh)
314 call ax%compute(p_res%x, p%x, c_xh, p%msh, p%Xh)
331 call bc_sym_surface%apply_surfvec_dev(wa1%x_d, wa2%x_d, wa3%x_d, &
332 ta1%x_d , ta2%x_d, ta3%x_d)
339 call bc_prs_surface%apply_surfvec_dev(ta1%x_d, ta2%x_d, ta3%x_d, &
344 wa1%x_d, wa2%x_d, wa3%x_d, dtbd, n)
346 call neko_error(
'No device backend configured')
354 w_res, p, f_x, f_y, f_z, c_Xh, msh, Xh, mu, rho, bd, dt, n)
355 class(
ax_t),
intent(in) :: Ax
356 type(
mesh_t),
intent(inout) :: msh
357 type(
space_t),
intent(inout) :: Xh
358 type(
field_t),
intent(inout) :: p, u, v, w
359 type(
field_t),
intent(inout) :: u_res, v_res, w_res
360 type(
field_t),
intent(inout) :: f_x, f_y, f_z
361 type(
coef_t),
intent(inout) :: c_Xh
362 type(
field_t),
intent(in) :: mu
363 type(
field_t),
intent(in) :: rho
364 real(kind=
rp),
intent(in) :: bd
365 real(kind=
rp),
intent(in) :: dt
366 real(kind=
rp) :: bddt
367 integer :: temp_indices(3)
368 type(
field_t),
pointer :: ta1, ta2, ta3
369 integer,
intent(in) :: n
381 call ax%compute_vector(u_res%x, v_res%x, w_res%x, u%x, v%x, w%x, c_xh,&
389 call opgrad(ta1%x, ta2%x, ta3%x, p%x, c_xh)
394 ta1%x_d, ta2%x_d, ta3%x_d, f_x%x_d, f_y%x_d, f_z%x_d, n)
397 ta1%x_d, ta2%x_d, ta3%x_d, f_x%x_d, f_y%x_d, f_z%x_d, n)
400 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_col2(a_d, b_d, n)
Vector multiplication .
subroutine, public device_rzero(a_d, n)
Zero a real vector.
subroutine, public device_invcol1(a_d, n)
Invert a vector .
subroutine, public device_vdot3(dot_d, u1_d, u2_d, u3_d, v1_d, v2_d, v3_d, n)
Compute a dot product (3-d version) assuming vector components etc.
subroutine, public device_cmult(a_d, c, n)
Multiplication by constant c .
subroutine, public device_copy(a_d, b_d, n)
Copy a vector .
subroutine, public device_sub2(a_d, b_d, n)
Vector substraction .
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 opgrad(ux, uy, uz, u, coef, es, ee)
Compute the weak gradient of a scalar field, i.e. the gradient multiplied by the mass matrix.
subroutine, public curl(w1, w2, w3, u1, u2, u3, work1, work2, coef)
subroutine, public strain_rate(s11, s22, s33, s12, s13, s23, u, v, w, coef)
Compute the strain rate tensor, i.e 0.5 * du_i/dx_j + du_j/dx_i.
subroutine, public cdtp(dtx, x, dr, ds, dt, coef, es, ee)
Apply D^T to a scalar field, where D is the derivative matrix.
subroutine, public dudxyz(du, u, dr, ds, dt, coef)
Compute derivative of a scalar field along a single direction.
Residuals in the Pn-Pn formulation (device version)
subroutine pnpn_vel_res_stress_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_stress_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_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_stress_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, void *rho, int *n)
void pnpn_prs_stress_res_part3_cuda(void *p_res, void *ta1, void *ta2, void *ta3, void *wa1, void *wa2, void *wa3, real *dtbd, 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.
Device implementation of the pressure residual for the PnPn fluid with full viscous stress formulatio...
Device implementation of the velocity residual for the PnPn fluid with full viscous stress formulatio...
Abstract type to compute pressure residual.
Abstract type to compute velocity residual.
The function space for the SEM solution fields.