41 use,
intrinsic :: iso_c_binding, only : c_ptr
65 type(
coef_t),
target,
intent(in) :: coef
68 real(kind=
rp) :: sx, sy, sz
69 real(kind=
rp),
parameter :: tol = 1d-3
76 call this%bc_x%init_base(this%coef)
77 call this%bc_y%init_base(this%coef)
78 call this%bc_z%init_base(this%coef)
80 associate(c => this%coef, nx => this%coef%nx, ny => this%coef%ny, &
82 bfp => this%marked_facet%array()
83 do i = 1, this%marked_facet%size()
94 sx = sx + abs(abs(nx(l, j, facet, el)) - 1d0)
95 sy = sy + abs(abs(ny(l, j, facet, el)) - 1d0)
96 sz = sz + abs(abs(nz(l, j, facet, el)) - 1d0)
100 do l = 2, c%Xh%lx - 1
101 do j = 2, c%Xh%lx - 1
102 sx = sx + abs(abs(nx(l, j, facet, el)) - 1d0)
103 sy = sy + abs(abs(ny(l, j, facet, el)) - 1d0)
104 sz = sz + abs(abs(nz(l, j, facet, el)) - 1d0)
108 do l = 2, c%Xh%lx - 1
109 do j = 2, c%Xh%lx - 1
110 sx = sx + abs(abs(nx(l, j, facet, el)) - 1d0)
111 sy = sy + abs(abs(ny(l, j, facet, el)) - 1d0)
112 sz = sz + abs(abs(nz(l, j, facet, el)) - 1d0)
116 sx = sx / (c%Xh%lx - 2)**2
117 sy = sy / (c%Xh%lx - 2)**2
118 sz = sz / (c%Xh%lx - 2)**2
120 if (sx .lt. tol)
then
121 call this%bc_x%mark_facet(facet, el)
124 if (sy .lt. tol)
then
125 call this%bc_y%mark_facet(facet, el)
128 if (sz .lt. tol)
then
129 call this%bc_z%mark_facet(facet, el)
133 call this%bc_x%finalize()
134 call this%bc_x%set_g(0.0_rp)
135 call this%bc_y%finalize()
136 call this%bc_y%set_g(0.0_rp)
137 call this%bc_z%finalize()
138 call this%bc_z%set_g(0.0_rp)
145 integer,
intent(in) :: n
146 real(kind=rp),
intent(inout),
dimension(n) :: x
147 real(kind=rp),
intent(in),
optional :: t
148 integer,
intent(in),
optional :: tstep
154 integer,
intent(in) :: n
155 real(kind=rp),
intent(inout),
dimension(n) :: x
156 real(kind=rp),
intent(inout),
dimension(n) :: y
157 real(kind=rp),
intent(inout),
dimension(n) :: z
158 real(kind=rp),
intent(in),
optional :: t
159 integer,
intent(in),
optional :: tstep
161 call this%bc_x%apply_scalar(x, n)
162 call this%bc_y%apply_scalar(y, n)
163 call this%bc_z%apply_scalar(z, n)
169 class(
symmetry_t),
intent(inout),
target :: this
171 real(kind=rp),
intent(in),
optional :: t
172 integer,
intent(in),
optional :: tstep
177 class(
symmetry_t),
intent(inout),
target :: this
181 real(kind=rp),
intent(in),
optional :: t
182 integer,
intent(in),
optional :: tstep
184 call device_symmetry_apply_vector(this%bc_x%msk_d, this%bc_y%msk_d, &
185 this%bc_z%msk_d, x_d, y_d, z_d, &
194 class(
symmetry_t),
target,
intent(inout) :: this
196 call this%free_base()
197 call this%bc_x%free()
198 call this%bc_y%free()
199 call this%bc_z%free()
Defines a boundary condition.
subroutine device_symmetry_apply_vector(xmsk, ymsk, zmsk, x, y, z, m, n, l)
Defines a dirichlet boundary condition.
integer, parameter, public rp
Global precision used in computations.
Mixed Dirichlet-Neumann axis aligned symmetry plane.
subroutine symmetry_init(this, coef)
Initialize symmetry mask for each axis.
subroutine symmetry_apply_vector(this, x, y, z, n, t, tstep)
Apply symmetry conditions (axis aligned)
subroutine symmetry_apply_scalar(this, x, n, t, tstep)
No-op scalar apply.
subroutine symmetry_apply_scalar_dev(this, x_d, t, tstep)
No-op scalar apply (device version)
subroutine symmetry_free(this)
Destructor.
subroutine symmetry_apply_vector_dev(this, x_d, y_d, z_d, t, tstep)
Apply symmetry conditions (axis aligned) (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 .
Mixed Dirichlet-Neumann symmetry plane condition.