Fixes for scaling in Krylov

kinsol-stop
sfilippone 6 days ago
parent c1fa595f0c
commit 6c0a8831d8

@ -34,8 +34,9 @@ extern "C" {
psb_i_t amg_c_ddescr(amg_c_dprec *ph);
psb_i_t amg_c_dkrylov(const char *method, psb_c_dspmat *ah, amg_c_dprec *ph,
psb_c_dvector *bh, psb_c_dvector *xh,
psb_c_descriptor *cdh, psb_c_SolverOptions *opt);
psb_c_dvector *bh, psb_c_dvector *xh,
psb_c_descriptor *cdh, psb_c_dvector *s1,
psb_c_dvector *s2, psb_c_SolverOptions *opt);
#ifdef __cplusplus

@ -35,8 +35,9 @@ extern "C" {
psb_i_t amg_c_zdescr(amg_c_zprec *ph);
psb_i_t amg_c_zkrylov(const char *method, psb_c_zspmat *ah, amg_c_zprec *ph,
psb_c_zvector *bh, psb_c_zvector *xh,
psb_c_descriptor *cdh, psb_c_SolverOptions *opt);
psb_c_zvector *bh, psb_c_zvector *xh,
psb_c_descriptor *cdh, psb_c_zvector *s1,
psb_c_zvector *s2, psb_c_SolverOptions *opt);
#ifdef __cplusplus

@ -3,6 +3,7 @@ module amg_dprec_cbind_mod
use iso_c_binding
use amg_prec_mod
use psb_base_cbind_mod
use psb_dlinsolve_cbind_mod
type, bind(c) :: amg_c_dprec
type(c_ptr) :: item = c_null_ptr
@ -18,7 +19,7 @@ contains
#define AMGC_ERROR(MSG)
#endif
#define amg_success_ 0
!#define AMGC_ERR_FILTER(INFO) min(0,INFO)
!#define AMGC_ERR_FILTER(INFO) min(0,INFO)
#define AMGC_ERR_FILTER(INFO) (INFO)
#define AMGC_ERR_HANDLE(INFO) if(INFO/=amg_success_)AMGC_ERROR("ERROR!")
function amg_c_dprecinit(cctxt,ph,ptype) bind(c) result(res)
@ -172,15 +173,17 @@ contains
end function amg_c_dprecbld
function amg_c_dhierarchy_build(ah,cdh,ph) bind(c) result(res)
use psb_base_mod
implicit none
integer(psb_c_ipk_) :: res
type(psb_c_object_type) :: ph,ah,cdh
integer(psb_ipk_) :: iret
type(amg_dprec_type), pointer :: precp
type(psb_dspmat_type), pointer :: ap
type(psb_desc_type), pointer :: descp
character(len=80) :: fptype
integer(psb_ipk_) :: iret, act
res = -1
@ -204,11 +207,16 @@ contains
res = AMGC_ERR_FILTER(iret)
AMGC_ERR_HANDLE(res)
if (res /=0) then
act = psb_act_abort_
call psb_error_handler(act)
end if
return
end function amg_c_dhierarchy_build
function amg_c_dsmoothers_build(ah,cdh,ph) bind(c) result(res)
use psb_base_mod
implicit none
integer(psb_c_ipk_) :: res
@ -217,7 +225,7 @@ contains
type(psb_dspmat_type), pointer :: ap
type(psb_desc_type), pointer :: descp
character(len=80) :: fptype
integer(psb_ipk_) :: iret
integer(psb_ipk_) :: iret, act
res = -1
@ -241,12 +249,16 @@ contains
res = AMGC_ERR_FILTER(iret)
AMGC_ERR_HANDLE(res)
if (res /=0) then
act = psb_act_abort_
call psb_error_handler(act)
end if
return
end function amg_c_dsmoothers_build
function amg_c_dkrylov(methd,&
& ah,ph,bh,xh,cdh,options) bind(c) result(res)
& ah,ph,bh,xh,cdh,s1,s2,options) bind(c) result(res)
use psb_base_mod
use psb_prec_mod
use psb_linsolve_mod
@ -254,20 +266,22 @@ contains
use psb_dlinsolve_cbind_mod
implicit none
integer(psb_c_ipk_) :: res
type(psb_c_object_type) :: ah,cdh,ph,bh,xh
type(psb_c_object_type) :: ah,cdh,ph,bh,xh,s1,s2
character(c_char) :: methd(*)
type(solveroptions) :: options
res= amg_c_dkrylov_opt(methd, ah, ph, bh, xh, options%eps,cdh, &
& itmax=options%itmax, iter=options%iter,&
& itrace=options%itrace, istop=options%istop,&
& irst=options%irst, err=options%err)
& irst=options%irst, err=options%err, s1=s1,s2=s2)
end function amg_c_dkrylov
function amg_c_dkrylov_opt(methd,&
& ah,ph,bh,xh,eps,cdh,itmax,iter,err,itrace,irst,istop) bind(c) result(res)
& ah,ph,bh,xh,eps,cdh,itmax,iter,err,&
& itrace,irst,istop,s1,s2) bind(c) result(res)
use psb_base_mod
use psb_prec_mod
use psb_linsolve_mod
@ -276,16 +290,17 @@ contains
use psb_base_string_cbind_mod
implicit none
integer(psb_c_ipk_) :: res
type(psb_c_object_type) :: ah,cdh,ph,bh,xh
type(psb_c_object_type) :: ah,cdh,ph,bh,xh,s1, s2
integer(psb_c_ipk_), value :: itmax,itrace,irst,istop
real(c_double), value :: eps
integer(psb_c_ipk_) :: iter
real(c_double) :: err
character(c_char) :: methd(*)
type(psb_desc_type), pointer :: descp
type(psb_dspmat_type), pointer :: ap
type(amg_dprec_type), pointer :: precp
type(psb_d_vect_type), pointer :: xp, bp
type(psb_d_vect_type), pointer :: xp, bp, s1p, s2p
integer(psb_ipk_) :: iret,fitmax,fitrace,first,fistop,fiter
character(len=20) :: fmethd
@ -317,7 +332,16 @@ contains
else
return
end if
if (c_associated(s1%item)) then
call c_f_pointer(s1%item,s1p)
else
nullify(s1p)
end if
if (c_associated(s2%item)) then
call c_f_pointer(s2%item,s2p)
else
nullify(s2p)
end if
call stringc2f(methd,fmethd)
feps = eps
@ -325,11 +349,40 @@ contains
fitrace = itrace
first = irst
fistop = istop
!!$ flush(0)
!!$
!!$ if (associated(s1p)) write(0,*) 'amg_c_dkrylov_opt S1 ',s1p%get_nrows(), &
!!$ & s1p%nrm2(s1p%get_nrows())
!!$ if (associated(s2p)) write(0,*) 'amg_c_dkrylov_opt S2 ',s2p%get_nrows(), &
!!$ & s2p%nrm2(s2p%get_nrows())
!!$ flush(0)
if (associated(s1p).and.associated(s2p)) then
call psb_krylov(fmethd, ap, precp, bp, xp, feps, &
& descp, iret,&
& itmax=fitmax,iter=fiter,itrace=fitrace,istop=fistop,&
& irst=first, err=ferr,s1=s1p,s2=s2p)
else if (associated(s1p)) then
call psb_krylov(fmethd, ap, precp, bp, xp, feps, &
& descp, iret,&
& itmax=fitmax,iter=fiter,itrace=fitrace,istop=fistop,&
& irst=first, err=ferr,s1=s1p)
else if (associated(s2p)) then
call psb_krylov(fmethd, ap, precp, bp, xp, feps, &
& descp, iret,&
& itmax=fitmax,iter=fiter,itrace=fitrace,istop=fistop,&
& irst=first, err=ferr,s2=s2p)
else
call psb_krylov(fmethd, ap, precp, bp, xp, feps, &
& descp, iret,&
& itmax=fitmax,iter=fiter,itrace=fitrace,istop=fistop,&
& irst=first, err=ferr)
end if
!!$ call psb_krylov(fmethd, ap, precp, bp, xp, feps, &
!!$ & descp, iret,&
!!$ & itmax=fitmax,iter=fiter,itrace=fitrace,istop=fistop,&
!!$ & irst=first, err=ferr)
call psb_krylov(fmethd, ap, precp, bp, xp, feps, &
& descp, iret,&
& itmax=fitmax,iter=fiter,itrace=fitrace,istop=fistop,&
& irst=first, err=ferr)
iter = fiter
err = ferr
res = min(iret,0)
@ -384,7 +437,7 @@ contains
res = AMGC_ERR_FILTER(info)
AMGC_ERR_HANDLE(res)
return
end function amg_c_dprecapply
end function amg_c_dprecapply
function amg_c_dprecapply_opt(ph,bc,xc,cdh,ctrans) bind(c,name="amg_c_dprecapply_opt") result(res)
use psb_base_mod
@ -432,7 +485,7 @@ end function amg_c_dprecapply
! Convert transpose flag
call stringc2f(ctrans,ftrans)
! Apply preconditioner
call precp%apply(bp,xp,descp,info,trans=ftrans)
@ -440,7 +493,7 @@ end function amg_c_dprecapply
res = AMGC_ERR_FILTER(info)
AMGC_ERR_HANDLE(res)
return
end function amg_c_dprecapply_opt
end function amg_c_dprecapply_opt
function amg_c_dprecfree(ph) bind(c) result(res)
implicit none

@ -3,6 +3,7 @@ module amg_zprec_cbind_mod
use iso_c_binding
use amg_prec_mod
use psb_base_cbind_mod
use psb_zlinsolve_cbind_mod
type, bind(c) :: amg_c_zprec
type(c_ptr) :: item = c_null_ptr
@ -140,11 +141,11 @@ contains
integer(psb_c_ipk_) :: res
type(psb_c_object_type) :: ph,ah,cdh
integer(psb_ipk_) :: iret
type(amg_zprec_type), pointer :: precp
type(psb_zspmat_type), pointer :: ap
type(psb_desc_type), pointer :: descp
character(len=80) :: fptype
integer(psb_ipk_) :: iret
res = -1
@ -173,15 +174,17 @@ contains
end function amg_c_zprecbld
function amg_c_zhierarchy_build(ah,cdh,ph) bind(c) result(res)
use psb_base_mod
implicit none
integer(psb_c_ipk_) :: res
type(psb_c_object_type) :: ph,ah,cdh
integer(psb_ipk_) :: iret
type(amg_zprec_type), pointer :: precp
type(psb_zspmat_type), pointer :: ap
type(psb_desc_type), pointer :: descp
character(len=80) :: fptype
integer(psb_ipk_) :: iret, act
res = -1
@ -205,20 +208,25 @@ contains
res = AMGC_ERR_FILTER(iret)
AMGC_ERR_HANDLE(res)
if (res /=0) then
act = psb_act_abort_
call psb_error_handler(act)
end if
return
end function amg_c_zhierarchy_build
function amg_c_zsmoothers_build(ah,cdh,ph) bind(c) result(res)
use psb_base_mod
implicit none
integer(psb_c_ipk_) :: res
type(psb_c_object_type) :: ph,ah,cdh
integer(psb_ipk_) :: iret
type(amg_zprec_type), pointer :: precp
type(psb_zspmat_type), pointer :: ap
type(psb_desc_type), pointer :: descp
character(len=80) :: fptype
integer(psb_ipk_) :: iret, act
res = -1
@ -242,37 +250,47 @@ contains
res = AMGC_ERR_FILTER(iret)
AMGC_ERR_HANDLE(res)
if (res /=0) then
act = psb_act_abort_
call psb_error_handler(act)
end if
return
end function amg_c_zsmoothers_build
function amg_c_zkrylov(methd,&
& ah,ph,bh,xh,cdh,options) bind(c) result(res)
& ah,ph,bh,xh,cdh,s1,s2,options) bind(c) result(res)
use psb_base_mod
use psb_prec_mod
use psb_linsolve_mod
use psb_prec_cbind_mod
use psb_zlinsolve_cbind_mod
implicit none
integer(psb_c_ipk_) :: res
type(psb_c_object_type) :: ah,cdh,ph,bh,xh
type(psb_c_object_type) :: ah,cdh,ph,bh,xh,s1,s2
character(c_char) :: methd(*)
type(solveroptions) :: options
res= amg_c_zkrylov_opt(methd, ah, ph, bh, xh, options%eps,cdh, &
& itmax=options%itmax, iter=options%iter,&
& itrace=options%itrace, istop=options%istop,&
& irst=options%irst, err=options%err)
& irst=options%irst, err=options%err, s1=s1,s2=s2)
end function amg_c_zkrylov
function amg_c_zkrylov_opt(methd,&
& ah,ph,bh,xh,eps,cdh,itmax,iter,err,itrace,irst,istop) bind(c) result(res)
& ah,ph,bh,xh,eps,cdh,itmax,iter,err,&
& itrace,irst,istop,s1,s2) bind(c) result(res)
use psb_base_mod
use psb_prec_mod
use psb_linsolve_mod
use psb_objhandle_mod
use psb_prec_cbind_mod
use psb_base_string_cbind_mod
implicit none
integer(psb_c_ipk_) :: res
type(psb_c_object_type) :: ah,cdh,ph,bh,xh
type(psb_c_object_type) :: ah,cdh,ph,bh,xh,s1, s2
integer(psb_c_ipk_), value :: itmax,itrace,irst,istop
real(c_double), value :: eps
integer(psb_c_ipk_) :: iter
@ -281,7 +299,7 @@ contains
type(psb_desc_type), pointer :: descp
type(psb_zspmat_type), pointer :: ap
type(amg_zprec_type), pointer :: precp
type(psb_z_vect_type), pointer :: xp, bp
type(psb_z_vect_type), pointer :: xp, bp, s1p, s2p
integer(psb_ipk_) :: iret,fitmax,fitrace,first,fistop,fiter
character(len=20) :: fmethd
@ -313,7 +331,16 @@ contains
else
return
end if
if (c_associated(s1%item)) then
call c_f_pointer(s1%item,s1p)
else
nullify(s1p)
end if
if (c_associated(s2%item)) then
call c_f_pointer(s2%item,s2p)
else
nullify(s2p)
end if
call stringc2f(methd,fmethd)
feps = eps
@ -322,10 +349,28 @@ contains
first = irst
fistop = istop
call psb_krylov(fmethd, ap, precp, bp, xp, feps, &
& descp, iret,&
& itmax=fitmax,iter=fiter,itrace=fitrace,istop=fistop,&
& irst=first, err=ferr)
if (associated(s1p).and.associated(s2p)) then
call psb_krylov(fmethd, ap, precp, bp, xp, feps, &
& descp, iret,&
& itmax=fitmax,iter=fiter,itrace=fitrace,istop=fistop,&
& irst=first, err=ferr,s1=s1p,s2=s2p)
else if (associated(s1p)) then
call psb_krylov(fmethd, ap, precp, bp, xp, feps, &
& descp, iret,&
& itmax=fitmax,iter=fiter,itrace=fitrace,istop=fistop,&
& irst=first, err=ferr,s1=s1p)
else if (associated(s2p)) then
call psb_krylov(fmethd, ap, precp, bp, xp, feps, &
& descp, iret,&
& itmax=fitmax,iter=fiter,itrace=fitrace,istop=fistop,&
& irst=first, err=ferr,s2=s2p)
else
call psb_krylov(fmethd, ap, precp, bp, xp, feps, &
& descp, iret,&
& itmax=fitmax,iter=fiter,itrace=fitrace,istop=fistop,&
& irst=first, err=ferr)
end if
iter = fiter
err = ferr
res = min(iret,0)

Loading…
Cancel
Save