Added C inteface to print ML preconditioner info

pizdaint-runs
Cirdans-Home 5 years ago
parent 8e2af97a35
commit a704873923

@ -29,7 +29,7 @@ extern "C" {
psb_i_t mld_c_dprecfree(mld_c_dprec *ph);
psb_i_t mld_c_dprecbld_opt(psb_c_dspmat *ah, psb_c_descriptor *cdh,
mld_c_dprec *ph, const char *afmt);
psb_i_t mld_c_ddescr(mld_c_dprec *ph);
psb_i_t mld_c_dkrylov(const char *method, psb_c_dspmat *ah, mld_c_dprec *ph,
psb_c_dvector *bh, psb_c_dvector *xh,

@ -30,6 +30,7 @@ extern "C" {
psb_i_t mld_c_zprecbld_opt(psb_c_zspmat *ah, psb_c_descriptor *cdh,
mld_c_zprec *ph, const char *afmt);
psb_i_t mld_c_zdescr(mld_c_zprec *ph);
psb_i_t mld_c_zkrylov(const char *method, psb_c_zspmat *ah, mld_c_zprec *ph,
psb_c_zvector *bh, psb_c_zvector *xh,

@ -377,4 +377,32 @@ contains
return
end function mld_c_dprecfree
function mld_c_ddescr(ph) bind(c) result(res)
use psb_base_mod
use mld_prec_mod
implicit none
integer(psb_c_ipk_) :: res
type(psb_c_object_type) :: ph
integer :: info
type(mld_dprec_type), pointer :: precp
res = -1
info = -1
if (c_associated(ph%item)) then
call c_f_pointer(ph%item,precp)
else
return
end if
call precp%descr()
call flush(output_unit)
info = 0
res = MLDC_ERR_FILTER(info)
MLDC_ERR_HANDLE(res)
return
end function mld_c_ddescr
end module mld_dprec_cbind_mod

@ -377,4 +377,32 @@ contains
return
end function mld_c_zprecfree
function mld_c_zdescr(ph) bind(c) result(res)
use psb_base_mod
use mld_prec_mod
implicit none
integer(psb_c_ipk_) :: res
type(psb_c_object_type) :: ph
integer :: info
type(mld_zprec_type), pointer :: precp
res = -1
info = -1
if (c_associated(ph%item)) then
call c_f_pointer(ph%item,precp)
else
return
end if
call precp%descr()
call flush(output_unit)
info = 0
res = MLDC_ERR_FILTER(info)
MLDC_ERR_HANDLE(res)
return
end function mld_c_zdescr
end module mld_zprec_cbind_mod

Loading…
Cancel
Save