65 use,
intrinsic :: iso_c_binding
70 real(kind=
rp),
allocatable :: r(:)
71 real(kind=
rp),
allocatable :: rc(:)
72 real(kind=
rp),
allocatable :: tmp(:)
73 type(c_ptr) :: r_d = c_null_ptr
74 type(c_ptr) :: rc_d = c_null_ptr
75 type(c_ptr) :: tmp_d = c_null_ptr
103 max_iter, cheby_degree)
105 class(
ax_t),
target,
intent(in) :: ax
106 type(
space_t),
target,
intent(in) :: Xh
107 type(
coef_t),
target,
intent(in) :: coef
108 type(
mesh_t),
target,
intent(in) :: msh
109 type(
gs_t),
target,
intent(in) :: gs_h
110 type(
bc_list_t),
target,
intent(in) :: blst
111 integer,
intent(in) :: nlvls
112 integer,
intent(in) :: max_iter
113 integer,
intent(in) :: cheby_degree
114 integer :: lvl, n, mlvl, target_num_aggs
115 integer,
allocatable :: agg_nhbr(:,:), asdf(:,:)
116 character(len=LOG_SIZE) :: log_buf
120 write(log_buf,
'(A28,I2,A8)')
'Creating AMG hierarchy with', &
125 call this%amg%init(ax, xh, coef, msh, gs_h, nlvls, blst)
131 allocate( agg_nhbr, source = msh%facet_neigh )
133 target_num_aggs = this%amg%lvl(mlvl-1)%nnodes / 8
135 if ( target_num_aggs .lt. 4 )
then
137 "TAMG: Too many levels. Not enough DOFs for coarsest grid.")
143 deallocate( agg_nhbr )
148 this%max_iter = max_iter
150 this%nlvls = this%amg%nlvls
151 if (this%nlvls .gt. this%amg%nlvls)
then
153 "Requested number multigrid levels &
154 & is greater than the initialized AMG levels")
157 allocate(this%smoo(0:(nlvls)))
159 n = this%amg%lvl(lvl+1)%fine_lvl_dofs
160 call this%smoo(lvl)%init(n, lvl, cheby_degree)
164 allocate(this%wrk(nlvls))
166 n = this%amg%lvl(lvl)%fine_lvl_dofs
167 allocate( this%wrk(lvl)%r(n) )
168 allocate( this%wrk(lvl)%rc(n) )
169 allocate( this%wrk(lvl)%tmp(n) )
171 call device_map( this%wrk(lvl)%r, this%wrk(lvl)%r_d, n)
172 call device_map( this%wrk(lvl)%rc, this%wrk(lvl)%rc_d, n)
173 call device_map( this%wrk(lvl)%tmp, this%wrk(lvl)%tmp_d, n)
195 integer,
intent(in) :: n
197 real(kind=
rp),
dimension(n),
intent(inout) :: z
198 real(kind=
rp),
dimension(n),
intent(inout) :: r
201 integer :: iter, max_iter
202 logical :: zero_initial_guess
204 max_iter = this%max_iter
211 zero_initial_guess = .true.
213 do iter = 1, max_iter
216 zero_initial_guess = .false.
221 zero_initial_guess = .true.
223 do iter = 1, max_iter
226 zero_initial_guess = .false.
241 integer,
intent(in) :: n
242 real(kind=
rp),
intent(inout) :: x(n)
243 real(kind=
rp),
intent(inout) :: b(n)
246 logical,
intent(in) :: zero_initial_guess
247 integer,
intent(in) :: lvl
248 real(kind=
rp) :: r(n)
249 real(kind=
rp) :: rc(n)
250 real(kind=
rp) :: tmp(n)
251 integer :: iter, num_iter
254 max_lvl = mgstuff%nlvls-1
258 call mgstuff%smoo(lvl)%solve(x, b, n, amg, &
260 if (lvl .eq. max_lvl)
then
270 call amg%interp_f2c(rc, r, lvl+1)
275 call tamg_mg_cycle(tmp, rc, amg%lvl(lvl+1)%nnodes, lvl+1, amg, mgstuff, &
280 call amg%interp_c2f(r, tmp, lvl+1)
288 call mgstuff%smoo(lvl)%solve(x,b, n, amg)
300 integer,
intent(in) :: n
301 real(kind=
rp),
intent(inout) :: x(n)
302 real(kind=
rp),
intent(inout) :: b(n)
307 logical,
intent(in) :: zero_initial_guess
308 integer,
intent(in) :: lvl
309 integer :: iter, num_iter
312 max_lvl = mgstuff%nlvls-1
316 call mgstuff%smoo(lvl)%device_solve(x, b, x_d, b_d, n, amg, &
318 if (lvl .eq. max_lvl)
then
321 associate( r => mgstuff%wrk(lvl+1)%r, r_d => mgstuff%wrk(lvl+1)%r_d, &
322 rc => mgstuff%wrk(lvl+1)%rc, rc_d => mgstuff%wrk(lvl+1)%rc_d, &
323 tmp => mgstuff%wrk(lvl+1)%tmp, tmp_d => mgstuff%wrk(lvl+1)%tmp_d )
327 call amg%device_matvec(r, x, r_d, x_d, lvl)
332 call amg%interp_f2c_d(rc_d, r_d, lvl+1)
338 amg%lvl(lvl+1)%nnodes, lvl+1, amg, mgstuff, .true.)
342 call amg%interp_c2f_d(r_d, tmp_d, lvl+1, r)
350 call mgstuff%smoo(lvl)%device_solve(x, b, x_d, b_d, n, amg)
363 integer,
intent(in) :: n
364 real(kind=rp),
intent(inout) :: r(n)
365 real(kind=rp),
intent(inout) :: x(n)
366 real(kind=rp),
intent(inout) :: b(n)
367 type(tamg_hierarchy_t),
intent(inout) :: amg
368 integer,
intent(in) :: lvl
370 call amg%matvec(r, x, lvl)
371 call sub3(r, b, r, n)
376 integer,
intent(in) :: lvl, nagg
377 character(len=LOG_SIZE) :: log_buf
379 write(log_buf,
'(A8,I2,A31)')
'-- level', lvl, &
380 '-- Calling Greedy Aggregation'
381 call neko_log%message(log_buf)
382 write(log_buf,
'(A33,I6)')
'Target Aggregates:', nagg
383 call neko_log%message(log_buf)
387 integer,
intent(in) :: lvl, n
388 real(kind=rp),
intent(inout) :: r(n)
389 real(kind=rp),
intent(inout) :: x(n)
390 real(kind=rp),
intent(inout) :: b(n)
394 type(tamg_hierarchy_t),
intent(inout) :: amg
396 character(len=LOG_SIZE) :: log_buf
398 call amg%device_matvec(r, x, r_d, x_d, lvl)
399 call device_sub3(r_d, b_d, r_d, n)
400 val = device_glsc2(r_d, r_d, n)
402 write(log_buf,
'(A33,I6,F12.6)')
'tAMG resid:', lvl, val
403 call neko_log%message(log_buf)
407 type(tamg_hierarchy_t),
intent(inout) :: amg
408 integer :: i, j, k, l, nid, n
409 do j = 1, amg%lvl(1)%nnodes
410 do k = 1, amg%lvl(1)%nodes(j)%ndofs
411 nid = amg%lvl(1)%nodes(j)%dofs(k)
412 amg%lvl(1)%map_finest2lvl(nid) = amg%lvl(1)%nodes(j)%gid
415 n =
size(amg%lvl(1)%map_finest2lvl)
418 nid = amg%lvl(l-1)%map_finest2lvl(i)
419 do j = 1, amg%lvl(l)%nnodes
420 do k = 1, amg%lvl(l)%nodes(j)%ndofs
421 if (nid .eq. amg%lvl(l)%nodes(j)%dofs(k))
then
422 amg%lvl(l)%map_finest2lvl(i) = amg%lvl(l)%nodes(j)%gid
428 if (neko_bcknd_device .eq. 1)
then
430 amg%lvl(l)%map_finest2lvl(0) = n
431 call device_memcpy( amg%lvl(l)%map_finest2lvl, &
432 amg%lvl(l)%map_finest2lvl_d, n, &
433 host_to_device, .true.)
434 call device_memcpy( amg%lvl(l)%map_f2c, &
435 amg%lvl(l)%map_f2c_d, amg%lvl(l)%fine_lvl_dofs+1, &
436 host_to_device, .true.)
__device__ T solve(const T u, const T y, const T guess, const T nu, const T kappa, const T B)
Return the device pointer for an associated Fortran array.
Map a Fortran array to a device (allocate and associate)
Copy data between host and device (or device and device)
Defines a Matrix-vector product.
subroutine, public device_sub3(a_d, b_d, c_d, n, strm)
Vector subtraction .
subroutine, public device_add2(a_d, b_d, n, strm)
Vector addition .
subroutine, public device_rzero(a_d, n, strm)
Zero a real vector.
subroutine, public device_col2(a_d, b_d, n, strm)
Vector multiplication .
real(kind=rp) function, public device_glsc2(a_d, b_d, n, strm)
Weighted inner product .
Device abstraction, common interface for various accelerators.
integer, parameter, public host_to_device
subroutine, public device_free(x_d)
Deallocate memory on the device.
type(log_t), public neko_log
Global log stream.
integer, parameter, public log_size
real(kind=rp) function, public glsc2(a, b, n)
Weighted inner product .
subroutine, public sub3(a, b, c, n)
Vector subtraction .
subroutine, public add2(a, b, n)
Vector addition .
subroutine, public col2(a, b, n)
Vector multiplication .
subroutine, public rzero(a, n)
Zero a real vector.
integer, parameter neko_bcknd_device
integer, parameter, public rp
Global precision used in computations.
Defines a function space.
Implements an aggregation for TreeAMG hierarchy structure.
subroutine aggregate_end(tamg, lvl_id)
Aggregate all dofs to a single point to form a tree-like structure.
subroutine aggregate_greedy(tamg, lvl_id, max_aggs, facet_neigh, agg_nhbr)
Aggregates dofs based on adjacent dofs.
subroutine aggregate_finest_level(tamg, lx, ly, lz, ne)
Aggregaiton on finest level Aggregates all dofs in an element into a single aggregate.
Implements multigrid using the TreeAMG hierarchy structure. USE:
recursive subroutine tamg_mg_cycle(x, b, n, lvl, amg, mgstuff, zero_initial_guess)
Recrsive multigrid cycle for the TreeAMG solver object.
subroutine calc_resid(r, x, b, amg, lvl, n)
Wrapper function to calculate residyal.
subroutine print_preagg_info(lvl, nagg)
recursive subroutine tamg_mg_cycle_d(x, b, x_d, b_d, n, lvl, amg, mgstuff, zero_initial_guess)
Recrsive multigrid cycle for the TreeAMG solver object on device.
subroutine fill_lvl_map(amg)
subroutine print_resid_info(r, x, b, r_d, x_d, b_d, amg, lvl, n)
subroutine tamg_mg_solve(this, z, r, n)
Solver function for the TreeAMG solver object.
subroutine tamg_mg_init(this, ax, xh, coef, msh, gs_h, nlvls, blst, max_iter, cheby_degree)
Initialization of the TreeAMG multigrid solver.
Implements smoothers for use with TreeAMG matrix vector product.
Implements the base type for TreeAMG hierarchy structure.
subroutine, public tamg_node_init(node, gid, ndofs)
Initialization of a TreeAMG tree node.
subroutine, public tamg_lvl_init(tamg_lvl, lvl, nnodes, ndofs)
Initialization of a TreeAMG level.
Base type for a matrix-vector product providing .
A list of allocatable `bc_t`. Follows the standard interface of lists.
Coefficients defined on a given (mesh, ) tuple. Arrays use indices (i,j,k,e): element e,...
The function space for the SEM solution fields.
Type for a TreeAMG hierarchy.
Type for the TreeAMG solver.
Type for Chebyshev iteration using TreeAMG matvec.