Skip to content

updates, mostly formatting #28

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
merged 1 commit into from
Nov 25, 2023
Merged

updates, mostly formatting #28

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
36 changes: 20 additions & 16 deletions base/PyNucleus_base/CSR_LinearOperator_{SCALAR}.pxi
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -11,9 +11,7 @@ cdef class {SCALAR_label}CSR_LinearOperator({SCALAR_label}LinearOperator):
INDEX_t[::1] indptr,
{SCALAR}_t[::1] data,
int NoThreads=1):
{SCALAR_label}LinearOperator.__init__(self,
indptr.shape[0]-1,
indptr.shape[0]-1)
{SCALAR_label}LinearOperator.__init__(self, indptr.shape[0]-1, indptr.shape[0]-1)
self.indices = indices
self.indptr = indptr
self.data = data
Expand Down Expand Up @@ -42,6 +40,18 @@ cdef class {SCALAR_label}CSR_LinearOperator({SCALAR_label}LinearOperator):
spmv(self.indptr, self.indices, self.data, x, y, overwrite=False)
return 0

cdef INDEX_t matvecTrans({SCALAR_label}CSR_LinearOperator self,
{SCALAR}_t[::1] x,
{SCALAR}_t[::1] y) except -1:
spmv(self.indptr, self.indices, self.data, x, y, overwrite=True, trans=True)
return 0

cdef INDEX_t matvecTrans_no_overwrite({SCALAR_label}CSR_LinearOperator self,
{SCALAR}_t[::1] x,
{SCALAR}_t[::1] y) except -1:
spmv(self.indptr, self.indices, self.data, x, y, overwrite=False, trans=True)
return 0

cdef INDEX_t matvec_multi({SCALAR_label}CSR_LinearOperator self,
{SCALAR}_t[:, ::1] x,
{SCALAR}_t[:, ::1] y) except -1:
Expand Down Expand Up @@ -220,8 +230,7 @@ cdef class {SCALAR_label}CSR_LinearOperator({SCALAR_label}LinearOperator):
def get_diagonal(self):
cdef:
INDEX_t i, jj
np.ndarray[{SCALAR}_t, ndim=1] diag_mem = np.zeros((self.num_rows),
dtype={SCALAR})
np.ndarray[{SCALAR}_t, ndim=1] diag_mem = np.zeros((self.num_rows), dtype={SCALAR})
{SCALAR}_t[::1] d = diag_mem

for i in range(self.num_rows):
Expand Down Expand Up @@ -273,8 +282,8 @@ cdef class {SCALAR_label}CSR_LinearOperator({SCALAR_label}LinearOperator):
@staticmethod
def HDF5read(node):
B = {SCALAR_label}CSR_LinearOperator(np.array(node['indices'], dtype=INDEX),
np.array(node['indptr'], dtype=INDEX),
np.array(node['data'], dtype={SCALAR}))
np.array(node['indptr'], dtype=INDEX),
np.array(node['data'], dtype={SCALAR}))
B.num_rows = node.attrs['num_rows']
B.num_columns = node.attrs['num_columns']
assert B.indptr.shape[0]-1 == B.num_rows
Expand Down Expand Up @@ -312,7 +321,7 @@ cdef class {SCALAR_label}CSR_LinearOperator({SCALAR_label}LinearOperator):
nnz = self.indptr[self.indptr.shape[0]-1]
for i in range(self.indptr.shape[0]-1):
s = self.indptr[i]
if s == nnz:
if s == nnz:
continue
p = self.indices[s]
for q in self.indices[self.indptr[i]+1:self.indptr[i+1]]:
Expand Down Expand Up @@ -440,9 +449,7 @@ cdef class {SCALAR_label}CSR_LinearOperator({SCALAR_label}LinearOperator):
return blockA


cdef BOOL_t sort_indices{SCALAR_label}(INDEX_t[::1] indptr,
INDEX_t[::1] indices,
{SCALAR}_t[::1] data):
cdef BOOL_t sort_indices{SCALAR_label}(INDEX_t[::1] indptr, INDEX_t[::1] indices, {SCALAR}_t[::1] data):
cdef:
INDEX_t n, i, jj, j, kk
{SCALAR}_t d
Expand Down Expand Up @@ -479,10 +486,7 @@ cdef class {SCALAR_label}CSR_VectorLinearOperator({SCALAR_label}VectorLinearOper
INDEX_t[::1] indices,
INDEX_t[::1] indptr,
{SCALAR}_t[:, ::1] data):
{SCALAR_label}VectorLinearOperator.__init__(self,
indptr.shape[0]-1,
indptr.shape[0]-1,
data.shape[1])
{SCALAR_label}VectorLinearOperator.__init__(self, indptr.shape[0]-1, indptr.shape[0]-1, data.shape[1])
self.indices = indices
self.indptr = indptr
self.data = data
Expand Down Expand Up @@ -659,7 +663,7 @@ cdef class {SCALAR_label}CSR_VectorLinearOperator({SCALAR_label}VectorLinearOper
nnz = self.indptr[self.indptr.shape[0]-1]
for i in range(self.indptr.shape[0]-1):
s = self.indptr[i]
if s == nnz:
if s == nnz:
continue
p = self.indices[s]
for q in self.indices[self.indptr[i]+1:self.indptr[i+1]]:
Expand Down
37 changes: 21 additions & 16 deletions base/PyNucleus_base/DenseLinearOperator_{SCALAR}.pxi
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -8,9 +8,7 @@
cdef class {SCALAR_label}Dense_LinearOperator({SCALAR_label}LinearOperator):
def __init__(self,
{SCALAR}_t[:, ::1] data):
{SCALAR_label}LinearOperator.__init__(self,
data.shape[0],
data.shape[1])
{SCALAR_label}LinearOperator.__init__(self, data.shape[0], data.shape[1])
self.data = data

cdef INDEX_t matvec(self,
Expand Down Expand Up @@ -42,11 +40,23 @@ cdef class {SCALAR_label}Dense_LinearOperator({SCALAR_label}LinearOperator):
y[i, k] = temp[k]
return 0

cdef INDEX_t matvecTrans(self,
{SCALAR}_t[::1] x,
{SCALAR}_t[::1] y) except -1:
gemvT(self.data, x, y)
return 0

cdef INDEX_t matvecTrans_no_overwrite(self,
{SCALAR}_t[::1] x,
{SCALAR}_t[::1] y) except -1:
gemvT(self.data, x, y, 1.)
return 0

property diagonal:
def __get__(self):
cdef INDEX_t i
diag = uninitialized((min(self.num_rows, self.num_columns)),
dtype={SCALAR})
dtype={SCALAR})
for i in range(min(self.num_rows, self.num_columns)):
diag[i] = self.data[i, i]
return diag
Expand Down Expand Up @@ -105,13 +115,13 @@ cdef class {SCALAR_label}Dense_LinearOperator({SCALAR_label}LinearOperator):
sizeInMB = self.getMemorySize() >> 20
if sizeInMB > 100:
return '<%dx%d %s, %d MB>' % (self.num_rows,
self.num_columns,
self.__class__.__name__,
sizeInMB)
self.num_columns,
self.__class__.__name__,
sizeInMB)
else:
return '<%dx%d %s>' % (self.num_rows,
self.num_columns,
self.__class__.__name__)
self.num_columns,
self.__class__.__name__)


cdef class {SCALAR_label}Dense_SubBlock_LinearOperator({SCALAR_label}LinearOperator):
Expand All @@ -121,9 +131,7 @@ cdef class {SCALAR_label}Dense_SubBlock_LinearOperator({SCALAR_label}LinearOpera
if mem is None:
mem = np.zeros((I.shape[0], J.shape[0]), dtype={SCALAR})
self.data = mem
{SCALAR_label}LinearOperator.__init__(self,
num_rows,
num_columns)
{SCALAR_label}LinearOperator.__init__(self, num_rows, num_columns)
self.lookupI = {}
self.lookupJ = {}
for i in range(I.shape[0]):
Expand Down Expand Up @@ -161,10 +169,7 @@ cdef class {SCALAR_label}Dense_SubBlock_LinearOperator({SCALAR_label}LinearOpera
cdef class {SCALAR_label}Dense_VectorLinearOperator({SCALAR_label}VectorLinearOperator):
def __init__(self,
{SCALAR}_t[:, :, ::1] data):
{SCALAR_label}VectorLinearOperator.__init__(self,
data.shape[0],
data.shape[1],
data.shape[2])
{SCALAR_label}VectorLinearOperator.__init__(self, data.shape[0], data.shape[1], data.shape[2])
self.data = data

cdef INDEX_t matvec(self,
Expand Down
7 changes: 2 additions & 5 deletions base/PyNucleus_base/DiagonalLinearOperator_{SCALAR}.pxi
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -7,9 +7,7 @@

cdef class {SCALAR_label}diagonalOperator({SCALAR_label}LinearOperator):
def __init__(self, {SCALAR}_t[::1] diagonal):
{SCALAR_label}LinearOperator.__init__(self,
diagonal.shape[0],
diagonal.shape[0])
{SCALAR_label}LinearOperator.__init__(self, diagonal.shape[0], diagonal.shape[0])
self.data = diagonal

cdef INDEX_t matvec(self,
Expand Down Expand Up @@ -87,5 +85,4 @@ cdef class {SCALAR_label}diagonalOperator({SCALAR_label}LinearOperator):

cdef class {SCALAR_label}invDiagonal({SCALAR_label}diagonalOperator):
def __init__(self, {SCALAR_label}LinearOperator A):
{SCALAR_label}diagonalOperator.__init__(self,
1./np.array(A.diagonal))
{SCALAR_label}diagonalOperator.__init__(self, 1./np.array(A.diagonal))
19 changes: 19 additions & 0 deletions base/PyNucleus_base/LinearOperator_decl_{SCALAR}.pxi
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -22,6 +22,12 @@ cdef class {SCALAR_label}LinearOperator:
cdef INDEX_t matvec_multi(self,
{SCALAR}_t[:, ::1] x,
{SCALAR}_t[:, ::1] y) except -1
cdef INDEX_t matvecTrans(self,
{SCALAR}_t[::1] x,
{SCALAR}_t[::1] y) except -1
cdef INDEX_t matvecTrans_no_overwrite(self,
{SCALAR}_t[::1] x,
{SCALAR}_t[::1] y) except -1
cdef void residual(self,
{SCALAR}_t[::1] x,
{SCALAR}_t[::1] rhs,
Expand Down Expand Up @@ -107,3 +113,16 @@ cdef class {SCALAR_label}VectorLinearOperator:
cdef void addToEntry(self, INDEX_t I, INDEX_t J, {SCALAR}_t[::1] val)
cdef void getEntry(self, INDEX_t I, INDEX_t J, {SCALAR}_t[::1] val)
cdef void setEntry(self, INDEX_t I, INDEX_t J, {SCALAR}_t[::1] val)


cdef class {SCALAR_label}Transpose_Linear_Operator({SCALAR_label}LinearOperator):
cdef:
public {SCALAR_label}LinearOperator A
cdef INDEX_t matvec(self,
{SCALAR}_t[::1] x,
{SCALAR}_t[::1] y) except -1
cdef void _residual(self,
{SCALAR}_t[::1] x,
{SCALAR}_t[::1] rhs,
{SCALAR}_t[::1] result,
BOOL_t simpleResidual=*)
93 changes: 78 additions & 15 deletions base/PyNucleus_base/LinearOperator_{SCALAR}.pxi
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -25,19 +25,32 @@ cdef class {SCALAR_label}LinearOperator:
{SCALAR}_t[:, ::1] y) except -1:
return -1

cdef INDEX_t matvecTrans(self,
{SCALAR}_t[::1] x,
{SCALAR}_t[::1] y) except -1:
return -1

cdef INDEX_t matvecTrans_no_overwrite(self,
{SCALAR}_t[::1] x,
{SCALAR}_t[::1] y) except -1:
return -1

def __call__(self,
{SCALAR}_t[::1] x,
{SCALAR}_t[::1] y,
BOOL_t no_overwrite=False):
if no_overwrite:
self.matvec_no_overwrite(x, y)
BOOL_t no_overwrite=False,
BOOL_t trans=False):
if not trans:
if no_overwrite:
self.matvec_no_overwrite(x, y)
else:
self.matvec(x, y)
else:
self.matvec(x, y)
self.matvecTrans(x, y)

def dot(self, {SCALAR}_t[::1] x):
cdef:
np.ndarray[{SCALAR}_t, ndim=1] y = np.zeros(self.num_rows,
dtype={SCALAR})
np.ndarray[{SCALAR}_t, ndim=1] y = np.zeros(self.num_rows, dtype={SCALAR})
self(x, y)
return y

Expand Down Expand Up @@ -274,6 +287,10 @@ cdef class {SCALAR_label}LinearOperator:
def getMemorySize(self):
return -1

@property
def T(self):
return {SCALAR_label}Transpose_Linear_Operator(self)


cdef class {SCALAR_label}TimeStepperLinearOperator({SCALAR_label}LinearOperator):
def __init__(self,
Expand All @@ -283,8 +300,7 @@ cdef class {SCALAR_label}TimeStepperLinearOperator({SCALAR_label}LinearOperator)
{SCALAR}_t facM=1.0):
assert M.num_columns == S.num_columns
assert M.num_rows == S.num_rows
super({SCALAR_label}TimeStepperLinearOperator, self).__init__(M.num_rows,
M.num_columns)
super({SCALAR_label}TimeStepperLinearOperator, self).__init__(M.num_rows, M.num_columns)
self.M = M
self.S = S
self.facM = facM
Expand All @@ -306,8 +322,8 @@ cdef class {SCALAR_label}TimeStepperLinearOperator({SCALAR_label}LinearOperator)
return 0

cdef INDEX_t matvec_no_overwrite(self,
{SCALAR}_t[::1] x,
{SCALAR}_t[::1] y) except -1:
{SCALAR}_t[::1] x,
{SCALAR}_t[::1] y) except -1:
if self.facS == 1.0:
self.S.matvec_no_overwrite(x, y)
elif self.facS != 0.:
Expand Down Expand Up @@ -374,8 +390,7 @@ cdef class {SCALAR_label}Multiply_Linear_Operator({SCALAR_label}LinearOperator):
def __init__(self,
{SCALAR_label}LinearOperator A,
{SCALAR}_t factor):
super({SCALAR_label}Multiply_Linear_Operator, self).__init__(A.num_rows,
A.num_columns)
super({SCALAR_label}Multiply_Linear_Operator, self).__init__(A.num_rows, A.num_columns)
self.A = A
self.factor = factor

Expand Down Expand Up @@ -420,7 +435,7 @@ cdef class {SCALAR_label}Multiply_Linear_Operator({SCALAR_label}LinearOperator):
elif isinstance(x, {SCALAR_label}Multiply_Linear_Operator) and isinstance(self, ({SCALAR}, float)):
return {SCALAR_label}Multiply_Linear_Operator(x.A, x.factor*self)
elif isinstance(x, COMPLEX):
return ComplexMultiply_Linear_Operator(wrapRealToComplex(self.A), self.factor*x)
return ComplexMultiply_Linear_Operator(wrapRealToComplex(self.A), self.factor*x)
else:
return super({SCALAR_label}Multiply_Linear_Operator, self).__mul__(x)

Expand All @@ -442,8 +457,7 @@ cdef class {SCALAR_label}Product_Linear_Operator({SCALAR_label}LinearOperator):
{SCALAR_label}LinearOperator B,
{SCALAR}_t[::1] temporaryMemory=None):
assert A.num_columns == B.num_rows, '{} and {} are not compatible'.format(A.num_columns, B.num_rows)
super({SCALAR_label}Product_Linear_Operator, self).__init__(A.num_rows,
B.num_columns)
super({SCALAR_label}Product_Linear_Operator, self).__init__(A.num_rows, B.num_columns)
self.A = A
self.B = B
if temporaryMemory is not None:
Expand Down Expand Up @@ -626,3 +640,52 @@ cdef class {SCALAR_label}VectorLinearOperator:

def getEntry_py(self, INDEX_t I, INDEX_t J, {SCALAR}_t[::1] val):
return self.getEntry(I, J, val)


cdef class {SCALAR_label}Transpose_Linear_Operator({SCALAR_label}LinearOperator):
def __init__(self,
{SCALAR_label}LinearOperator A):
super({SCALAR_label}Transpose_Linear_Operator, self).__init__(A.num_columns, A.num_rows)
self.A = A

cdef INDEX_t matvec(self,
{SCALAR}_t[::1] x,
{SCALAR}_t[::1] y) except -1:
self.A.matvecTrans(x, y)
return 0

cdef INDEX_t matvec_no_overwrite(self,
{SCALAR}_t[::1] x,
{SCALAR}_t[::1] y) except -1:
self.A.matvecTrans_no_overwrite(x, y)
return 0

def isSparse(self):
return self.A.isSparse()

def to_csr(self):
return self.A.to_csr().T

def to_csr_linear_operator(self):
if isinstance(self.A, {SCALAR_label}Dense_LinearOperator):
return {SCALAR_label}Dense_LinearOperator(self.A.transpose())
else:
B = self.A.transpose()
Bcsr = {SCALAR_label}CSR_LinearOperator(B.indices, B.indptr, B.data)
Bcsr.num_rows = B.shape[0]
Bcsr.num_columns = B.shape[1]
return Bcsr

def toarray(self):
return self.A.transpose().toarray()

def get_diagonal(self):
return np.array(self.A.diagonal, copy=False)

diagonal = property(fget=get_diagonal)

def __repr__(self):
return 'transpose({})'.format(self.A)

def getMemorySize(self):
return self.A.getMemorySize()
Loading