From f9e52bed20f23f79c84355842151db5ce304b488 Mon Sep 17 00:00:00 2001 From: Mathieu Faverge <mathieu.faverge@inria.fr> Date: Mon, 14 May 2018 19:24:07 +0200 Subject: [PATCH] Add the python wrapper --- CMakeLists.txt | 4 + wrappers/CMakeLists.txt | 20 +++ wrappers/python/CMakeLists.txt | 50 +++++++ wrappers/python/spm/__init__.py | 48 +++++++ wrappers/python/spm/__spm__.py | 200 ++++++++++++++++++++++++++++ wrappers/python/spm/enum.py.in | 155 ++++++++++++++++++++++ wrappers/python/spm/spm.py | 227 ++++++++++++++++++++++++++++++++ wrappers/python/spm_driver.py | 42 ++++++ wrappers/python/spm_scipy.py | 36 +++++ 9 files changed, 782 insertions(+) create mode 100644 wrappers/CMakeLists.txt create mode 100644 wrappers/python/CMakeLists.txt create mode 100644 wrappers/python/spm/__init__.py create mode 100644 wrappers/python/spm/__spm__.py create mode 100644 wrappers/python/spm/enum.py.in create mode 100644 wrappers/python/spm/spm.py create mode 100755 wrappers/python/spm_driver.py create mode 100755 wrappers/python/spm_scipy.py diff --git a/CMakeLists.txt b/CMakeLists.txt index fdcd443e..577e96c4 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -267,3 +267,7 @@ add_documented_files( # Testing executables add_subdirectory(tests) + +### Wrappers +add_subdirectory(wrappers) + diff --git a/wrappers/CMakeLists.txt b/wrappers/CMakeLists.txt new file mode 100644 index 00000000..cf763161 --- /dev/null +++ b/wrappers/CMakeLists.txt @@ -0,0 +1,20 @@ +### +# +# @copyright 2017-2018 Bordeaux INP, CNRS (LaBRI UMR 5800), Inria, +# Univ. Bordeaux. All rights reserved. +# +# @version 6.0.0 +# @author Mathieu Faverge +# @date 2017-05-22 +# +### + +if (SPM_WITH_FORTRAN) + add_subdirectory( fortran90 ) +endif() + +if (BUILD_SHARED_LIBS) + add_subdirectory( python ) +else() + message(STATUS "--- Python wrapper is disabled with static libraries") +endif() diff --git a/wrappers/python/CMakeLists.txt b/wrappers/python/CMakeLists.txt new file mode 100644 index 00000000..a4b59201 --- /dev/null +++ b/wrappers/python/CMakeLists.txt @@ -0,0 +1,50 @@ +### +# +# @copyright 2017 Bordeaux INP, CNRS (LaBRI UMR 5800), Inria, +# Univ. Bordeaux. All rights reserved. +# +# @version 6.0.0 +# @author Mathieu Faverge +# @date 2017-05-22 +# +### + +# Configure enum.py +if (SPM_INT64) + set(SPM_PYTHON_INTEGER c_int64) +else() + set(SPM_PYTHON_INTEGER c_int) +endif() + +configure_file( + "${CMAKE_CURRENT_SOURCE_DIR}/spm/enum.py.in" + "${CMAKE_CURRENT_SOURCE_DIR}/spm/enum.py" @ONLY) + +# Install python package +install(FILES + ${CMAKE_CURRENT_SOURCE_DIR}/spm/__init__.py + ${CMAKE_CURRENT_SOURCE_DIR}/spm/__spm__.py + ${CMAKE_CURRENT_SOURCE_DIR}/spm/spm.py + ${CMAKE_CURRENT_SOURCE_DIR}/spm/enum.py + DESTINATION lib/python/spm ) + +# # Install python examples +# install(FILES +# ${CMAKE_CURRENT_SOURCE_DIR}/examples/simple.py +# ${CMAKE_CURRENT_SOURCE_DIR}/examples/simple_obj.py +# ${CMAKE_CURRENT_SOURCE_DIR}/examples/schur.py +# ${CMAKE_CURRENT_SOURCE_DIR}/examples/schur_obj.py +# DESTINATION examples +# ) + +# ## CTest execution +# find_package(PythonInterp QUIET) +# if (PYTHONINTERP_FOUND) +# set( PYTHON_TESTS +# simple step-by-step schur simple_obj schur_obj) + +# foreach(example ${PYTHON_TESTS} ) +# add_test(python_${example} ${PYTHON_EXECUTABLE} ${CMAKE_CURRENT_SOURCE_DIR}/examples/${example}.py) +# endforeach() +# endif() + diff --git a/wrappers/python/spm/__init__.py b/wrappers/python/spm/__init__.py new file mode 100644 index 00000000..caf07917 --- /dev/null +++ b/wrappers/python/spm/__init__.py @@ -0,0 +1,48 @@ +# +# @file __init__.py +# +# SParse Matrix package python module intialization +# +# @copyright 2017-2018 Bordeaux INP, CNRS (LaBRI UMR 5800), Inria, +# Univ. Bordeaux. All rights reserved. +# +# @version 6.0.0 +# @author Pierre Ramet +# @author Mathieu Faverge +# @author Louis Poirel +# @date 2017-05-04 +# +""" +PySpm +===== + +Provides + 1. A sparse matrix structure + 2. Mathematical operations over this sparse matrix structure + 3. Driver to read from different file formats and to convert from Scipy package + +""" +import ctypes +import ctypes.util + +# Load the SPM library +libspm_name = ctypes.util.find_library('spm') +if libspm_name == None: + raise EnvironmentError("Could not find shared library: spm." + "The path to libpastix_spm.so should be in " + "$LIBRARY_PATH") + +try: + libspm = ctypes.cdll.LoadLibrary(libspm_name) +except: + raise EnvironmentError("Could not load shared library: spm." + "The path to libpastix_spm.so should be in " + "$LD_LIBRARY_PATH or $DYLD_LIBRARY_PATH on MacOS"); + +__all__ = [ 'libspm' ] + +from .enum import * +from .spm import * + +#__all__.extend(enum.__all__) +#__all__.extend(spm.__all__) diff --git a/wrappers/python/spm/__spm__.py b/wrappers/python/spm/__spm__.py new file mode 100644 index 00000000..02e151fa --- /dev/null +++ b/wrappers/python/spm/__spm__.py @@ -0,0 +1,200 @@ +""" + + @file __spm__.py + + SPM python wrapper + + @copyright 2017-2017 Bordeaux INP, CNRS (LaBRI UMR 5800), Inria, + Univ. Bordeaux. All rights reserved. + + @version 6.0.0 + @author Pierre Ramet + @author Mathieu Faverge + @author Louis Poirel + @date 2017-05-04 + +This file has been automatically generated with gen_wrappers.py + +""" +from ctypes import * +import numpy as np + +from . import libspm +from .enum import __spm_int__ + +class pyspm_spmatrix_t(Structure): + _fields_ = [("mtxtype", c_int ), + ("flttype", c_int ), + ("fmttype", c_int ), + ("gN", __spm_int__ ), + ("n", __spm_int__ ), + ("gnnz", __spm_int__ ), + ("nnz", __spm_int__ ), + ("gNexp", __spm_int__ ), + ("nexp", __spm_int__ ), + ("gnnzexp", __spm_int__ ), + ("nnzexp", __spm_int__ ), + ("dof", __spm_int__ ), + ("dofs", POINTER(__spm_int__)), + ("layout", c_int ), + ("colptr", POINTER(__spm_int__)), + ("rowptr", POINTER(__spm_int__)), + ("loc2glob", POINTER(__spm_int__)), + ("values", c_void_p ) ] + +def pyspm_spmInit( spm ): + libspm.spmInit.argtypes = [ POINTER(pyspm_spmatrix_t) ] + libspm.spmInit( spm ) + +def pyspm_spmExit( spm ): + libspm.spmExit.argtypes = [ POINTER(pyspm_spmatrix_t) ] + libspm.spmExit( spm ) + +def pyspm_spmCopy( spm ): + libspm.spmCopy.argtypes = [ POINTER(pyspm_spmatrix_t) ] + libspm.spmCopy.restype = POINTER(pyspm_spmatrix_t) + return libspm.spmCopy( spm ) + +def pyspm_spmBase( spm, baseval ): + libspm.spmBase.argtypes = [ POINTER(pyspm_spmatrix_t), c_int ] + libspm.spmBase( spm, baseval ) + +def pyspm_spmFindBase( spm ): + libspm.spmFindBase.argtypes = [ POINTER(pyspm_spmatrix_t) ] + libspm.spmFindBase.restype = __spm_int__ + return libspm.spmFindBase( spm ) + +def pyspm_spmConvert( ofmttype, ospm ): + libspm.spmConvert.argtypes = [ c_int, POINTER(pyspm_spmatrix_t) ] + libspm.spmConvert.restype = c_int + return libspm.spmConvert( ofmttype, ospm ) + +def pyspm_spmUpdateComputedFields( spm ): + libspm.spmUpdateComputedFields.argtypes = [ POINTER(pyspm_spmatrix_t) ] + libspm.spmUpdateComputedFields( spm ) + +def pyspm_spmGenFakeValues( spm ): + libspm.spmGenFakeValues.argtypes = [ POINTER(pyspm_spmatrix_t) ] + libspm.spmGenFakeValues( spm ) + +def pyspm_spmNorm( ntype, spm ): + libspm.spmNorm.argtypes = [ c_int, POINTER(pyspm_spmatrix_t) ] + libspm.spmNorm.restype = c_double + return libspm.spmNorm( ntype, spm ) + +def pyspm_spmMatVec( trans, alpha, spm, x, beta, y ): + libspm.spmMatVec.argtypes = [ c_int, c_void_p, POINTER(pyspm_spmatrix_t), + c_void_p, c_void_p, c_void_p ] + libspm.spmMatVec.restype = c_int + return libspm.spmMatVec( trans, alpha, spm, x, beta, y ) + +def pyspm_spmMatMat( trans, n, alpha, A, B, ldb, beta, C, ldc ): + libspm.spmMatMat.argtypes = [ c_int, __spm_int__, c_void_p, + POINTER(pyspm_spmatrix_t), c_void_p, + __spm_int__, c_void_p, c_void_p, __spm_int__ ] + libspm.spmMatMat.restype = c_int + return libspm.spmMatMat( trans, n, alpha, A, B, ldb, beta, C, ldc ) + +def pyspm_spmScalMatrix( alpha, spm ): + libspm.spmScalMatrix.argtypes = [ c_double, POINTER(pyspm_spmatrix_t) ] + libspm.spmScalMatrix( alpha, spm ) + +def pyspm_spmScalVector( flt, alpha, n, x, incx ): + libspm.spmScalVector.argtypes = [ c_int, c_double, __spm_int__, c_void_p, + __spm_int__ ] + libspm.spmScalVector( flt, alpha, n, x, incx ) + +def pyspm_spmSort( spm ): + libspm.spmSort.argtypes = [ POINTER(pyspm_spmatrix_t) ] + libspm.spmSort.restype = c_int + return libspm.spmSort( spm ) + +def pyspm_spmMergeDuplicate( spm ): + libspm.spmMergeDuplicate.argtypes = [ POINTER(pyspm_spmatrix_t) ] + libspm.spmMergeDuplicate.restype = __spm_int__ + return libspm.spmMergeDuplicate( spm ) + +def pyspm_spmSymmetrize( spm ): + libspm.spmSymmetrize.argtypes = [ POINTER(pyspm_spmatrix_t) ] + libspm.spmSymmetrize.restype = __spm_int__ + return libspm.spmSymmetrize( spm ) + +def pyspm_spmCheckAndCorrect( spm ): + libspm.spmCheckAndCorrect.argtypes = [ POINTER(pyspm_spmatrix_t) ] + libspm.spmCheckAndCorrect.restype = POINTER(pyspm_spmatrix_t) + return libspm.spmCheckAndCorrect( spm ) + +def pyspm_spmGenRHS( type, nrhs, spm, x, ldx, b, ldb ): + libspm.spmGenRHS.argtypes = [ c_int, __spm_int__, POINTER(pyspm_spmatrix_t), + c_void_p, __spm_int__, c_void_p, __spm_int__ ] + libspm.spmGenRHS.restype = c_int + return libspm.spmGenRHS( type, nrhs, spm, x, ldx, b, ldb ) + +def pyspm_spmCheckAxb( eps, nrhs, spm, x0, ldx0, b, ldb, x, ldx ): + libspm.spmCheckAxb.argtypes = [ c_double, __spm_int__, + POINTER(pyspm_spmatrix_t), c_void_p, + __spm_int__, c_void_p, __spm_int__, + c_void_p, __spm_int__ ] + libspm.spmCheckAxb.restype = c_int + return libspm.spmCheckAxb( eps, nrhs, spm, x0, ldx0, b, ldb, x, ldx ) + +def pyspm_spmIntConvert( n, input ): + libspm.spmIntConvert.argtypes = [ __spm_int__, c_int_p ] + libspm.spmIntConvert.restype = POINTER(__spm_int__) + return libspm.spmIntConvert( n, input ) + +def pyspm_spmIntMSortIntAsc( pbase, n ): + libspm.spmIntMSortIntAsc.argtypes = [ c_void_p, __spm_int__ ] + libspm.spmIntMSortIntAsc( pointer( pbase ), n ) + +def pyspm_spmLoad( spm ): + libspm.spmLoad.argtypes = [ POINTER(pyspm_spmatrix_t), c_void_p ] + libspm.spmLoad.restype = c_int + return libspm.spmLoad( spm, None ) + +def pyspm_spmSave( spm ): + libspm.spmSave.argtypes = [ POINTER(pyspm_spmatrix_t), c_void_p ] + libspm.spmSave.restype = c_int + return libspm.spmSave( spm, None ) + +def pyspm_spmReadDriver( driver, filename, spm ): + libspm.spmReadDriver.argtypes = [ c_int, c_char_p, + POINTER(pyspm_spmatrix_t) ] + libspm.spmReadDriver.restype = c_int + return libspm.spmReadDriver( driver, filename, spm ) + +def pyspm_spmParseLaplacianInfo( filename, flttype, dim1, dim2, dim3, alpha, + beta ): + libspm.spmParseLaplacianInfo.argtypes = [ c_char_p, c_int_p, + POINTER(__spm_int__), + POINTER(__spm_int__), + POINTER(__spm_int__), + POINTER(c_double), + POINTER(c_double) ] + libspm.spmParseLaplacianInfo.restype = c_int + return libspm.spmParseLaplacianInfo( filename, flttype, dim1, dim2, dim3, + alpha, beta ) + +def pyspm_spm2Dense( spm ): + libspm.spm2Dense.argtypes = [ POINTER(pyspm_spmatrix_t) ] + libspm.spm2Dense( spm ) + +def pyspm_spmPrint( spm ): + libspm.spmPrint.argtypes = [ POINTER(pyspm_spmatrix_t), c_void_p ] + libspm.spmPrint( spm, None ) + +def pyspm_spmPrintInfo( spm ): + libspm.spmPrintInfo.argtypes = [ POINTER(pyspm_spmatrix_t), c_void_p ] + libspm.spmPrintInfo( spm, None ) + +def pyspm_spmExpand( spm ): + libspm.spmExpand.argtypes = [ POINTER(pyspm_spmatrix_t) ] + libspm.spmExpand.restype = POINTER(pyspm_spmatrix_t) + return libspm.spmExpand( spm ) + +def pyspm_spmDofExtend( spm, type, dof ): + libspm.spmDofExtend.argtypes = [ POINTER(pyspm_spmatrix_t), c_int, c_int ] + libspm.spmDofExtend.restype = POINTER(pyspm_spmatrix_t) + return libspm.spmDofExtend( spm, type, dof ) + + diff --git a/wrappers/python/spm/enum.py.in b/wrappers/python/spm/enum.py.in new file mode 100644 index 00000000..b73be1c7 --- /dev/null +++ b/wrappers/python/spm/enum.py.in @@ -0,0 +1,155 @@ +""" + + @file enum.py + + SPM python wrapper to define enums and datatypes + + @copyright 2017-2017 Bordeaux INP, CNRS (LaBRI UMR 5800), Inria, + Univ. Bordeaux. All rights reserved. + + @version 6.0.0 + @author Pierre Ramet + @author Mathieu Faverge + @author Louis Poirel + @date 2017-05-04 + +This file has been automatically generated with gen_wrappers.py + +""" +from ctypes import * +import numpy as np + +# Start with __ to prevent broadcast to file importing enum +__spm_int__ = @SPM_PYTHON_INTEGER@ + +class verbose: + Not = 0 + No = 1 + Yes = 2 + +class coeftype: + Pattern = 0 + Float = 2 + Double = 3 + Complex32 = 4 + Complex64 = 5 + + @staticmethod + def getptype ( dtype ): + np_dict = { + np.dtype('float32') : coeftype.Float, + np.dtype('float64') : coeftype.Double, + np.dtype('complex64') : coeftype.Complex32, + np.dtype('complex128') : coeftype.Complex64, + } + if dtype in np_dict: + return np_dict[dtype] + else: + return -1 + + @staticmethod + def getctype ( flttype ): + class c_float_complex(Structure): + _fields_ = [("r",c_float),("i", c_float)] + class c_double_complex(Structure): + _fields_ = [("r",c_double),("i", c_double)] + + np_dict = { + coeftype.Float : c_float, + coeftype.Double : c_double, + coeftype.Complex32 : c_float_complex, + coeftype.Complex64 : c_double_complex + } + if flttype in np_dict: + return np_dict[flttype] + else: + return -1 + + @staticmethod + def getnptype ( flttype ): + np_dict = { + coeftype.Float : np.dtype('float32'), + coeftype.Double : np.dtype('float64'), + coeftype.Complex32 : np.dtype('complex64'), + coeftype.Complex64 : np.dtype('complex128') + } + if flttype in np_dict: + return np_dict[flttype] + else: + return -1 + +class fmttype: + CSC = 0 + CSR = 1 + IJV = 2 + +class error: + SUCCESS = 0 + UNKNOWN = 1 + ALLOC = 2 + NOTIMPLEMENTED = 3 + OUTOFMEMORY = 4 + THREAD = 5 + INTERNAL = 6 + BADPARAMETER = 7 + FILE = 8 + INTEGER_TYPE = 9 + IO = 10 + MPI = 11 + +class driver: + RSA = 0 + HB = 1 + IJV = 2 + MM = 3 + Laplacian = 4 + XLaplacian = 5 + Graph = 6 + SPM = 7 + +class rhstype: + One = 0 + I = 1 + RndX = 2 + RndB = 3 + +class layout: + RowMajor = 101 + ColMajor = 102 + +class trans: + NoTrans = 111 + Trans = 112 + ConjTrans = 113 + +class mtxtype: + General = trans.NoTrans + Symmetric = trans.Trans + Hermitian = trans.ConjTrans + SymPosDef = trans.ConjTrans + 1 + HerPosDef = trans.ConjTrans + 2 + +class uplo: + Upper = 121 + Lower = 122 + UpperLower = 123 + +class diag: + NonUnit = 131 + Unit = 132 + +class side: + Left = 141 + Right = 142 + +class normtype: + One = 171 + Frobenius = 174 + Inf = 175 + Max = 177 + +class dir: + Forward = 391 + Backward = 392 + + diff --git a/wrappers/python/spm/spm.py b/wrappers/python/spm/spm.py new file mode 100644 index 00000000..875d9763 --- /dev/null +++ b/wrappers/python/spm/spm.py @@ -0,0 +1,227 @@ +""" + @file spm.py + + SPM python wrapper + + @copyright 2017-2018 Bordeaux INP, CNRS (LaBRI UMR 5800), Inria, + Univ. Bordeaux. All rights reserved. + + @version 6.0.0 + @author Pierre Ramet + @author Mathieu Faverge + @date 2017-05-04 + +""" +from ctypes import * +import numpy as np + +__use_sps__ = True +try: + import scipy.sparse as sps +except ImportError: + __use_sps__ = False + +from .enum import * +from .enum import __spm_int__ +from .__spm__ import * + +class spmatrix(): + + dtype = None + + def __init__( self, A=None, mtxtype_=mtxtype.General, driver=None, filename="" ): + """ + Initialize the SPM wrapper by loading the libraries + """ + if mtxtype_ == mtxtype.SymPosDef: + mtxtype_ = mtxtype.Symmetric + if mtxtype_ == mtxtype.HerPosDef: + mtxtype_ = mtxtype.Hermitian + + self.spm_c = pyspm_spmatrix_t( mtxtype_, + coeftype.Double, + fmttype.CSC, + 0, 0, 0, 0, 0, 0, 0, 0, + 1, None, + layout.ColMajor, + None, None, None, None ) + self.id_ptr = pointer( self.spm_c ) + if A is not None: + self.fromsps( A, mtxtype_ ) + elif driver is not None: + self.fromdriver( driver, filename ) + + if __use_sps__: + def fromsps( self, A, mtxtype_=mtxtype.General ): + """ + Initialize the SPM wrapper by loading the libraries + """ + + if not sps.isspmatrix(A): + raise TypeError( "A must be of type scipy.sparse matrix" ) + + # Assume A is already in Scipy sparse format + self.dtype = A.dtype + flttype = coeftype.getptype( A.dtype ) + if flttype == -1: + raise TypeError( "Invalid data type. Must be part of (f4, f8, c8 or c16)" ) + + A = sps.csc_matrix( A ) + A.sort_indices() + + # Pointer variables + self.py_colptr = np.array( A.indptr[:], dtype=__spm_int__ ) + self.py_rowptr = np.array( A.indices[:], dtype=__spm_int__ ) + self.py_values = np.array( A.data[:] ) + + if mtxtype_ == mtxtype.SymPosDef: + mtxtype_ = mtxtype.Symmetric + if mtxtype_ == mtxtype.HerPosDef: + mtxtype_ = mtxtype.Hermitian + + self.spm_c.mtxtype = mtxtype_ + self.spm_c.flttype = flttype + self.spm_c.fmttype = fmttype.CSC + self.spm_c.n = A.shape[0] + self.spm_c.nnz = A.getnnz() + self.spm_c.dof = 1 + self.spm_c.dofs = None + self.spm_c.layout = layout.ColMajor + self.spm_c.colptr = self.py_colptr.ctypes.data_as(POINTER(__spm_int__)) + self.spm_c.rowptr = self.py_rowptr.ctypes.data_as(POINTER(__spm_int__)) + self.spm_c.loc2glob = None + self.spm_c.values = self.py_values.ctypes.data_as(c_void_p) + + self.id_ptr = pointer( self.spm_c ) + + self.updateComputedField() + self.checkAndCorrect() + + def tosps( self ): + """ + Return a Scipy sparse matrix + """ + n = int( self.spm_c.n ) + nnz = int( self.spm_c.nnz ) + cflt = coeftype.getctype( self.spm_c.flttype ) + nflt = coeftype.getnptype( self.spm_c.flttype ) + colptr = self.py_colptr.copy() + rowptr = self.py_rowptr.copy() + values = self.py_values.copy() + + baseval = colptr[0] + colptr = colptr - baseval + rowptr = rowptr - baseval + + return sps.csc_matrix((values, rowptr, colptr), shape=(n, n)) + + def fromdriver( self, driver=driver.Laplacian, filename="10:10:10" ): + """ + Initialize the SPM wrapper by loading the libraries + """ + if filename == "": + raise ValueError("filename must be prodived") + + pyspm_spmReadDriver( driver, filename.encode('utf-8'), self.id_ptr ) + + self.dtype = coeftype.getnptype( self.spm_c.flttype ) + self.checkAndCorrect() + + def scale( self, alpha ): + return pyspm_spmScalMatrix( float(alpha), self.id_ptr ) + + def norm( self, ntype=normtype.Frobenius ): + return pyspm_spmNorm( ntype, self.id_ptr ) + + def base( self, baseval ): + pyspm_spmBase( self.id_ptr, baseval ) + + def findBase( self ): + return pyspm_spmFindBase( self.id_ptr ) + + def updateComputedField( self ): + pyspm_spmUpdateComputedFields( self.id_ptr ) + + def printInfo( self ): + pyspm_spmPrintInfo( self.id_ptr ) + + def printSpm( self ): + pyspm_spmPrint( self.id_ptr ) + + def checkAndCorrect( self ): + spm1 = self.id_ptr + spm2 = pyspm_spmCheckAndCorrect( self.id_ptr ) + if (( spm1.contents.fmttype == spm2.contents.fmttype ) and + ( spm1.contents.nnzexp == spm2.contents.nnzexp ) ): + return + self.spm_c = cast( spm2, POINTER(pyspm_spmatrix_t) ).contents + self.id_ptr = pointer( self.spm_c ) + + self.py_colptr = np.frombuffer( (__spm_int__ * (n+1)).from_address( cast(self.spm_c.colptr, c_void_p).value ), __spm_int__ ).copy() + self.py_rowptr = np.frombuffer( (__spm_int__ * nnz ).from_address( cast(self.spm_c.rowptr, c_void_p).value ), __spm_int__ ).copy() + self.py_values = np.frombuffer( (cflt * nnz ).from_address( self.spm_c.values ), nflt ).copy() + + def __checkVector( self, n, nrhs, x ): + if x.dtype != self.dtype: + raise TypeError( "Vectors must use the same arithmetic as the spm" ) + + if x.ndim > 2: + raise TypeError( "Vectors must be of dimension 1 or 2" ) + + if x.shape[0] < n: + raise TypeError( "Vectors must be of size at least ", n ) + + if (x.ndim == 1 and nrhs > 1) or (x.ndim>1 and x.shape[1] < nrhs): + raise TypeError( "At least nrhs vectors must be stored in the vector" ) + + def checkAxb( self, x0, b, x, nrhs=-1 ): + # if libspm == None: + # raise EnvironmentError( "SPM Instance badly instanciated" ) + + b = np.array(b, self.dtype) + x = np.asarray(x, self.dtype) + + n = self.spm_c.n + if nrhs == -1: + if x.ndim == 1: + nrhs = 1 + else: + nrhs = x.shape[1] + + if x0 is not None: + x0 = np.asarray(x0, self.dtype) + self.__checkVector( n, nrhs, x0 ) + ldx0 = x0.shape[0] + x0ptr = x0.ctypes.data_as(c_void_p) + else: + ldx0 = 1 + x0ptr = None + self.__checkVector( n, nrhs, b ) + self.__checkVector( n, nrhs, x ) + + ldb = b.shape[0] + ldx = x.shape[0] + + pyspm_spmCheckAxb( -1., nrhs, self.id_ptr, + x0ptr, ldx0, + b.ctypes.data_as(c_void_p), ldb, + x.ctypes.data_as(c_void_p), ldx ) + + def genRHS( self, rhstype=rhstype.One, nrhs=1 ): + # if libspm == None: + # raise EnvironmentError( "SPM Instance badly instanciated" ) + + n = self.spm_c.n + b = np.zeros((n, nrhs), self.dtype) + x = np.zeros((n, nrhs), self.dtype) + + ldb = b.shape[0] + ldx = x.shape[0] + + self.__checkVector( n, nrhs, x ) + self.__checkVector( n, nrhs, b ) + + pyspm_spmGenRHS( rhstype, nrhs, self.id_ptr, + x.ctypes.data_as(c_void_p), ldx, + b.ctypes.data_as(c_void_p), ldb ) + return x, b diff --git a/wrappers/python/spm_driver.py b/wrappers/python/spm_driver.py new file mode 100755 index 00000000..fa88e594 --- /dev/null +++ b/wrappers/python/spm_driver.py @@ -0,0 +1,42 @@ +#!/usr/bin/env python3 +""" + @file spm_driver.py + + SPM example to generate a sparse matrix from the spm drivers + + @copyright 2017-2018 Bordeaux INP, CNRS (LaBRI UMR 5800), Inria, + Univ. Bordeaux. All rights reserved. + + @version 6.0.0 + @author Pierre Ramet + @author Mathieu Faverge + @author Louis Poirel + @date 2017-05-04 + +""" +import spm +import numpy as np + +# Hack to make sure that the mkl is loaded +tmp = np.eye(2).dot(np.ones(2)) + +# Load a sparse matrix from the Laplacian driver +A = spm.spmatrix( None, driver=spm.driver.Laplacian, filename="10:10:10:2.:1." ) + +# Example from a RSA file +#A = spm( None, driver=driver.RSA, filename="$PASTIX_DIR/test/matrix/oilpan.rsa" ) + +A.printInfo() + +# Scale A for low-rank: A / ||A||_f +norm = A.norm() +A.scale( 1. / norm ) + +# Generate b and x0 vectors such that A * x0 = b +nrhs = 10 +x0, b = A.genRHS( spm.rhstype.RndX, nrhs ) + +# Check that A * x = b +x = x0.copy() +A.checkAxb( x0, b, x ) + diff --git a/wrappers/python/spm_scipy.py b/wrappers/python/spm_scipy.py new file mode 100755 index 00000000..6a574b6a --- /dev/null +++ b/wrappers/python/spm_scipy.py @@ -0,0 +1,36 @@ +#!/usr/bin/env python3 +""" + @file spm_scipy.py + + SPM example to geneate a sparse matrix from Scipy to SPM + + @copyright 2017-2018 Bordeaux INP, CNRS (LaBRI UMR 5800), Inria, + Univ. Bordeaux. All rights reserved. + + @version 6.0.0 + @author Pierre Ramet + @author Mathieu Faverge + @author Louis Poirel + @date 2017-05-04 + +""" +import spm +import scipy.sparse as sps +import numpy as np + +# Hack to make sure that the mkl is loaded +tmp = np.eye(2).dot(np.ones(2)) + +# Set the problem +n = 9 +A = sps.spdiags([np.ones(n)*i for i in [4, -1, -1, -1, -1]], + [0, 1, 3, -1, -3], n, n) +x0 = np.arange(n) +b = A.dot(x0) + +spmA = spm.spmatrix( A ) + +# Check that A * x = b +x = x0.copy() +spmA.checkAxb( x0, b, x ) + -- GitLab