Start refactoring of C interface
Salvatore Filippone
7 years ago
parent
c8852cc22e
commit
ecfe0c1d42
@ -0,0 +1,27 @@
|
||||
include ../Make.inc
|
||||
|
||||
HERE=.
|
||||
LIBDIR=../lib
|
||||
INCDIR=../include
|
||||
MODDIR=../modules/
|
||||
LIBNAME=$(CBINDLIBNAME)
|
||||
LIBNAME=libmld_cbind.a
|
||||
|
||||
lib: mlprecd
|
||||
# /bin/cp -p $(CPUPDFLAG) $(HERE)/$(LIBNAME) $(LIBDIR)#
|
||||
# /bin/cp -p $(CPUPDFLAG) *.h $(INCDIR)
|
||||
# /bin/cp -p $(CPUPDFLAG) *$(.mod) $(MODDIR)
|
||||
|
||||
|
||||
mlprecd:
|
||||
cd mlprec && $(MAKE) lib LIBNAME=$(LIBNAME)
|
||||
|
||||
|
||||
clean:
|
||||
cd mlprec && $(MAKE) clean
|
||||
cd krylov && $(MAKE) clean
|
||||
|
||||
|
||||
veryclean: clean
|
||||
cd test/pargen && $(MAKE) clean
|
||||
/bin/rm -f $(HERE)/$(LIBNAME) $(LIBMOD) *$(.mod) *.h
|
||||
@ -0,0 +1,34 @@
|
||||
include ../../Make.inc
|
||||
LIBDIR=../../lib
|
||||
INCDIR=../../include
|
||||
MODDIR=../../modules
|
||||
HERE=.
|
||||
|
||||
CINCLUDES=-I. -I$(LIBDIR) -I$(PSBLAS_INCDIR)
|
||||
FINCLUDES=$(FMFLAG)$(HERE) $(FMFLAG)$(INCDIR) $(FMFLAG)$(MODDIR) $(PSBLAS_INCLUDES)
|
||||
|
||||
|
||||
OBJS=mld_prec_cbind_mod.o mld_dprec_cbind_mod.o mld_c_dprec.o
|
||||
CMOD=mld_cbind.h mld_c_dprec.h mld_const.h
|
||||
|
||||
|
||||
LIBMOD=mld_prec_cbind_mod$(.mod)
|
||||
LOCAL_MODS=$(LIBMOD)
|
||||
#LIBNAME=$(CPRECLIBNAME)
|
||||
|
||||
|
||||
lib: $(OBJS) $(CMOD)
|
||||
$(AR) $(HERE)/$(LIBNAME) $(OBJS)
|
||||
$(RANLIB) $(HERE)/$(LIBNAME)
|
||||
/bin/cp -p $(HERE)/$(LIBNAME) $(LIBDIR)
|
||||
/bin/cp -p $(LIBMOD) $(CMOD) $(INCDIR)
|
||||
|
||||
mld_prec_cbind_mod.o: mld_dprec_cbind_mod.o
|
||||
#mld_prec_cbind_mod.o: psb_prec_cbind_mod.o
|
||||
veryclean: clean
|
||||
/bin/rm -f $(HERE)/$(LIBNAME)
|
||||
|
||||
clean:
|
||||
/bin/rm -f $(OBJS) $(LOCAL_MODS)
|
||||
|
||||
veryclean: clean
|
||||
@ -0,0 +1,21 @@
|
||||
#include <stdlib.h>
|
||||
#include "mld_c_dprec.h"
|
||||
|
||||
mld_c_dprec* mld_c_new_dprec()
|
||||
{
|
||||
mld_c_dprec* temp;
|
||||
|
||||
temp=(mld_c_dprec *) malloc(sizeof(mld_c_dprec));
|
||||
temp->dprec=NULL;
|
||||
return(temp);
|
||||
}
|
||||
|
||||
|
||||
int mld_c_delete_dprec(mld_c_dprec* p)
|
||||
{
|
||||
int iret;
|
||||
iret=mld_c_dprecfree(p);
|
||||
if (iret ==0) free(p);
|
||||
return(iret);
|
||||
}
|
||||
|
||||
@ -0,0 +1,43 @@
|
||||
#ifndef MLD_C_DPREC_
|
||||
#define MLD_C_DPREC_
|
||||
|
||||
#include "mld_const.h"
|
||||
#include "psb_base_cbind.h"
|
||||
#include "psb_prec_cbind.h"
|
||||
#include "psb_krylov_cbind.h"
|
||||
|
||||
/* Object handle related routines */
|
||||
/* Note: mld_get_XXX_handle returns: <= 0 unsuccessful */
|
||||
/* >0 valid handle */
|
||||
#ifdef __cplusplus
|
||||
extern "C" {
|
||||
#endif
|
||||
typedef struct MLD_C_DPREC {
|
||||
void *dprec;
|
||||
} mld_c_dprec;
|
||||
|
||||
mld_c_dprec* mld_c_new_dprec();
|
||||
psb_i_t mld_c_delete_dprec(mld_c_dprec* p);
|
||||
|
||||
psb_i_t mld_c_dprecinit(mld_c_dprec *ph, const char *ptype);
|
||||
psb_i_t mld_c_dprecseti(mld_c_dprec *ph, const char *what, psb_i_t val);
|
||||
psb_i_t mld_c_dprecsetc(mld_c_dprec *ph, const char *what, const char *val);
|
||||
psb_i_t mld_c_dprecsetr(mld_c_dprec *ph, const char *what, double val);
|
||||
psb_i_t mld_c_dprecbld(psb_c_dspmat *ah, psb_c_descriptor *cdh, mld_c_dprec *ph);
|
||||
psb_i_t mld_c_dhierarchy_build(psb_c_dspmat *ah, psb_c_descriptor *cdh, mld_c_dprec *ph);
|
||||
psb_i_t mld_c_dsmoothers_build(psb_c_dspmat *ah, psb_c_descriptor *cdh, mld_c_dprec *ph);
|
||||
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_dkrylov(const char *method, psb_c_dspmat *ah, mld_c_dprec *ph,
|
||||
psb_c_dvector *bh, psb_c_dvector *xh,
|
||||
psb_c_descriptor *cdh, psb_c_SolverOptions *opt);
|
||||
|
||||
|
||||
#ifdef __cplusplus
|
||||
}
|
||||
#endif
|
||||
|
||||
#endif
|
||||
@ -0,0 +1,10 @@
|
||||
#ifndef MLD_CBIND_
|
||||
#define MLD_CBIND_
|
||||
|
||||
#define MLD_VALID_PRECONDITIONER_STRINGS "NONE","DIAG","BJAC","ML","AS"
|
||||
#define MLD_VALID_PRECONDITIONER_STRING "NONE DIAG BJAC ML AS"
|
||||
|
||||
#include "mld_const.h"
|
||||
#include "mld_c_dprec.h"
|
||||
|
||||
#endif
|
||||
@ -0,0 +1,14 @@
|
||||
/* This file was generated by a script using the mld_base_prec_type.F90 file as a basis. */
|
||||
#ifndef MLD_CONST_H_
|
||||
#define MLD_CONST_H_
|
||||
#ifdef __cplusplus
|
||||
extern "C" {
|
||||
#endif
|
||||
#define MLD_VERSION_STRING_ ( "2.3.0" )
|
||||
#define MLD_VERSION_MAJOR_ ( 2 )
|
||||
#define MLD_VERSION_MINOR_ ( 3 )
|
||||
#define MLD_PATCHLEVEL_ ( 0 )
|
||||
#ifdef __cplusplus
|
||||
}
|
||||
#endif
|
||||
#endif
|
||||
@ -0,0 +1,381 @@
|
||||
module mld_dprec_cbind_mod
|
||||
|
||||
use iso_c_binding
|
||||
use mld_prec_mod
|
||||
use psb_base_cbind_mod
|
||||
|
||||
type, bind(c) :: mld_c_dprec
|
||||
type(c_ptr) :: item = c_null_ptr
|
||||
end type mld_c_dprec
|
||||
|
||||
contains
|
||||
|
||||
#if 1
|
||||
#define MLDC_DEBUG(MSG) write(*,*) __FILE__,':',__LINE__,':',MSG
|
||||
#define MLDC_ERROR(MSG) write(*,*) __FILE__,':',__LINE__,':'," ERROR: ",MSG
|
||||
#else
|
||||
#define MLDC_DEBUG(MSG)
|
||||
#define MLDC_ERROR(MSG)
|
||||
#endif
|
||||
#define mld_success_ 0
|
||||
!#define MLDC_ERR_FILTER(INFO) min(0,INFO)
|
||||
#define MLDC_ERR_FILTER(INFO) (INFO)
|
||||
#define MLDC_ERR_HANDLE(INFO) if(INFO/=mld_success_)MLDC_ERROR("ERROR!")
|
||||
|
||||
function mld_c_dprecinit(ictxt,ph,ptype) bind(c) result(res)
|
||||
use psb_base_mod
|
||||
use mld_prec_mod
|
||||
implicit none
|
||||
|
||||
integer(psb_c_ipk_) :: res
|
||||
type(mld_c_dprec) :: ph
|
||||
integer(psb_c_ipk_), value :: ictxt
|
||||
character(c_char) :: ptype(*)
|
||||
integer :: info
|
||||
type(mld_dprec_type), pointer :: precp
|
||||
character(len=80) :: fptype
|
||||
|
||||
res = -1
|
||||
res = -1
|
||||
if (c_associated(ph%item)) then
|
||||
return
|
||||
end if
|
||||
|
||||
allocate(precp,stat=info)
|
||||
if (info /= 0) return
|
||||
|
||||
ph%item = c_loc(precp)
|
||||
|
||||
call stringc2f(ptype,fptype)
|
||||
|
||||
call precp%init(ictxt,fptype,info)
|
||||
|
||||
res = MLDC_ERR_FILTER(info)
|
||||
MLDC_ERR_HANDLE(res)
|
||||
return
|
||||
end function mld_c_dprecinit
|
||||
|
||||
function mld_c_dprecseti(ph,what,val) 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
|
||||
character(c_char) :: what(*)
|
||||
integer(psb_c_ipk_), value :: val
|
||||
integer :: info
|
||||
character(len=80) :: fwhat
|
||||
type(mld_dprec_type), pointer :: precp
|
||||
|
||||
res = -1
|
||||
if (c_associated(ph%item)) then
|
||||
call c_f_pointer(ph%item,precp)
|
||||
else
|
||||
return
|
||||
end if
|
||||
|
||||
call stringc2f(what,fwhat)
|
||||
|
||||
call mld_precset(precp,fwhat,val,info)
|
||||
|
||||
res = MLDC_ERR_FILTER(info)
|
||||
MLDC_ERR_HANDLE(res)
|
||||
return
|
||||
end function mld_c_dprecseti
|
||||
|
||||
|
||||
function mld_c_dprecsetr(ph,what,val) 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
|
||||
character(c_char) :: what(*)
|
||||
real(c_double), value :: val
|
||||
integer :: info
|
||||
character(len=80) :: fwhat
|
||||
type(mld_dprec_type), pointer :: precp
|
||||
|
||||
res = -1
|
||||
if (c_associated(ph%item)) then
|
||||
call c_f_pointer(ph%item,precp)
|
||||
else
|
||||
return
|
||||
end if
|
||||
|
||||
call stringc2f(what,fwhat)
|
||||
|
||||
call mld_precset(precp,fwhat,val,info)
|
||||
|
||||
res = MLDC_ERR_FILTER(info)
|
||||
MLDC_ERR_HANDLE(res)
|
||||
return
|
||||
end function mld_c_dprecsetr
|
||||
|
||||
function mld_c_dprecsetc(ph,what,val) 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
|
||||
character(c_char) :: what(*), val(*)
|
||||
integer :: info
|
||||
character(len=80) :: fwhat,fval
|
||||
type(mld_dprec_type), pointer :: precp
|
||||
|
||||
res = -1
|
||||
if (c_associated(ph%item)) then
|
||||
call c_f_pointer(ph%item,precp)
|
||||
else
|
||||
return
|
||||
end if
|
||||
|
||||
call stringc2f(what,fwhat)
|
||||
call stringc2f(val,fval)
|
||||
|
||||
call mld_precset(precp,fwhat,fval,info)
|
||||
|
||||
res = MLDC_ERR_FILTER(info)
|
||||
MLDC_ERR_HANDLE(res)
|
||||
return
|
||||
end function mld_c_dprecsetc
|
||||
|
||||
function mld_c_dprecbld(ah,cdh,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,ah,cdh
|
||||
integer :: info
|
||||
type(mld_dprec_type), pointer :: precp
|
||||
type(psb_dspmat_type), pointer :: ap
|
||||
type(psb_desc_type), pointer :: descp
|
||||
character(len=80) :: fptype
|
||||
|
||||
res = -1
|
||||
|
||||
if (c_associated(cdh%item)) then
|
||||
call c_f_pointer(cdh%item,descp)
|
||||
else
|
||||
return
|
||||
end if
|
||||
if (c_associated(ah%item)) then
|
||||
call c_f_pointer(ah%item,ap)
|
||||
else
|
||||
return
|
||||
end if
|
||||
if (c_associated(ph%item)) then
|
||||
call c_f_pointer(ph%item,precp)
|
||||
else
|
||||
return
|
||||
end if
|
||||
|
||||
call mld_precbld(ap,descp,precp,info)
|
||||
|
||||
res = MLDC_ERR_FILTER(info)
|
||||
MLDC_ERR_HANDLE(res)
|
||||
|
||||
return
|
||||
end function mld_c_dprecbld
|
||||
|
||||
function mld_c_dhierarchy_build(ah,cdh,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,ah,cdh
|
||||
integer :: info
|
||||
type(mld_dprec_type), pointer :: precp
|
||||
type(psb_dspmat_type), pointer :: ap
|
||||
type(psb_desc_type), pointer :: descp
|
||||
character(len=80) :: fptype
|
||||
|
||||
res = -1
|
||||
|
||||
if (c_associated(cdh%item)) then
|
||||
call c_f_pointer(cdh%item,descp)
|
||||
else
|
||||
return
|
||||
end if
|
||||
if (c_associated(ah%item)) then
|
||||
call c_f_pointer(ah%item,ap)
|
||||
else
|
||||
return
|
||||
end if
|
||||
if (c_associated(ph%item)) then
|
||||
call c_f_pointer(ph%item,precp)
|
||||
else
|
||||
return
|
||||
end if
|
||||
|
||||
call precp%hierarchy_build(ap,descp,info)
|
||||
|
||||
res = MLDC_ERR_FILTER(info)
|
||||
MLDC_ERR_HANDLE(res)
|
||||
|
||||
return
|
||||
end function mld_c_dhierarchy_build
|
||||
|
||||
function mld_c_dsmoothers_build(ah,cdh,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,ah,cdh
|
||||
integer :: info
|
||||
type(mld_dprec_type), pointer :: precp
|
||||
type(psb_dspmat_type), pointer :: ap
|
||||
type(psb_desc_type), pointer :: descp
|
||||
character(len=80) :: fptype
|
||||
|
||||
res = -1
|
||||
|
||||
if (c_associated(cdh%item)) then
|
||||
call c_f_pointer(cdh%item,descp)
|
||||
else
|
||||
return
|
||||
end if
|
||||
if (c_associated(ah%item)) then
|
||||
call c_f_pointer(ah%item,ap)
|
||||
else
|
||||
return
|
||||
end if
|
||||
if (c_associated(ph%item)) then
|
||||
call c_f_pointer(ph%item,precp)
|
||||
else
|
||||
return
|
||||
end if
|
||||
|
||||
call precp%smoothers_build(ap,descp,info)
|
||||
|
||||
res = MLDC_ERR_FILTER(info)
|
||||
MLDC_ERR_HANDLE(res)
|
||||
|
||||
return
|
||||
end function mld_c_dsmoothers_build
|
||||
|
||||
function mld_c_dkrylov(methd,&
|
||||
& ah,ph,bh,xh,cdh,options) bind(c) result(res)
|
||||
use psb_base_mod
|
||||
use psb_prec_mod
|
||||
use psb_krylov_mod
|
||||
use psb_prec_cbind_mod
|
||||
use psb_dkrylov_cbind_mod
|
||||
implicit none
|
||||
integer(psb_c_ipk_) :: res
|
||||
type(psb_c_object_type) :: ah,cdh,ph,bh,xh
|
||||
character(c_char) :: methd(*)
|
||||
type(solveroptions) :: options
|
||||
|
||||
res= mld_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)
|
||||
|
||||
end function mld_c_dkrylov
|
||||
|
||||
|
||||
function mld_c_dkrylov_opt(methd,&
|
||||
& ah,ph,bh,xh,eps,cdh,itmax,iter,err,itrace,irst,istop) bind(c) result(res)
|
||||
use psb_base_mod
|
||||
use psb_prec_mod
|
||||
use psb_krylov_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
|
||||
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(mld_dprec_type), pointer :: precp
|
||||
type(psb_d_vect_type), pointer :: xp, bp
|
||||
|
||||
integer :: info,fitmax,fitrace,first,fistop,fiter
|
||||
character(len=20) :: fmethd
|
||||
real(kind(1.d0)) :: feps,ferr
|
||||
|
||||
res = -1
|
||||
if (c_associated(cdh%item)) then
|
||||
call c_f_pointer(cdh%item,descp)
|
||||
else
|
||||
return
|
||||
end if
|
||||
if (c_associated(xh%item)) then
|
||||
call c_f_pointer(xh%item,xp)
|
||||
else
|
||||
return
|
||||
end if
|
||||
if (c_associated(bh%item)) then
|
||||
call c_f_pointer(bh%item,bp)
|
||||
else
|
||||
return
|
||||
end if
|
||||
if (c_associated(ah%item)) then
|
||||
call c_f_pointer(ah%item,ap)
|
||||
else
|
||||
return
|
||||
end if
|
||||
if (c_associated(ph%item)) then
|
||||
call c_f_pointer(ph%item,precp)
|
||||
else
|
||||
return
|
||||
end if
|
||||
|
||||
|
||||
call stringc2f(methd,fmethd)
|
||||
feps = eps
|
||||
fitmax = itmax
|
||||
fitrace = itrace
|
||||
first = irst
|
||||
fistop = istop
|
||||
|
||||
call psb_krylov(fmethd, ap, precp, bp, xp, feps, &
|
||||
& descp, info,&
|
||||
& itmax=fitmax,iter=fiter,itrace=fitrace,istop=fistop,&
|
||||
& irst=first, err=ferr)
|
||||
iter = fiter
|
||||
err = ferr
|
||||
res = min(info,0)
|
||||
|
||||
end function mld_c_dkrylov_opt
|
||||
|
||||
function mld_c_dprecfree(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
|
||||
character(len=80) :: fptype
|
||||
|
||||
res = -1
|
||||
if (c_associated(ph%item)) then
|
||||
call c_f_pointer(ph%item,precp)
|
||||
else
|
||||
return
|
||||
end if
|
||||
|
||||
|
||||
call precp%free(info)
|
||||
|
||||
res = MLDC_ERR_FILTER(info)
|
||||
MLDC_ERR_HANDLE(res)
|
||||
return
|
||||
end function mld_c_dprecfree
|
||||
|
||||
end module mld_dprec_cbind_mod
|
||||
|
||||
@ -0,0 +1,8 @@
|
||||
module mld_prec_cbind_mod
|
||||
|
||||
use iso_c_binding
|
||||
use psb_base_cbind_mod
|
||||
use psb_prec_cbind_mod
|
||||
use mld_dprec_cbind_mod
|
||||
|
||||
end module mld_prec_cbind_mod
|
||||
@ -0,0 +1,53 @@
|
||||
include ../../Make.inc
|
||||
#
|
||||
# Libraries used
|
||||
#
|
||||
LIBDIR=../../lib/
|
||||
INCLUDEDIR=../../include/
|
||||
HERE=.
|
||||
|
||||
FINCLUDES=$(FMFLAG). $(FMFLAG)$(LIBDIR) $(FMFLAG)$(PSBLAS_INCDIR)
|
||||
#PSBLAS_LIBS= -L$(PSBLAS_LIBDIR) -L$(LIBDIR) $(CPSBLAS_LIB) $(PSBLAS_LIB)
|
||||
# -lpsb_krylov_cbind -lpsb_prec_cbind -lpsb_base_cbind
|
||||
PSBC_LIBS= -L$(PSBLAS_CBIND_LIBDIR) \
|
||||
-lpsb_krylov_cbind -lpsb_prec_cbind -lpsb_base_cbind
|
||||
MLDC_LIBS=-L$(LIBDIR) -lmld_prec_cbind -lmld_prec -L$(MLDLIBDIR)
|
||||
#
|
||||
# Compilers and such
|
||||
#
|
||||
CCOPT= -g
|
||||
FINCLUDES=$(FMFLAG)$(LIBDIR) $(FMFLAG).
|
||||
CINCLUDES=-I$(LIBDIR) -I$(PSBLAS_CBIND_INCDIR) -I$(INCLUDEDIR)
|
||||
|
||||
EXEDIR=./runs
|
||||
# 20110404 specifying UMFLIBS here is not portable!
|
||||
#UMFLIBS=-lumfpack -lamd -lcholmod -lcolamd -lcamd -lccolamd -L/usr/include/suitesparse
|
||||
#UMFFLAGS=-DHave_UMF_ -I/usr/include/suitesparse
|
||||
|
||||
all: mldec
|
||||
|
||||
mldec: mldec.o
|
||||
$(MPF90) mldec.o -o mldec $(MLDC_LIBS) $(PSBC_LIBS) $(PSBCLDLIBS) $(PSBLAS_LIBS) \
|
||||
$(UMFLIBS) $(PSBLDLIBS) -lm -lgfortran
|
||||
# \
|
||||
# -lifcore -lifcoremt -lguide -limf -lirc -lintlc -lcxaguard -L/opt/intel/fc/10.0.023/lib/ -lm
|
||||
|
||||
/bin/mv mldec $(EXEDIR)
|
||||
|
||||
.f90.o:
|
||||
$(MPF90) $(F90COPT) $(FINCLUDES) $(FDEFINES) -c $<
|
||||
.c.o:
|
||||
$(MPCC) $(CCOPT) $(CINCLUDES) $(CDEFINES) -c $<
|
||||
|
||||
|
||||
clean:
|
||||
/bin/rm -f mldec.o $(EXEDIR)/mldec
|
||||
verycleanlib:
|
||||
(cd ../..; make veryclean)
|
||||
lib:
|
||||
(cd ../../; make library)
|
||||
|
||||
tests: all
|
||||
cd runs ; ./mldec < mlde.inp
|
||||
|
||||
|
||||
@ -0,0 +1,421 @@
|
||||
/*----------------------------------------------------------------------------------*/
|
||||
/* Parallel Sparse BLAS v2.2 */
|
||||
/* (C) Copyright 2007 Salvatore Filippone University of Rome Tor Vergata */
|
||||
/* */
|
||||
/* Redistribution and use in source and binary forms, with or without */
|
||||
/* modification, are permitted provided that the following conditions */
|
||||
/* are met: */
|
||||
/* 1. Redistributions of source code must retain the above copyright */
|
||||
/* notice, this list of conditions and the following disclaimer. */
|
||||
/* 2. Redistributions in binary form must reproduce the above copyright */
|
||||
/* notice, this list of conditions, and the following disclaimer in the */
|
||||
/* documentation and/or other materials provided with the distribution. */
|
||||
/* 3. The name of the PSBLAS group or the names of its contributors may */
|
||||
/* not be used to endorse or promote products derived from this */
|
||||
/* software without specific written permission. */
|
||||
/* */
|
||||
/* THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS */
|
||||
/* ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED */
|
||||
/* TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR */
|
||||
/* PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE PSBLAS GROUP OR ITS CONTRIBUTORS */
|
||||
/* BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR */
|
||||
/* CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF */
|
||||
/* SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS */
|
||||
/* INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN */
|
||||
/* CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) */
|
||||
/* ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE */
|
||||
/* POSSIBILITY OF SUCH DAMAGE. */
|
||||
/* */
|
||||
/* */
|
||||
/* File: ppdec.c */
|
||||
/* */
|
||||
/* Program: ppdec */
|
||||
/* This sample program shows how to build and solve a sparse linear */
|
||||
/* */
|
||||
/* The program solves a linear system based on the partial differential */
|
||||
/* equation */
|
||||
/* */
|
||||
/* */
|
||||
/* */
|
||||
/* The equation generated is */
|
||||
/* */
|
||||
/* b1 d d (u) b2 d d (u) a1 d (u)) a2 d (u))) */
|
||||
/* - ------ - ------ + ----- + ------ + a3 u = 0 */
|
||||
/* dx dx dy dy dx dy */
|
||||
/* */
|
||||
/* */
|
||||
/* with Dirichlet boundary conditions on the unit cube */
|
||||
/* */
|
||||
/* 0<=x,y,z<=1 */
|
||||
/* */
|
||||
/* The equation is discretized with finite differences and uniform stepsize; */
|
||||
/* the resulting discrete equation is */
|
||||
/* */
|
||||
/* ( u(x,y,z)(2b1+2b2+a1+a2)+u(x-1,y)(-b1-a1)+u(x,y-1)(-b2-a2)+ */
|
||||
/* -u(x+1,y)b1-u(x,y+1)b2)*(1/h**2) */
|
||||
/* */
|
||||
/* Example adapted from: C.T.Kelley */
|
||||
/* Iterative Methods for Linear and Nonlinear Equations */
|
||||
/* SIAM 1995 */
|
||||
/* */
|
||||
/* */
|
||||
/* In this sample program the index space of the discretized */
|
||||
/* computational domain is first numbered sequentially in a standard way, */
|
||||
/* then the corresponding vector is distributed according to an HPF BLOCK */
|
||||
/* distribution directive. */
|
||||
/* */
|
||||
/* Boundary conditions are set in a very simple way, by adding */
|
||||
/* equations of the form */
|
||||
/* */
|
||||
/* u(x,y) = rhs(x,y) */
|
||||
/* */
|
||||
/*----------------------------------------------------------------------------------*/
|
||||
|
||||
#include <stdio.h>
|
||||
#include <stdlib.h>
|
||||
#include <string.h>
|
||||
#include <math.h>
|
||||
|
||||
#include "psb_cbind.h"
|
||||
#include "mld_cbind.h"
|
||||
|
||||
#define LINEBUFSIZE 1024
|
||||
#define NBMAX 20
|
||||
|
||||
double a1(double x, double y, double z)
|
||||
{
|
||||
return(1.0);
|
||||
}
|
||||
double a2(double x, double y, double z)
|
||||
{
|
||||
return(20.0*y);
|
||||
}
|
||||
double a3(double x, double y, double z)
|
||||
{
|
||||
return(1.0);
|
||||
}
|
||||
double a4(double x, double y, double z)
|
||||
{
|
||||
return(1.0);
|
||||
}
|
||||
double b1(double x, double y, double z)
|
||||
{
|
||||
return(1.0);
|
||||
}
|
||||
double b2(double x, double y, double z)
|
||||
{
|
||||
return(1.0);
|
||||
}
|
||||
double b3(double x, double y, double z)
|
||||
{
|
||||
return(1.0);
|
||||
}
|
||||
|
||||
int matgen(int ictxt, int ng,int idim,int vg[],psb_c_dspmat *ah,psb_c_descriptor *cdh,
|
||||
psb_c_dvector *xh, psb_c_dvector *bh, psb_c_dvector *rh)
|
||||
{
|
||||
int iam, np;
|
||||
int x, y, z, el,glob_row,i,info,ret;
|
||||
double gx, gy, gz, deltah;
|
||||
double val[10*NBMAX], zt[NBMAX];
|
||||
int irow[10*NBMAX], icol[10*NBMAX];
|
||||
|
||||
info = 0;
|
||||
psb_c_info(ictxt,&iam,&np);
|
||||
deltah = (double) 1.0/(idim-1);
|
||||
|
||||
for (glob_row=1; glob_row<=ng; glob_row++) {
|
||||
|
||||
/* Check if I have to do something about this entry */
|
||||
if (vg[glob_row-1] == iam) {
|
||||
el=0;
|
||||
if ( (glob_row%(idim*idim)) == 0) {
|
||||
x = glob_row/(idim*idim);
|
||||
} else {
|
||||
x = glob_row/(idim*idim)+1;
|
||||
}
|
||||
if (((glob_row-(x-1)*idim*idim)%idim) == 0) {
|
||||
y = (glob_row-(x-1)*idim*idim)/idim;
|
||||
} else {
|
||||
y = (glob_row-(x-1)*idim*idim)/idim+1;
|
||||
}
|
||||
z = glob_row-(x-1)*idim*idim-(y-1)*idim;
|
||||
gx=x*deltah;
|
||||
gy=y*deltah;
|
||||
gz=z*deltah;
|
||||
zt[0] = 0.0;
|
||||
/* internal point: build discretization */
|
||||
/* term depending on (x-1,y,z) */
|
||||
|
||||
if (x==1) {
|
||||
val[el] = -b1(gx,gy,gz)-a1(gx,gy,gz);
|
||||
val[el] /= deltah*deltah;
|
||||
zt[0] = exp(-gy*gy-gz*gz)*(-val[el]);
|
||||
} else {
|
||||
val[el]=-b1(gx,gy,gz) -a1(gx,gy,gz);
|
||||
val[el] = val[el]/(deltah*deltah);
|
||||
icol[el]=(x-2)*idim*idim+(y-1)*idim+(z);
|
||||
el=el+1;
|
||||
}
|
||||
/* term depending on (x,y-1,z) */
|
||||
if (y==1) {
|
||||
val[el]=-b2(gx,gy,gz)-a2(gx,gy,gz);
|
||||
val[el] = val[el]/(deltah*deltah);
|
||||
zt[0] = exp(-gy*gy-gz*gz)*exp(-gx)*(-val[el]);
|
||||
} else {
|
||||
val[el]=-b2(gx,gy,gz)-a2(gx,gy,gz);
|
||||
val[el] = val[el]/(deltah*deltah);
|
||||
icol[el]=(x-1)*idim*idim+(y-2)*idim+(z);
|
||||
el=el+1;
|
||||
}
|
||||
/* term depending on (x,y,z-1)*/
|
||||
if (z==1) {
|
||||
val[el]=-b3(gx,gy,gz)-a3(gx,gy,gz);
|
||||
val[el] = val[el]/(deltah*deltah);
|
||||
zt[0] = exp(-gy*gy-gz*gz)*exp(-gx)*(-val[el]);
|
||||
} else {
|
||||
val[el]=-b3(gx,gy,gz)-a3(gx,gy,gz);
|
||||
val[el] = val[el]/(deltah*deltah);
|
||||
icol[el]=(x-1)*idim*idim+(y-1)*idim+(z-1);
|
||||
el=el+1;
|
||||
}
|
||||
/* term depending on (x,y,z)*/
|
||||
val[el]=2*b1(gx,gy,gz)+2*b2(gx,gy,gz)+2*b3(gx,gy,gz)
|
||||
+a1(gx,gy,gz)+a2(gx,gy,gz)+a3(gx,gy,gz);
|
||||
val[el] = val[el]/(deltah*deltah);
|
||||
icol[el]=(x-1)*idim*idim+(y-1)*idim+(z);
|
||||
el=el+1;
|
||||
/* term depending on (x,y,z+1) */
|
||||
if (z==idim) {
|
||||
val[el]=-b1(gx,gy,gz);
|
||||
val[el] = val[el]/(deltah*deltah);
|
||||
zt[0] = exp(-gy*gy-gz*gz)*exp(-gx)*(-val[el]);
|
||||
} else {
|
||||
val[el]=-b1(gx,gy,gz);
|
||||
val[el] = val[el]/(deltah*deltah);
|
||||
icol[el]=(x-1)*idim*idim+(y-1)*idim+(z+1);
|
||||
el=el+1;
|
||||
}
|
||||
/* term depending on (x,y+1,z) */
|
||||
if (y==idim) {
|
||||
val[el]=-b2(gx,gy,gz);
|
||||
val[el] = val[el]/(deltah*deltah);
|
||||
zt[0] = exp(-gy*gy-gz*gz)*exp(-gx)*(-val[el]);
|
||||
} else {
|
||||
val[el]=-b2(gx,gy,gz);
|
||||
val[el] = val[el]/(deltah*deltah);
|
||||
icol[el]=(x-1)*idim*idim+(y)*idim+(z);
|
||||
el=el+1;
|
||||
}
|
||||
/* term depending on (x+1,y,z) */
|
||||
if (x<idim) {
|
||||
val[el]=-b3(gx,gy,gz);
|
||||
val[el] = val[el]/(deltah*deltah);
|
||||
icol[el]=(x)*idim*idim+(y-1)*idim+(z);
|
||||
el=el+1;
|
||||
}
|
||||
for (i=0; i<el; i++) irow[i]=glob_row;
|
||||
if ((ret=psb_c_dspins(el,irow,icol,val,ah,cdh))!=0)
|
||||
fprintf(stderr,"From psb_c_dspins: %d\n",ret);
|
||||
psb_c_dgeins(1,irow,zt,bh,cdh);
|
||||
zt[0]=0.0;
|
||||
psb_c_dgeins(1,irow,zt,xh,cdh);
|
||||
}
|
||||
}
|
||||
|
||||
if ((info=psb_c_cdasb(cdh))!=0) return(info);
|
||||
|
||||
if ((info=psb_c_dspasb(ah,cdh))!=0) return(info);
|
||||
|
||||
if ((info=psb_c_dgeasb(xh,cdh))!=0) return(info);
|
||||
if ((info=psb_c_dgeasb(bh,cdh))!=0) return(info);
|
||||
if ((info=psb_c_dgeasb(rh,cdh))!=0) return(info);
|
||||
return(info);
|
||||
|
||||
}
|
||||
|
||||
int main(int argc, char *argv[])
|
||||
{
|
||||
int ictxt, iam, np;
|
||||
char methd[40], ptype[20], afmt[8], buffer[LINEBUFSIZE+1];
|
||||
int nparms;
|
||||
int idim,info,istop,itmax,itrace,irst,i,iter,ret;
|
||||
mld_c_dprec *ph;
|
||||
psb_c_dspmat *ah;
|
||||
psb_c_dvector *bh, *xh, *rh;
|
||||
int *vg, ng, nb,nlr;
|
||||
double t1,t2,eps,err;
|
||||
double *xv, *bv, *rv;
|
||||
double one=1.0, zero=0.0, res2;
|
||||
psb_c_SolverOptions options;
|
||||
psb_c_descriptor *cdh;
|
||||
|
||||
|
||||
ictxt = psb_c_init();
|
||||
psb_c_info(ictxt,&iam,&np);
|
||||
fprintf(stdout,"Initialization: am %d of %d\n",iam,np);
|
||||
|
||||
fflush(stdout);
|
||||
psb_c_barrier(ictxt);
|
||||
if (iam == 0) {
|
||||
fgets(buffer,LINEBUFSIZE,stdin);
|
||||
sscanf(buffer,"%d ",&nparms);
|
||||
fgets(buffer,LINEBUFSIZE,stdin);
|
||||
sscanf(buffer,"%s",methd);
|
||||
fgets(buffer,LINEBUFSIZE,stdin);
|
||||
sscanf(buffer,"%s",ptype);
|
||||
fgets(buffer,LINEBUFSIZE,stdin);
|
||||
sscanf(buffer,"%s",afmt);
|
||||
fgets(buffer,LINEBUFSIZE,stdin);
|
||||
sscanf(buffer,"%d",&idim);
|
||||
fgets(buffer,LINEBUFSIZE,stdin);
|
||||
sscanf(buffer,"%d",&istop);
|
||||
fgets(buffer,LINEBUFSIZE,stdin);
|
||||
sscanf(buffer,"%d",&itmax);
|
||||
fgets(buffer,LINEBUFSIZE,stdin);
|
||||
sscanf(buffer,"%d",&itrace);
|
||||
fgets(buffer,LINEBUFSIZE,stdin);
|
||||
sscanf(buffer,"%d",&irst);
|
||||
}
|
||||
/* Now broadcast the values, and check they're OK */
|
||||
psb_c_ibcast(ictxt,1,&nparms,0);
|
||||
psb_c_hbcast(ictxt,methd,0);
|
||||
psb_c_hbcast(ictxt,ptype,0);
|
||||
psb_c_hbcast(ictxt,afmt,0);
|
||||
psb_c_ibcast(ictxt,1,&idim,0);
|
||||
psb_c_ibcast(ictxt,1,&istop,0);
|
||||
psb_c_ibcast(ictxt,1,&itmax,0);
|
||||
psb_c_ibcast(ictxt,1,&itrace,0);
|
||||
psb_c_ibcast(ictxt,1,&irst,0);
|
||||
|
||||
psb_c_barrier(ictxt);
|
||||
|
||||
cdh=psb_c_new_descriptor();
|
||||
|
||||
/* Simple minded BLOCK data distribution */
|
||||
ng = idim*idim*idim;
|
||||
nb = (ng+np-1)/np;
|
||||
if ((vg=malloc(ng*sizeof(int)))==NULL) {
|
||||
fprintf(stderr,"On %d: malloc failure\n",iam);
|
||||
psb_c_abort(ictxt);
|
||||
}
|
||||
for (i=0; i<ng; i++) {
|
||||
vg[i] = i/nb;
|
||||
}
|
||||
if ((info=psb_c_cdall_vg(ng,vg,ictxt,cdh))!=0) {
|
||||
fprintf(stderr,"From cdall: %d\nBailing out\n",info);
|
||||
psb_c_abort(ictxt);
|
||||
}
|
||||
|
||||
bh = psb_c_new_dvector();
|
||||
xh = psb_c_new_dvector();
|
||||
rh = psb_c_new_dvector();
|
||||
ah = psb_c_new_dspmat();
|
||||
|
||||
/* Allocate mem space for sparse matrix and vectors */
|
||||
ret=psb_c_dspall(ah,cdh);
|
||||
psb_c_dgeall(bh,cdh);
|
||||
psb_c_dgeall(xh,cdh);
|
||||
psb_c_dgeall(rh,cdh);
|
||||
/* Matrix generation */
|
||||
if (matgen(ictxt,ng,idim,vg,ah,cdh,xh,bh,rh) != 0) {
|
||||
fprintf(stderr,"Error during matrix build loop\n");
|
||||
psb_c_abort(ictxt);
|
||||
}
|
||||
psb_c_barrier(ictxt);
|
||||
/* Set up the preconditioner */
|
||||
ph = mld_c_new_dprec();
|
||||
mld_c_dprecinit(ph,ptype,3);
|
||||
mld_c_dprecseti(ph,"SMOOTHER_SWEEPS",2);
|
||||
mld_c_dprecseti(ph,"SUB_FILLIN",1);
|
||||
mld_c_dprecsetc(ph,"COARSE_SOLVE","BJAC");
|
||||
mld_c_dprecsetc(ph,"COARSE_SUBSOLVE","ILU");
|
||||
mld_c_dprecseti(ph,"COARSE_FILLIN",1);
|
||||
mld_c_dprecbld(ah,cdh,ph);
|
||||
|
||||
psb_c_barrier(ictxt);
|
||||
/* Set up the solver options */
|
||||
psb_c_DefaultSolverOptions(&options);
|
||||
options.eps = 1.e-9;
|
||||
options.itmax = itmax;
|
||||
options.irst = irst;
|
||||
options.itrace = 1;
|
||||
options.istop = istop;
|
||||
|
||||
t1=psb_c_wtime();
|
||||
ret=mld_c_dkrylov(methd,ah,ph,bh,xh,cdh,&options);
|
||||
t2=psb_c_wtime();
|
||||
iter = options.iter;
|
||||
err = options.err;
|
||||
fprintf(stderr,"From krylov: %d %lf, %d\n",iter,err,ret);
|
||||
|
||||
|
||||
/* Check 2-norm of residual on exit */
|
||||
psb_c_dgeaxpby(one,bh,zero,rh,cdh);
|
||||
psb_c_dspmm(-one,ah,xh,one,rh,cdh);
|
||||
res2=psb_c_dgenrm2(rh,cdh);
|
||||
|
||||
if (iam==0) {
|
||||
fprintf(stdout,"Time: %lf\n",(t2-t1));
|
||||
fprintf(stdout,"Iter: %d\n",iter);
|
||||
fprintf(stdout,"Err: %lg\n",err);
|
||||
fprintf(stdout,"||r||_2: %lg\n",res2);
|
||||
|
||||
}
|
||||
|
||||
#if 0
|
||||
bv = psb_c_dvect_get_cpy(bh);
|
||||
nlr=psb_c_cd_get_local_rows(cdh);
|
||||
for (i=0;i<nlr; i++)
|
||||
fprintf(stdout,"RHS: %d %d %lf\n",iam,i,bv[i]);
|
||||
|
||||
|
||||
xv = psb_c_dvect_get_cpy(xh);
|
||||
nlr=psb_c_cd_get_local_rows(cdh);
|
||||
for (i=0;i<nlr; i++)
|
||||
fprintf(stdout,"SOL: %d %d %lf\n",iam,i,xv[i]);
|
||||
|
||||
|
||||
rv = psb_c_dvect_get_cpy(rh);
|
||||
nlr=psb_c_cd_get_local_rows(cdh);
|
||||
for (i=0;i<nlr; i++)
|
||||
fprintf(stdout,"RES: %d %d %lf\n",iam,i,rv[i]);
|
||||
|
||||
#endif
|
||||
|
||||
/* Clean up memory */
|
||||
if ((info=psb_c_dgefree(xh,cdh))!=0) {
|
||||
fprintf(stderr,"From dgefree: %d\nBailing out\n",info);
|
||||
psb_c_abort(ictxt);
|
||||
}
|
||||
if ((info=psb_c_dgefree(bh,cdh))!=0) {
|
||||
fprintf(stderr,"From dgefree: %d\nBailing out\n",info);
|
||||
psb_c_abort(ictxt);
|
||||
}
|
||||
if ((info=psb_c_dgefree(rh,cdh))!=0) {
|
||||
fprintf(stderr,"From dgefree: %d\nBailing out\n",info);
|
||||
psb_c_abort(ictxt);
|
||||
}
|
||||
|
||||
if ((info=psb_c_cdfree(cdh))!=0) {
|
||||
fprintf(stderr,"From cdfree: %d\nBailing out\n",info);
|
||||
psb_c_abort(ictxt);
|
||||
}
|
||||
fprintf(stderr,"pointer from cdfree: %p\n",cdh->descriptor);
|
||||
|
||||
/* Clean up object handles */
|
||||
free(ph);
|
||||
free(xh);
|
||||
free(bh);
|
||||
free(ah);
|
||||
free(cdh);
|
||||
|
||||
|
||||
fprintf(stderr,"program completed successfully\n");
|
||||
|
||||
psb_c_barrier(ictxt);
|
||||
psb_c_exit(ictxt);
|
||||
}
|
||||
|
||||
|
||||
|
||||
@ -0,0 +1,11 @@
|
||||
7 Number of entries below this
|
||||
BICGSTAB Iterative method BICGSTAB CGS BICG BICGSTABL RGMRES
|
||||
ML Preconditioner NONE DIAG BJAC
|
||||
CSR A Storage format CSR COO JAD
|
||||
60 Domain size (acutal system is this**3)
|
||||
1 Stopping criterion
|
||||
80 MAXIT
|
||||
01 ITRACE
|
||||
20 IRST restart for RGMRES and BiCGSTABL
|
||||
|
||||
|
||||
Loading…
Reference in New Issue