31 f_y, f_z, c_Xh, gs_Xh, bc_prs_surface, bc_sym_surface, Ax, bd, dt, mu, &
33 type(
field_t),
intent(inout) :: p, u, v, w
34 type(
field_t),
intent(inout) :: u_e, v_e, w_e
35 type(
field_t),
intent(inout) :: p_res
36 type(
field_t),
intent(inout) :: f_x, f_y, f_z
37 type(
coef_t),
intent(inout) :: c_Xh
38 type(
gs_t),
intent(inout) :: gs_Xh
41 class(
ax_t),
intent(inout) :: Ax
42 real(kind=
rp),
intent(inout) :: bd
43 real(kind=
rp),
intent(in) :: dt
45 type(
field_t),
intent(in) :: rho
46 real(kind=
rp) :: dtbd, rho_val, mu_val
49 type(
field_t),
pointer :: ta1, ta2, ta3, wa1, wa2, wa3, work1, work2
50 integer :: temp_indices(8)
64 rho_val = rho%x(1,1,1,1)
65 mu_val = mu%x(1,1,1,1)
67 c_xh%h1(i,1,1,1) = 1.0_rp / rho_val
68 c_xh%h2(i,1,1,1) = 0.0_rp
72 call curl(ta1, ta2, ta3, u_e, v_e, w_e, work1, work2, c_xh)
73 call curl(wa1, wa2, wa3, ta1, ta2, ta3, work1, work2, c_xh)
76 wa1%x(i,1,1,1) = (wa1%x(i,1,1,1) * mu_val / rho_val) * c_xh%B(i,1,1,1)
77 wa2%x(i,1,1,1) = (wa2%x(i,1,1,1) * mu_val / rho_val) * c_xh%B(i,1,1,1)
78 wa3%x(i,1,1,1) = (wa3%x(i,1,1,1) * mu_val / rho_val) * c_xh%B(i,1,1,1)
82 ta1%x(i,1,1,1) = f_x%x(i,1,1,1) / rho_val - wa1%x(i,1,1,1)
83 ta2%x(i,1,1,1) = f_y%x(i,1,1,1) / rho_val - wa2%x(i,1,1,1)
84 ta3%x(i,1,1,1) = f_z%x(i,1,1,1) / rho_val - wa3%x(i,1,1,1)
87 call gs_xh%op(ta1, gs_op_add)
88 call gs_xh%op(ta2, gs_op_add)
89 call gs_xh%op(ta3, gs_op_add)
92 ta1%x(i,1,1,1) = ta1%x(i,1,1,1) * c_xh%Binv(i,1,1,1)
93 ta2%x(i,1,1,1) = ta2%x(i,1,1,1) * c_xh%Binv(i,1,1,1)
94 ta3%x(i,1,1,1) = ta3%x(i,1,1,1) * c_xh%Binv(i,1,1,1)
97 call ax%compute(p_res%x, p%x, c_xh, p%msh, p%Xh)
99 call cdtp(wa1%x, ta1%x, c_xh%drdx, c_xh%dsdx, c_xh%dtdx, c_xh)
100 call cdtp(wa2%x, ta2%x, c_xh%drdy, c_xh%dsdy, c_xh%dtdy, c_xh)
101 call cdtp(wa3%x, ta3%x, c_xh%drdz, c_xh%dsdz, c_xh%dtdz, c_xh)
103 p_res%x(i,1,1,1) = (-p_res%x(i,1,1,1)) &
104 + wa1%x(i,1,1,1) + wa2%x(i,1,1,1) + wa3%x(i,1,1,1)
111 wa1%x(i,1,1,1) = 0.0_rp
112 wa2%x(i,1,1,1) = 0.0_rp
113 wa3%x(i,1,1,1) = 0.0_rp
116 call bc_sym_surface%apply_surfvec(wa1%x, wa2%x, wa3%x, ta1%x, ta2%x, ta3%x,&
121 ta1%x(i,1,1,1) = 0.0_rp
122 ta2%x(i,1,1,1) = 0.0_rp
123 ta3%x(i,1,1,1) = 0.0_rp
126 call bc_prs_surface%apply_surfvec(ta1%x, ta2%x, ta3%x, u%x, v%x, w%x, n)
129 p_res%x(i,1,1,1) = p_res%x(i,1,1,1) &
130 - (dtbd * (ta1%x(i,1,1,1) + ta2%x(i,1,1,1) + ta3%x(i,1,1,1)))&
131 - (wa1%x(i,1,1,1) + wa2%x(i,1,1,1) + wa3%x(i,1,1,1))
138 p, f_x, f_y, f_z, c_Xh, msh, Xh, mu, rho, bd, dt, n)
139 class(
ax_t),
intent(in) :: Ax
140 type(
mesh_t),
intent(inout) :: msh
141 type(
space_t),
intent(inout) :: Xh
142 type(
field_t),
intent(inout) :: p, u, v, w
143 type(
field_t),
intent(inout) :: u_res, v_res, w_res
144 type(
field_t),
intent(inout) :: f_x, f_y, f_z
145 type(
coef_t),
intent(inout) :: c_Xh
146 type(
field_t),
intent(in) :: mu
147 type(
field_t),
intent(in) :: rho
148 real(kind=
rp),
intent(in) :: bd
149 real(kind=
rp),
intent(in) :: dt
150 integer,
intent(in) :: n
151 integer :: temp_indices(3)
152 type(
field_t),
pointer :: ta1, ta2, ta3
154 real(kind=
rp) :: rho_val, mu_val
156 rho_val = rho%x(1,1,1,1)
157 mu_val = mu%x(1,1,1,1)
160 c_xh%h1(i,1,1,1) = mu_val
161 c_xh%h2(i,1,1,1) = rho_val * bd / dt
166 call ax%compute(u_res%x, u%x, c_xh, msh, xh)
167 call ax%compute(v_res%x, v%x, c_xh, msh, xh)
168 call ax%compute(w_res%x, w%x, c_xh, msh, xh)
174 call opgrad(ta1%x, ta2%x, ta3%x, p%x, c_xh)
177 u_res%x(i,1,1,1) = (-u_res%x(i,1,1,1)) - ta1%x(i,1,1,1) + f_x%x(i,1,1,1)
178 v_res%x(i,1,1,1) = (-v_res%x(i,1,1,1)) - ta2%x(i,1,1,1) + f_y%x(i,1,1,1)
179 w_res%x(i,1,1,1) = (-w_res%x(i,1,1,1)) - ta3%x(i,1,1,1) + f_z%x(i,1,1,1)
Defines a Matrix-vector product.
Dirichlet condition applied in the facet normal direction.
subroutine, public cmult2(a, b, c, n)
Multiplication by constant c .
subroutine, public invers2(a, b, n)
Compute inverted vector .
subroutine, public copy(a, b, n)
Copy a vector .
subroutine, public rzero(a, n)
Zero a real vector.
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 cdtp(dtx, x, dr, ds, dt, coef, es, ee)
Apply D^T to a scalar field, where D is the derivative matrix.
Residuals in the Pn-Pn formulation (SX version)
subroutine pnpn_vel_res_sx_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_sx_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.
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.