37 use,
intrinsic :: iso_c_binding, only : c_ptr
52 real(kind=
rp),
allocatable,
private :: tau1_(:)
54 real(kind=
rp),
allocatable,
private :: tau2_(:)
66 procedure, pass(this) :: init_shear_stress => &
80 integer,
intent(in) :: n
81 real(kind=
rp),
intent(inout),
dimension(n) :: x
82 real(kind=
rp),
intent(in),
optional :: t
83 integer,
intent(in),
optional :: tstep
84 integer :: i, m, k, facet
88 call neko_error(
"The shear stress bc is not applicable to scalar fields.")
96 integer,
intent(in) :: n
97 real(kind=
rp),
intent(inout),
dimension(n) :: x
98 real(kind=
rp),
intent(inout),
dimension(n) :: y
99 real(kind=
rp),
intent(inout),
dimension(n) :: z
100 real(kind=
rp),
intent(in),
optional :: t
101 integer,
intent(in),
optional :: tstep
103 call this%neumann1%apply_scalar(x, n, t, tstep)
104 call this%neumann2%apply_scalar(z, n, t, tstep)
105 call this%dirichlet%apply_scalar(y, n, t, tstep)
114 real(kind=
rp),
intent(in),
optional :: t
115 integer,
intent(in),
optional :: tstep
117 call neko_error(
"shear_stress bc not implemented on the device")
128 real(kind=
rp),
intent(in),
optional :: t
129 integer,
intent(in),
optional :: tstep
131 call neko_error(
"shear_stress bc not implemented on the device")
141 type(
coef_t),
target,
intent(in) :: coef
143 call this%init_base(coef)
145 call this%dirichlet%init_base(coef)
146 call this%neumann1%init_base(coef)
147 call this%neumann2%init_base(coef)
149 call this%dirichlet%set_g(0.0_rp)
158 real(kind=
rp),
intent(in) :: tau1
159 real(kind=
rp),
intent(in) :: tau2
161 call this%neumann1%mark_facets(this%marked_facet)
162 call this%neumann2%mark_facets(this%marked_facet)
163 call this%dirichlet%mark_facets(this%marked_facet)
165 call this%dirichlet%finalize()
166 call this%neumann1%finalize_neumann(tau1)
167 call this%neumann2%finalize_neumann(tau2)
174 real(kind=
rp) :: tau1(this%msk(0))
182 real(kind=
rp) :: tau2(this%msk(0))
190 real(kind=
rp),
intent(in) :: tau1(this%msk(0))
191 real(kind=
rp),
intent(in) :: tau2(this%msk(0))
193 call copy(this%tau1_, tau1, this%msk(0))
194 call copy(this%tau2_, tau2, this%msk(0))
202 call this%dirichlet%free
203 call this%neumann1%free
204 call this%neumann2%free
__device__ void nonlinear_index(const int idx, const int lx, int *index)
Defines a boundary condition.
subroutine, public device_copy(a_d, b_d, n)
Copy a vector .
Defines a dirichlet boundary condition.
subroutine, public copy(a, b, n)
Copy a vector .
Defines a Neumann boundary condition.
integer, parameter, public rp
Global precision used in computations.
Defines a shear stress boundary condition for a vector field.
subroutine shear_stress_apply_vector_dev(this, x_d, y_d, z_d, t, tstep)
Boundary condition apply for a generic shear_stress condition to vectors x, y and z (device version)
subroutine shear_stress_finalize(this, tau1, tau2)
Finalize by allocating the stress arrays and marking the facets for the bc components.
subroutine shear_stress_free(this)
Destructor.
pure real(kind=rp) function, dimension(this%msk(0)) shear_stress_tau2(this)
Get the stress in the 2nd wall-parallel direction.
subroutine shear_stress_set_stress(this, tau1, tau2)
Set the shear stress components.
subroutine shear_stress_init_shear_stress(this, coef)
Constructor.
subroutine shear_stress_apply_scalar(this, x, n, t, tstep)
Apply shear stress for a scalar field x.
subroutine shear_stress_apply_vector(this, x, y, z, n, t, tstep)
Boundary condition apply for a generic shear_stress condition to vectors x, y and z.
pure real(kind=rp) function, dimension(this%msk(0)) shear_stress_tau1(this)
Get the stress in the 1st wall-parallel direction.
subroutine shear_stress_apply_scalar_dev(this, x_d, t, tstep)
Boundary condition apply for a generic shear_stress condition to a vector x (device version)
Base type for a boundary condition.
Coefficients defined on a given (mesh, ) tuple. Arrays use indices (i,j,k,e): element e,...
Generic Dirichlet boundary condition on .
Generic Neumann boundary condition. This sets the flux of the field to the chosen value.
A shear stress boundary condition.