35 f_x, f_y, f_z, c_Xh, gs_Xh, bc_prs_surface, bc_sym_surface, Ax, bd, dt,&
37 type(
field_t),
intent(inout) :: p, u, v, w
38 type(
field_t),
intent(in) :: u_e, v_e, w_e
39 type(
field_t),
intent(inout) :: p_res
40 type(
field_t),
intent(in) :: f_x, f_y, f_z
41 type(
coef_t),
intent(inout) :: c_Xh
42 type(
gs_t),
intent(inout) :: gs_Xh
45 class(
ax_t),
intent(inout) :: Ax
46 real(kind=
rp),
intent(in) :: bd
47 real(kind=
rp),
intent(in) :: dt
49 type(
field_t),
intent(in) :: rho
51 integer :: n, nelv, lxyz
54 type(
field_t),
pointer :: ta1, ta2, ta3, wa1, wa2, wa3, work1, work2, work3
55 type(
field_t),
pointer :: s11, s22, s33, s12, s13, s23
56 integer :: temp_indices(15)
82 call rzero(c_xh%h2, n)
86 call curl(ta1, ta2, ta3, u_e, v_e, w_e, work1, work2, c_xh)
87 call curl(wa1, wa2, wa3, ta1, ta2, ta3, work1, work2, c_xh)
89 call col2(wa1%x, mu%x, n)
90 call col2(wa2%x, mu%x, n)
91 call col2(wa3%x, mu%x, n)
95 call strain_rate(s11%x, s22%x, s33%x, s12%x, s13%x, s23%x, &
100 call dudxyz(ta1%x, mu%x, c_xh%drdx, c_xh%dsdx, c_xh%dtdx, c_xh)
101 call dudxyz(ta2%x, mu%x, c_xh%drdy, c_xh%dsdy, c_xh%dtdy, c_xh)
102 call dudxyz(ta3%x, mu%x, c_xh%drdz, c_xh%dsdz, c_xh%dtdz, c_xh)
104 call cmult(ta1%x, 2.0_rp, n)
105 call cmult(ta2%x, 2.0_rp, n)
106 call cmult(ta3%x, 2.0_rp, n)
110 call vdot3(work1%x(:, :, :, e), &
111 ta1%x(:, :, :, e), ta2%x(:, :, :, e), ta3%x(:, :, :, e), &
112 s11%x(:, :, :, e), s12%x(:, :, :, e), s13%x(:, :, :, e), &
115 call vdot3 (work2%x(:, :, :, e), &
116 ta1%x(:, :, :, e), ta2%x(:, :, :, e), ta3%x(:, :, :, e), &
117 s12%x(:, :, :, e), s22%x(:, :, :, e), s23%x(:, :, :, e), &
120 call vdot3 (work3%x(:, :, :, e), &
121 ta1%x(:, :, :, e), ta2%x(:, :, :, e), ta3%x(:, :, :, e), &
122 s13%x(:, :, :, e), s23%x(:, :, :, e), s33%x(:, :, :, e), &
130 call sub2(wa1%x, work1%x, n)
131 call sub2(wa2%x, work2%x, n)
132 call sub2(wa3%x, work3%x, n)
134 do concurrent(i = 1:n)
135 ta1%x(i,1,1,1) = f_x%x(i,1,1,1) / rho%x(i,1,1,1) &
136 - ((wa1%x(i,1,1,1) / rho%x(i,1,1,1)) * c_xh%B(i,1,1,1))
137 ta2%x(i,1,1,1) = f_y%x(i,1,1,1) / rho%x(i,1,1,1) &
138 - ((wa2%x(i,1,1,1) / rho%x(i,1,1,1)) * c_xh%B(i,1,1,1))
139 ta3%x(i,1,1,1) = f_z%x(i,1,1,1) / rho%x(i,1,1,1) &
140 - ((wa3%x(i,1,1,1) / rho%x(i,1,1,1)) * c_xh%B(i,1,1,1))
143 call gs_xh%op(ta1, gs_op_add)
144 call gs_xh%op(ta2, gs_op_add)
145 call gs_xh%op(ta3, gs_op_add)
148 ta1%x(i,1,1,1) = ta1%x(i,1,1,1) * c_xh%Binv(i,1,1,1)
149 ta2%x(i,1,1,1) = ta2%x(i,1,1,1) * c_xh%Binv(i,1,1,1)
150 ta3%x(i,1,1,1) = ta3%x(i,1,1,1) * c_xh%Binv(i,1,1,1)
154 call cdtp(wa1%x, ta1%x, c_xh%drdx, c_xh%dsdx, c_xh%dtdx, c_xh)
155 call cdtp(wa2%x, ta2%x, c_xh%drdy, c_xh%dsdy, c_xh%dtdy, c_xh)
156 call cdtp(wa3%x, ta3%x, c_xh%drdz, c_xh%dsdz, c_xh%dtdz, c_xh)
159 call ax%compute(p_res%x, p%x, c_xh, p%msh, p%Xh)
162 p_res%x(i,1,1,1) = (-p_res%x(i,1,1,1)) &
163 + wa1%x(i,1,1,1) + wa2%x(i,1,1,1) + wa3%x(i,1,1,1)
170 wa1%x(i,1,1,1) = 0.0_rp
171 wa2%x(i,1,1,1) = 0.0_rp
172 wa3%x(i,1,1,1) = 0.0_rp
175 call bc_sym_surface%apply_surfvec(wa1%x, wa2%x, wa3%x, ta1%x, ta2%x, ta3%x,&
180 ta1%x(i,1,1,1) = 0.0_rp
181 ta2%x(i,1,1,1) = 0.0_rp
182 ta3%x(i,1,1,1) = 0.0_rp
185 call bc_prs_surface%apply_surfvec(ta1%x, ta2%x, ta3%x, u%x, v%x, w%x, n)
188 p_res%x(i,1,1,1) = p_res%x(i,1,1,1) &
189 - (dtbd * (ta1%x(i,1,1,1) + ta2%x(i,1,1,1) + ta3%x(i,1,1,1)))&
190 - (wa1%x(i,1,1,1) + wa2%x(i,1,1,1) + wa3%x(i,1,1,1))
198 p, f_x, f_y, f_z, c_Xh, msh, Xh, mu, rho, bd, dt, n)
199 class(
ax_t),
intent(in) :: Ax
200 type(
mesh_t),
intent(inout) :: msh
201 type(
space_t),
intent(inout) :: Xh
202 type(
field_t),
intent(inout) :: p, u, v, w
203 type(
field_t),
intent(inout) :: u_res, v_res, w_res
204 type(
field_t),
intent(in) :: f_x, f_y, f_z
205 type(
coef_t),
intent(inout) :: c_Xh
206 type(
field_t),
intent(in) :: mu
207 type(
field_t),
intent(in) :: rho
208 real(kind=
rp),
intent(in) :: bd
209 real(kind=
rp),
intent(in) :: dt
210 integer :: temp_indices(3)
211 type(
field_t),
pointer :: ta1, ta2, ta3
212 integer,
intent(in) :: n
215 call copy(c_xh%h1, mu%x, n)
216 call cmult2(c_xh%h2, rho%x, bd / dt, n)
220 call ax%compute_vector(u_res%x, v_res%x, w_res%x, u%x, v%x, w%x, c_xh,&
228 call opgrad(ta1%x, ta2%x, ta3%x, p%x, c_xh)
232 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)
233 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)
234 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)