137 subroutine hsmg_init(this, msh, Xh, coef, dof, gs_h, bclst, crs_pctype)
138 class(
hsmg_t),
intent(inout),
target :: this
139 type(
mesh_t),
intent(inout),
target :: msh
140 type(
space_t),
intent(inout),
target :: Xh
141 type(
coef_t),
intent(inout),
target :: coef
142 type(
dofmap_t),
intent(inout),
target :: dof
143 type(
gs_t),
intent(inout),
target :: gs_h
144 type(
bc_list_t),
intent(inout),
target :: bclst
145 character(len=*),
optional :: crs_pctype
147 integer :: lx_crs, lx_mid
152 if (xh%lx .lt. 5)
then
153 lx_mid =
max(xh%lx-1,3)
155 if(xh%lx .le. 2)
then
156 call neko_error(
'Polynomial order < 2 not supported for hsmg precon')
163 allocate(this%grids(this%nlvls))
164 allocate(this%w(dof%size()))
165 allocate(this%r(dof%size()))
169 call msh%all_deformed()
172 call this%e%init(dof,
'work array')
173 call this%wf%init(dof,
'work 2')
175 call this%Xh_crs%init(
gll, lx_crs, lx_crs, lx_crs)
176 call this%dm_crs%init(msh, this%Xh_crs)
177 call this%gs_crs%init(this%dm_crs)
178 call this%e_crs%init(this%dm_crs,
'work crs')
179 call this%c_crs%init(this%gs_crs)
181 call this%Xh_mg%init(
gll, lx_mid, lx_mid, lx_mid)
182 call this%dm_mg%init(msh, this%Xh_mg)
183 call this%gs_mg%init(this%dm_mg)
184 call this%e_mg%init(this%dm_mg,
'work midl')
185 call this%c_mg%init(this%gs_mg)
188 call ax_helm_factory(this%ax, full_formulation = .false.)
191 call precon_factory(this%pc_crs,
'jacobi')
194 if (
present(crs_pctype))
then
195 call krylov_solver_factory(this%crs_solver, &
196 this%dm_crs%size(), trim(crs_pctype),
ksp_max_iter, m = this%pc_crs)
198 call krylov_solver_factory(this%crs_solver, &
199 this%dm_crs%size(),
'cg',
ksp_max_iter, m = this%pc_crs)
202 call this%bc_crs%init_base(this%c_crs)
203 call this%bc_mg%init_base(this%c_mg)
204 call this%bc_reg%init_base(coef)
205 if (bclst%n .gt. 0)
then
207 call this%bc_reg%mark_facets(bclst%bc(i)%bcp%marked_facet)
208 call this%bc_crs%mark_facets(bclst%bc(i)%bcp%marked_facet)
209 call this%bc_mg%mark_facets(bclst%bc(i)%bcp%marked_facet)
212 call this%bc_reg%finalize()
213 call this%bc_reg%set_g(
real(0d0,
rp))
217 call this%bc_crs%finalize()
218 call this%bc_crs%set_g(
real(0d0,
rp))
223 call this%bc_mg%finalize()
224 call this%bc_mg%set_g(0.0_rp)
228 call this%schwarz%init(xh, dof, gs_h, this%bclst_reg, msh)
229 call this%schwarz_mg%init(this%Xh_mg, this%dm_mg, this%gs_mg,&
232 call this%interp_fine_mid%init(xh, this%Xh_mg)
233 call this%interp_mid_crs%init(this%Xh_mg, this%Xh_crs)
235 call hsmg_fill_grid(dof, gs_h, xh, coef, this%bclst_reg, this%schwarz, &
236 this%e, this%grids, 3)
237 call hsmg_fill_grid(this%dm_mg, this%gs_mg, this%Xh_mg, this%c_mg, &
238 this%bclst_mg, this%schwarz_mg, this%e_mg, &
241 this%c_crs, this%bclst_crs, this%schwarz_crs, &
242 this%e_crs, this%grids, 1)
250 select type(pc => this%pc_crs)
252 call pc%init(this%c_crs, this%dm_crs, this%gs_crs)
254 call pc%init(this%c_crs, this%dm_crs, this%gs_crs)
256 call pc%init(this%c_crs, this%dm_crs, this%gs_crs)
345 integer,
intent(in) :: n
346 class(
hsmg_t),
intent(inout) :: this
347 real(kind=
rp),
dimension(n),
intent(inout) :: z
348 real(kind=
rp),
dimension(n),
intent(inout) :: r
349 type(c_ptr) :: z_d, r_d
351 integer :: thrdid, nthrds
363 call device_col2(this%r_d, this%grids(3)%coef%mult_d, &
364 this%grids(3)%dof%size())
366 call this%interp_fine_mid%map(this%e%x, this%r, &
367 this%msh%nelv, this%grids(2)%Xh)
368 call this%grids(2)%gs_h%op(this%e%x, &
369 this%grids(2)%dof%size(), gs_op_add, this%gs_event)
373 call device_copy(this%w_d, this%e%x_d, this%grids(2)%dof%size())
375 this%grids(2)%dof%size())
377 call device_col2(this%w_d, this%grids(2)%coef%mult_d, &
378 this%grids(2)%dof%size())
380 call this%interp_mid_crs%map(this%wf%x, this%w, this%msh%nelv, &
383 call device_copy(this%w_d, this%e%x_d, this%grids(2)%dof%size())
385 this%grids(2)%dof%size())
394 if (thrdid .eq. 0)
then
396 call this%grids(3)%schwarz%compute(z, this%r)
397 call this%grids(2)%schwarz%compute(this%grids(2)%e%x, this%w)
400 if (nthrds .eq. 1 .or. thrdid .eq. 1)
then
402 call this%grids(1)%gs_h%op(this%wf%x, &
403 this%grids(1)%dof%size(), gs_op_add, this%gs_event)
406 this%grids(1)%dof%size())
408 crs_info = this%crs_solver%solve(this%Ax, this%grids(1)%e, &
410 this%grids(1)%dof%size(), &
411 this%grids(1)%coef, &
412 this%grids(1)%bclst, &
413 this%grids(1)%gs_h, this%niter)
416 this%grids(1)%dof%size())
421 call this%interp_mid_crs%map(this%w, this%grids(1)%e%x, &
422 this%msh%nelv, this%grids(2)%Xh)
423 call device_add2(this%grids(2)%e%x_d, this%w_d, this%grids(2)%dof%size())
425 call this%interp_fine_mid%map(this%w, this%grids(2)%e%x, &
426 this%msh%nelv, this%grids(3)%Xh)
427 call device_add2(z_d, this%w_d, this%grids(3)%dof%size())
428 call this%grids(3)%gs_h%op(z, this%grids(3)%dof%size(), &
429 gs_op_add, this%gs_event)
432 this%grids(3)%dof%size())
435 call copy(this%r, r, n)
438 call this%grids(3)%schwarz%compute(z, this%r)
440 call col2(this%r, this%grids(3)%coef%mult, &
441 this%grids(3)%dof%size())
443 call this%interp_fine_mid%map(this%w, this%r, &
444 this%msh%nelv, this%grids(2)%Xh)
445 call this%grids(2)%gs_h%op(this%w, this%grids(2)%dof%size(), gs_op_add)
447 call this%grids(2)%schwarz%compute(this%grids(2)%e%x, this%w)
448 call col2(this%w, this%grids(2)%coef%mult, this%grids(2)%dof%size())
450 call this%interp_mid_crs%map(this%r, this%w, &
451 this%msh%nelv, this%grids(1)%Xh)
454 call this%grids(1)%gs_h%op(this%r, this%grids(1)%dof%size(), gs_op_add)
456 this%grids(1)%dof%size())
458 crs_info = this%crs_solver%solve(this%Ax, this%grids(1)%e, this%r, &
459 this%grids(1)%dof%size(), &
460 this%grids(1)%coef, &
461 this%grids(1)%bclst, &
462 this%grids(1)%gs_h, this%niter)
465 this%grids(1)%dof%size())
468 call this%interp_mid_crs%map(this%w, this%grids(1)%e%x, &
469 this%msh%nelv, this%grids(2)%Xh)
470 call add2(this%grids(2)%e%x, this%w, this%grids(2)%dof%size())
472 call this%interp_fine_mid%map(this%w, this%grids(2)%e%x, &
473 this%msh%nelv, this%grids(3)%Xh)
474 call add2(z, this%w, this%grids(3)%dof%size())
475 call this%grids(3)%gs_h%op(z, this%grids(3)%dof%size(), gs_op_add)
476 call col2(z, this%grids(3)%coef%mult, this%grids(3)%dof%size())