43 use json_module,
only : json_file
51 use,
intrinsic :: iso_c_binding, only : c_ptr, c_null_ptr, c_associated
59 character(len=:),
allocatable :: filter_type
65 real(kind=
rp),
allocatable :: fh(:,:), fht(:,:)
66 type(c_ptr) :: fh_d = c_null_ptr
67 type(c_ptr) :: fht_d = c_null_ptr
69 real(kind=
rp),
allocatable :: transfer(:)
74 procedure, pass(this) :: init_from_components => &
88 type(json_file),
intent(inout) :: json
89 type(
coef_t),
intent(in) :: coef
90 real(kind=
rp),
allocatable :: transfer(:)
91 character(len=:),
allocatable :: filter_type
94 call this%init_base(json, coef)
96 call this%init_from_components(coef%dof%xh%lx)
99 this%filter_type,
"nonBoyd")
101 if (json%valid_path(
'filter.transfer_function'))
then
102 call json_get(json,
'filter.transfer_function', transfer)
103 if (
size(transfer) .eq. coef%dof%xh%lx)
then
104 this%transfer = transfer
106 call neko_error(
"The transfer function of the elementwise " // &
107 "filter must correspond the order of the polynomial")
123 allocate(this%fh(nx, nx))
124 allocate(this%fht(nx, nx))
125 allocate(this%transfer(nx))
127 call rzero(this%fh, nx*nx)
128 call rzero(this%fht, nx*nx)
129 call rone(this%transfer, nx)
132 call device_map(this%fh, this%fh_d, this%nx * this%nx)
133 call device_map(this%fht, this%fht_d, this%nx * this%nx)
135 call device_cfill(this%fht_d, 0.0_rp, this%nx * this%nx)
144 if (
allocated(this%fh))
then
148 if (
allocated(this%fht))
then
152 if (
allocated(this%transfer))
then
153 deallocate(this%transfer)
156 if (c_associated(this%fh_d))
then
160 if (c_associated(this%fht_d))
then
164 this%filter_type =
""
168 call this%free_base()
177 this%nx, this%filter_type)
190 type(
field_t),
intent(inout) :: F_out
191 type(
field_t),
intent(in) :: F_in
194 call tnsr3d(f_out%x, this%nx, f_in%x, this%nx, this%fh, this%fht, this%fht, &
207 integer,
intent(in) :: nx
208 real(kind=
rp),
intent(inout) :: fh(nx, nx), fht(nx, nx)
209 real(kind=
rp),
intent(in) :: transfer(nx)
210 real(kind=
rp) :: diag(nx, nx), rmult(nx), lj(nx), zpts(nx)
212 integer :: n, i, j, k
214 character(len=*),
intent(in) :: filter_type
216 call phi%init(nx, nx)
217 call pht%init(nx, nx)
219 call zwgll(zpts, rmult, nx)
225 select case (filter_type)
230 pht%x(k,j) = lj(k)-lj(k-2)
237 call trsp(phi%x, nx, pht%x, nx)
245 diag(i,i) = transfer(i)
248 call mxm (diag, nx, pht%x, nx, fh, nx)
249 call mxm (phi%x, nx, fh, nx, pht%x, nx)
251 call copy (fh, pht%x, nx*nx)
252 call trsp (fht, nx, fh, nx)
Map a Fortran array to a device (allocate and associate)
Copy data between host and device (or device and device)
Retrieves a parameter by name or assigns a provided default value. In the latter case also adds the m...
Retrieves a parameter by name or throws an error.
subroutine, public device_cfill(a_d, c, n, strm)
Set all elements to a constant c .
Device abstraction, common interface for various accelerators.
integer, parameter, public host_to_device
subroutine, public device_free(x_d)
Deallocate memory on the device.
Implements elementwise_filter_t.
subroutine build_1d_cpu(fh, fht, transfer, nx, filter_type)
Build the 1d filter for an element on the CPU. Suppose field x is filtered into x_hat by x_hat = fh*x...
subroutine elementwise_filter_init_from_json(this, json, coef)
Constructor.
subroutine elementwise_filter_free(this)
Destructor.
subroutine build_1d(this)
Build the 1d filter for an element.
subroutine elementwise_filter_init_from_components(this, nx)
Actual Constructor.
subroutine elementwise_field_filter_3d(this, f_out, f_in)
Filter a 3D field.
Filter to be applied to a scalar field.
Utilities for retrieving parameters from the case files.
subroutine, public rone(a, n)
Set all elements to one.
subroutine, public copy(a, b, n)
Copy a vector .
subroutine, public rzero(a, n)
Zero a real vector.
Wrapper for all matrix-matrix product implementations.
subroutine, public mxm(a, n1, b, n2, c, n3)
Compute matrix-matrix product for contiguously packed matrices A,B, and C.
integer, parameter neko_bcknd_device
integer, parameter, public rp
Global precision used in computations.
LIBRARY ROUTINES FOR SPECTRAL METHODS.
subroutine zwgll(z, w, np)
Generate NP Gauss-Lobatto Legendre points (Z) and weights (W) associated with Jacobi polynomial P(N)(...
subroutine legendre_poly(l, x, n)
Evaluate Legendre polynomials of degrees 0-N at point x and store in array L.
subroutine, public trsp(a, lda, b, ldb)
Transpose of a rectangular tensor .
subroutine, public tnsr3d(v, nv, u, nu, a, bt, ct, nelv)
Tensor product performed on nelv elements.
Coefficients defined on a given (mesh, ) tuple. Arrays use indices (i,j,k,e): element e,...
Implements the elementwise filter for SEM.
Base abstract class for filter.