2 #include <petsc/finclude/petscksp.h>
13 public :: kspsetconvergencetest
23 integer(I4B) :: icnvg_ims
24 integer(I4B) :: icnvgopt
27 integer(I4B) :: max_its
28 real(dp) :: rnorm_l2_init
30 real(dp) :: t_convergence_check
39 subroutine cnvgcheckfunc(ksp, n, rnorm, flag, context, ierr)
44 kspconvergedreason :: flag
46 petscerrorcode :: ierr
49 subroutine cnvgdestroyfunc(context, ierr)
52 petscerrorcode :: ierr
55 subroutine kspsetconvergencetest(ksp, check_convergence, context, &
59 procedure(cnvgcheckfunc) :: check_convergence
61 procedure(cnvgdestroyfunc) :: destroy
62 petscerrorcode :: ierr
68 subroutine create(this, mat, settings, summary)
74 petscerrorcode :: ierr
77 this%icnvgopt = settings%icnvgopt
78 this%dvclose = settings%dvclose
79 this%rclose = settings%rclose
80 this%max_its = settings%iter1
81 this%cnvg_summary => summary
82 call matcreatevecs(mat, this%x_old, petsc_null_vec, ierr)
84 call matcreatevecs(mat, this%delta_x, petsc_null_vec, ierr)
86 call matcreatevecs(mat, this%residual, petsc_null_vec, ierr)
97 kspconvergedreason :: flag
99 petscerrorcode :: ierr
101 petscreal,
parameter :: min_one = -1.0
102 petscreal,
dimension(:),
pointer :: local_dx, local_res
103 petscreal :: xnorm_inf_ims, rnorm_inf_ims, rnorm_l2_ims
104 petscreal :: dvmax_model, rmax_model
105 petscint :: idx_dv, idx_r
110 petscint :: i, j, istart, iend
112 summary => context%cnvg_summary
116 call kspbuildsolution(ksp, petsc_null_vec, x, ierr)
119 call kspgetrhs(ksp, rhs, ierr)
122 call kspgetoperators(ksp, amat, petsc_null_mat, ierr)
125 call matmult(amat, x, context%residual, ierr)
129 call vecaypx(context%residual, -1.0_dp, rhs, ierr)
132 call vecnorm(context%residual, norm_2, rnorm_l2_ims, ierr)
137 context%rnorm_L2_init = rnorm_l2_ims
140 flag = ksp_converged_happy_breakdown
142 call veccopy(x, context%x_old, ierr)
144 flag = ksp_converged_iterating
151 summary%iter_cnt = summary%iter_cnt + 1
152 iter_cnt = summary%iter_cnt
154 if (summary%nitermax > 1)
then
155 summary%itinner(iter_cnt) = n
156 do i = 1, summary%convnmod
157 summary%convdvmax(i, iter_cnt) =
dzero
158 summary%convlocdv(i, iter_cnt) = 0
159 summary%convrmax(i, iter_cnt) =
dzero
160 summary%convlocr(i, iter_cnt) = 0
164 call vecwaxpy(context%delta_x, min_one, context%x_old, x, ierr)
167 call vecnorm(context%delta_x, norm_infinity, xnorm_inf_ims, ierr)
171 if (context%icnvgopt == 0 .or. context%icnvgopt == 1)
then
172 call vecnorm(context%residual, norm_infinity, rnorm_inf_ims, ierr)
176 call veccopy(x, context%x_old, ierr)
180 call vecgetarrayreadf90(context%delta_x, local_dx, ierr)
182 call vecgetarrayreadf90(context%residual, local_res, ierr)
184 do i = 1, summary%convnmod
191 istart = summary%model_bounds(i)
192 iend = summary%model_bounds(i + 1) - 1
194 if (abs(local_dx(j)) > abs(dvmax_model))
then
195 dvmax_model = local_dx(j)
198 if (abs(local_res(j)) > abs(rmax_model))
then
199 rmax_model = local_res(j)
203 if (summary%nitermax > 1)
then
204 summary%convdvmax(i, iter_cnt) = dvmax_model
205 summary%convlocdv(i, iter_cnt) = idx_dv
206 summary%convrmax(i, iter_cnt) = rmax_model
207 summary%convlocr(i, iter_cnt) = idx_r
210 call vecrestorearrayf90(context%delta_x, local_dx, ierr)
212 call vecrestorearrayf90(context%residual, local_res, ierr)
217 flag = ksp_converged_happy_breakdown
220 flag =
apply_check(context, n, xnorm_inf_ims, rnorm_inf_ims, rnorm_l2_ims)
223 if (flag == ksp_converged_iterating)
then
225 if (n == context%max_its)
then
226 flag = ksp_diverged_its
234 function apply_check(ctx, nit, dvmax, rnorm_inf, rnorm_L2)
result(flag)
240 real(dp) :: rnorm_inf
242 kspconvergedreason :: flag
248 flag = ksp_converged_iterating
253 if (ctx%icnvgopt == 2 .or. &
254 ctx%icnvgopt == 3 .or. &
255 ctx%icnvgopt == 4)
then
261 dvmax, rcnvg, ctx%rnorm_L2_init, &
262 epfact, ctx%dvclose, ctx%rclose)
264 if (ctx%icnvg_ims /= 0)
then
266 flag = ksp_converged_happy_breakdown
276 call vecdestroy(this%x_old, ierr)
278 call vecdestroy(this%delta_x, ierr)
280 call vecdestroy(this%residual, ierr)
This module contains simulation constants.
real(dp), parameter dzero
real constant zero
real(dp), parameter dprec
real constant machine precision
subroutine destroy(this)
Cleanup.
This module contains the IMS linear accelerator subroutines.
subroutine ims_base_testcnvg(Icnvgopt, Icnvg, Iiter, Dvmax, Rmax, Rmax0, Epfact, Dvclose, Rclose)
@ brief Test for solver convergence
real(dp) function ims_base_epfact(icnvgopt, kstp)
Function returning EPFACT.
This module defines variable data types.
real(dp), parameter, private rnorm_l2_tol
Context for the custom convergence check.
function apply_check(ctx, nit, dvmax, rnorm_inf, rnorm_L2)
Apply the IMS convergence check.
subroutine, public petsc_check_convergence(ksp, n, rnorm_L2, flag, context, ierr)
Routine to check the convergence. This is called.
subroutine create(this, mat, settings, summary)
integer(i4b), pointer, public kstp
current time step number
This structure stores the generic convergence info for a solution.
x vector from the previous iteration