15 use,
intrinsic :: iso_c_binding
23 integer,
allocatable :: dir_el(:)
25 integer,
allocatable :: el_lvl(:)
27 integer,
allocatable :: pt_lvl(:,:,:,:)
33 type(
mesh_t),
pointer :: msh => null()
37 real(kind=
rp) :: tol = 1e-7
50 type(
dofmap_t),
target,
intent(in) :: dof
51 type(
gs_t),
intent(inout) :: gs
52 integer,
intent(in) :: dir
53 real(kind=
rp),
intent(in) :: tol
54 integer :: nelv, lx, n, i, e, lvl
55 real(kind=
rp),
contiguous,
pointer :: line(:,:,:,:)
56 real(kind=
rp),
allocatable :: min_vals(:,:,:,:)
57 type(c_ptr) :: min_vals_d = c_null_ptr
58 real(kind=
rp) :: el_dim(3,3), glb_min, glb_max, el_min
61 if (neko_bcknd_device .eq. 1)
then
63 call neko_warning(
'map_1d does not copy indices to device, but ok if used on cpu')
76 else if (dir .eq. 2)
then
78 else if(dir .eq. 3)
then
81 call neko_error(
'Invalid dir for geopmetric comm')
84 allocate(this%dir_el(nelv))
85 allocate(this%el_lvl(nelv))
86 allocate(min_vals(lx, lx, lx, nelv))
87 allocate(this%pt_lvl(lx, lx, lx, nelv))
88 if (neko_bcknd_device .eq. 1)
then
97 el_dim(1,:) = abs(this%msh%elements(i)%e%pts(1)%p%x - this%msh%elements(i)%e%pts(2)%p%x)
98 el_dim(1,:) = el_dim(1,:)/norm2(el_dim(1,:))
99 el_dim(2,:) = abs(this%msh%elements(i)%e%pts(1)%p%x - this%msh%elements(i)%e%pts(3)%p%x)
100 el_dim(2,:) = el_dim(2,:)/norm2(el_dim(2,:))
101 el_dim(3,:) = abs(this%msh%elements(i)%e%pts(1)%p%x - this%msh%elements(i)%e%pts(5)%p%x)
102 el_dim(3,:) = el_dim(3,:)/norm2(el_dim(3,:))
105 this%dir_el(i) = maxloc(el_dim(:,this%dir),dim=1)
107 glb_min =
glmin(line,n)
108 glb_max =
glmax(line,n)
114 el_min = minval(line(:,:,:,e))
115 min_vals(:,:,:,e) = el_min
117 if (
relcmp(el_min, glb_min, this%tol))
then
118 if(this%el_lvl(e) .eq. -1) this%el_lvl(e) = i
124 do while (.not.
relcmp(
glmax(min_vals,n), glb_min, this%tol))
128 if (this%dir_el(e) .eq. 1)
then
129 if (line(1,1,1,e) .gt. line(lx,1,1,e))
then
130 min_vals(lx,:,:,e) = glb_max
132 min_vals(1,:,:,e) = glb_max
135 if (this%dir_el(e) .eq. 2)
then
136 if (line(1,1,1,e) .gt. line(1,lx,1,e))
then
137 min_vals(:,lx,:,e) = glb_max
139 min_vals(:,1,:,e) = glb_max
142 if (this%dir_el(e) .eq. 3)
then
143 if (line(1,1,1,e) .gt. line(1,1,lx,e))
then
144 min_vals(:,:,lx,e) = glb_max
146 min_vals(:,:,1,e) = glb_max
150 if (neko_bcknd_device .eq. 1) &
154 call gs%op(min_vals,n,gs_op_min)
155 if (neko_bcknd_device .eq. 1) &
163 el_min = minval(min_vals(:,:,:,e))
164 min_vals(:,:,:,e) = el_min
165 if (
relcmp(el_min, glb_min, this%tol))
then
166 if (this%el_lvl(e) .eq. -1) this%el_lvl(e) = i
170 this%n_el_lvls =
glimax(this%el_lvl,nelv)
172 write(*,*)
'Number of element levels', this%n_el_lvls
178 lvl = lx * (this%el_lvl(e) - 1) + i
179 if (this%dir_el(e) .eq. 1)
then
180 if (line(1,1,1,e) .gt. line(lx,1,1,e))
then
181 this%pt_lvl(lx-i+1,:,:,e) = lvl
183 this%pt_lvl(i,:,:,e) = lvl
186 if (this%dir_el(e) .eq. 2)
then
187 if (line(1,1,1,e) .gt. line(1,lx,1,e))
then
188 this%pt_lvl(:,lx-i+1,:,e) = lvl
190 this%pt_lvl(:,i,:,e) = lvl
193 if (this%dir_el(e) .eq. 3)
then
194 if (line(1,1,1,e) .gt. line(1,1,lx,e))
then
195 this%pt_lvl(:,:,lx-i+1,e) = lvl
197 this%pt_lvl(:,:,i,e) = lvl
210 if(
allocated(this%dir_el))
deallocate(this%dir_el)
211 if(
allocated(this%el_lvl))
deallocate(this%el_lvl)
212 if(
allocated(this%pt_lvl))
deallocate(this%pt_lvl)
213 if(
associated(this%dof))
nullify(this%dof)
214 if(
associated(this%msh))
nullify(this%msh)
Deassociate a Fortran array from a device pointer.
Map a Fortran array to a device (allocate and associate)
Copy data between host and device (or device and device)
Device abstraction, common interface for various accelerators.
integer, parameter, public host_to_device
subroutine, public device_free(x_d)
Deallocate memory on the device.
integer, parameter, public device_to_host
Defines a mapping of the degrees of freedom.
type(log_t), public neko_log
Global log stream.
integer, parameter, public log_size
Creates a 1d GLL point map along a specified direction based on the connectivity in the mesh.
subroutine map_1d_free(this)
subroutine map_1d_init(this, dof, gs, dir, tol)
integer function, public glimin(a, n)
Min of an integer vector of length n.
real(kind=rp) function, public glmax(a, n)
Max of a vector of length n.
integer function, public glimax(a, n)
Max of an integer vector of length n.
real(kind=rp) function, public glmin(a, n)
Min of a vector of length n.
integer, parameter, public rp
Global precision used in computations.
Defines a function space.
subroutine neko_warning(warning_msg)
Type that encapsulates a mapping from each gll point in the mesh to its corresponding (global) GLL po...
The function space for the SEM solution fields.