Skip to content

Commit 8e41f36

Browse files
committed
updates, mostly formatting
1 parent a7e8350 commit 8e41f36

File tree

76 files changed

+1005
-859
lines changed

Some content is hidden

Large Commits have some content hidden by default. Use the searchbox below for content that may be hidden.

76 files changed

+1005
-859
lines changed

base/PyNucleus_base/CSR_LinearOperator_{SCALAR}.pxi

+20-16
Original file line numberDiff line numberDiff line change
@@ -11,9 +11,7 @@ cdef class {SCALAR_label}CSR_LinearOperator({SCALAR_label}LinearOperator):
1111
INDEX_t[::1] indptr,
1212
{SCALAR}_t[::1] data,
1313
int NoThreads=1):
14-
{SCALAR_label}LinearOperator.__init__(self,
15-
indptr.shape[0]-1,
16-
indptr.shape[0]-1)
14+
{SCALAR_label}LinearOperator.__init__(self, indptr.shape[0]-1, indptr.shape[0]-1)
1715
self.indices = indices
1816
self.indptr = indptr
1917
self.data = data
@@ -42,6 +40,18 @@ cdef class {SCALAR_label}CSR_LinearOperator({SCALAR_label}LinearOperator):
4240
spmv(self.indptr, self.indices, self.data, x, y, overwrite=False)
4341
return 0
4442

43+
cdef INDEX_t matvecTrans({SCALAR_label}CSR_LinearOperator self,
44+
{SCALAR}_t[::1] x,
45+
{SCALAR}_t[::1] y) except -1:
46+
spmv(self.indptr, self.indices, self.data, x, y, overwrite=True, trans=True)
47+
return 0
48+
49+
cdef INDEX_t matvecTrans_no_overwrite({SCALAR_label}CSR_LinearOperator self,
50+
{SCALAR}_t[::1] x,
51+
{SCALAR}_t[::1] y) except -1:
52+
spmv(self.indptr, self.indices, self.data, x, y, overwrite=False, trans=True)
53+
return 0
54+
4555
cdef INDEX_t matvec_multi({SCALAR_label}CSR_LinearOperator self,
4656
{SCALAR}_t[:, ::1] x,
4757
{SCALAR}_t[:, ::1] y) except -1:
@@ -220,8 +230,7 @@ cdef class {SCALAR_label}CSR_LinearOperator({SCALAR_label}LinearOperator):
220230
def get_diagonal(self):
221231
cdef:
222232
INDEX_t i, jj
223-
np.ndarray[{SCALAR}_t, ndim=1] diag_mem = np.zeros((self.num_rows),
224-
dtype={SCALAR})
233+
np.ndarray[{SCALAR}_t, ndim=1] diag_mem = np.zeros((self.num_rows), dtype={SCALAR})
225234
{SCALAR}_t[::1] d = diag_mem
226235

227236
for i in range(self.num_rows):
@@ -273,8 +282,8 @@ cdef class {SCALAR_label}CSR_LinearOperator({SCALAR_label}LinearOperator):
273282
@staticmethod
274283
def HDF5read(node):
275284
B = {SCALAR_label}CSR_LinearOperator(np.array(node['indices'], dtype=INDEX),
276-
np.array(node['indptr'], dtype=INDEX),
277-
np.array(node['data'], dtype={SCALAR}))
285+
np.array(node['indptr'], dtype=INDEX),
286+
np.array(node['data'], dtype={SCALAR}))
278287
B.num_rows = node.attrs['num_rows']
279288
B.num_columns = node.attrs['num_columns']
280289
assert B.indptr.shape[0]-1 == B.num_rows
@@ -312,7 +321,7 @@ cdef class {SCALAR_label}CSR_LinearOperator({SCALAR_label}LinearOperator):
312321
nnz = self.indptr[self.indptr.shape[0]-1]
313322
for i in range(self.indptr.shape[0]-1):
314323
s = self.indptr[i]
315-
if s == nnz:
324+
if s == nnz:
316325
continue
317326
p = self.indices[s]
318327
for q in self.indices[self.indptr[i]+1:self.indptr[i+1]]:
@@ -440,9 +449,7 @@ cdef class {SCALAR_label}CSR_LinearOperator({SCALAR_label}LinearOperator):
440449
return blockA
441450

442451

443-
cdef BOOL_t sort_indices{SCALAR_label}(INDEX_t[::1] indptr,
444-
INDEX_t[::1] indices,
445-
{SCALAR}_t[::1] data):
452+
cdef BOOL_t sort_indices{SCALAR_label}(INDEX_t[::1] indptr, INDEX_t[::1] indices, {SCALAR}_t[::1] data):
446453
cdef:
447454
INDEX_t n, i, jj, j, kk
448455
{SCALAR}_t d
@@ -479,10 +486,7 @@ cdef class {SCALAR_label}CSR_VectorLinearOperator({SCALAR_label}VectorLinearOper
479486
INDEX_t[::1] indices,
480487
INDEX_t[::1] indptr,
481488
{SCALAR}_t[:, ::1] data):
482-
{SCALAR_label}VectorLinearOperator.__init__(self,
483-
indptr.shape[0]-1,
484-
indptr.shape[0]-1,
485-
data.shape[1])
489+
{SCALAR_label}VectorLinearOperator.__init__(self, indptr.shape[0]-1, indptr.shape[0]-1, data.shape[1])
486490
self.indices = indices
487491
self.indptr = indptr
488492
self.data = data
@@ -659,7 +663,7 @@ cdef class {SCALAR_label}CSR_VectorLinearOperator({SCALAR_label}VectorLinearOper
659663
nnz = self.indptr[self.indptr.shape[0]-1]
660664
for i in range(self.indptr.shape[0]-1):
661665
s = self.indptr[i]
662-
if s == nnz:
666+
if s == nnz:
663667
continue
664668
p = self.indices[s]
665669
for q in self.indices[self.indptr[i]+1:self.indptr[i+1]]:

base/PyNucleus_base/DenseLinearOperator_{SCALAR}.pxi

+21-16
Original file line numberDiff line numberDiff line change
@@ -8,9 +8,7 @@
88
cdef class {SCALAR_label}Dense_LinearOperator({SCALAR_label}LinearOperator):
99
def __init__(self,
1010
{SCALAR}_t[:, ::1] data):
11-
{SCALAR_label}LinearOperator.__init__(self,
12-
data.shape[0],
13-
data.shape[1])
11+
{SCALAR_label}LinearOperator.__init__(self, data.shape[0], data.shape[1])
1412
self.data = data
1513

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

43+
cdef INDEX_t matvecTrans(self,
44+
{SCALAR}_t[::1] x,
45+
{SCALAR}_t[::1] y) except -1:
46+
gemvT(self.data, x, y)
47+
return 0
48+
49+
cdef INDEX_t matvecTrans_no_overwrite(self,
50+
{SCALAR}_t[::1] x,
51+
{SCALAR}_t[::1] y) except -1:
52+
gemvT(self.data, x, y, 1.)
53+
return 0
54+
4555
property diagonal:
4656
def __get__(self):
4757
cdef INDEX_t i
4858
diag = uninitialized((min(self.num_rows, self.num_columns)),
49-
dtype={SCALAR})
59+
dtype={SCALAR})
5060
for i in range(min(self.num_rows, self.num_columns)):
5161
diag[i] = self.data[i, i]
5262
return diag
@@ -105,13 +115,13 @@ cdef class {SCALAR_label}Dense_LinearOperator({SCALAR_label}LinearOperator):
105115
sizeInMB = self.getMemorySize() >> 20
106116
if sizeInMB > 100:
107117
return '<%dx%d %s, %d MB>' % (self.num_rows,
108-
self.num_columns,
109-
self.__class__.__name__,
110-
sizeInMB)
118+
self.num_columns,
119+
self.__class__.__name__,
120+
sizeInMB)
111121
else:
112122
return '<%dx%d %s>' % (self.num_rows,
113-
self.num_columns,
114-
self.__class__.__name__)
123+
self.num_columns,
124+
self.__class__.__name__)
115125

116126

117127
cdef class {SCALAR_label}Dense_SubBlock_LinearOperator({SCALAR_label}LinearOperator):
@@ -121,9 +131,7 @@ cdef class {SCALAR_label}Dense_SubBlock_LinearOperator({SCALAR_label}LinearOpera
121131
if mem is None:
122132
mem = np.zeros((I.shape[0], J.shape[0]), dtype={SCALAR})
123133
self.data = mem
124-
{SCALAR_label}LinearOperator.__init__(self,
125-
num_rows,
126-
num_columns)
134+
{SCALAR_label}LinearOperator.__init__(self, num_rows, num_columns)
127135
self.lookupI = {}
128136
self.lookupJ = {}
129137
for i in range(I.shape[0]):
@@ -161,10 +169,7 @@ cdef class {SCALAR_label}Dense_SubBlock_LinearOperator({SCALAR_label}LinearOpera
161169
cdef class {SCALAR_label}Dense_VectorLinearOperator({SCALAR_label}VectorLinearOperator):
162170
def __init__(self,
163171
{SCALAR}_t[:, :, ::1] data):
164-
{SCALAR_label}VectorLinearOperator.__init__(self,
165-
data.shape[0],
166-
data.shape[1],
167-
data.shape[2])
172+
{SCALAR_label}VectorLinearOperator.__init__(self, data.shape[0], data.shape[1], data.shape[2])
168173
self.data = data
169174

170175
cdef INDEX_t matvec(self,

base/PyNucleus_base/DiagonalLinearOperator_{SCALAR}.pxi

+2-5
Original file line numberDiff line numberDiff line change
@@ -7,9 +7,7 @@
77

88
cdef class {SCALAR_label}diagonalOperator({SCALAR_label}LinearOperator):
99
def __init__(self, {SCALAR}_t[::1] diagonal):
10-
{SCALAR_label}LinearOperator.__init__(self,
11-
diagonal.shape[0],
12-
diagonal.shape[0])
10+
{SCALAR_label}LinearOperator.__init__(self, diagonal.shape[0], diagonal.shape[0])
1311
self.data = diagonal
1412

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

8886
cdef class {SCALAR_label}invDiagonal({SCALAR_label}diagonalOperator):
8987
def __init__(self, {SCALAR_label}LinearOperator A):
90-
{SCALAR_label}diagonalOperator.__init__(self,
91-
1./np.array(A.diagonal))
88+
{SCALAR_label}diagonalOperator.__init__(self, 1./np.array(A.diagonal))

base/PyNucleus_base/LinearOperator_decl_{SCALAR}.pxi

+19
Original file line numberDiff line numberDiff line change
@@ -22,6 +22,12 @@ cdef class {SCALAR_label}LinearOperator:
2222
cdef INDEX_t matvec_multi(self,
2323
{SCALAR}_t[:, ::1] x,
2424
{SCALAR}_t[:, ::1] y) except -1
25+
cdef INDEX_t matvecTrans(self,
26+
{SCALAR}_t[::1] x,
27+
{SCALAR}_t[::1] y) except -1
28+
cdef INDEX_t matvecTrans_no_overwrite(self,
29+
{SCALAR}_t[::1] x,
30+
{SCALAR}_t[::1] y) except -1
2531
cdef void residual(self,
2632
{SCALAR}_t[::1] x,
2733
{SCALAR}_t[::1] rhs,
@@ -107,3 +113,16 @@ cdef class {SCALAR_label}VectorLinearOperator:
107113
cdef void addToEntry(self, INDEX_t I, INDEX_t J, {SCALAR}_t[::1] val)
108114
cdef void getEntry(self, INDEX_t I, INDEX_t J, {SCALAR}_t[::1] val)
109115
cdef void setEntry(self, INDEX_t I, INDEX_t J, {SCALAR}_t[::1] val)
116+
117+
118+
cdef class {SCALAR_label}Transpose_Linear_Operator({SCALAR_label}LinearOperator):
119+
cdef:
120+
public {SCALAR_label}LinearOperator A
121+
cdef INDEX_t matvec(self,
122+
{SCALAR}_t[::1] x,
123+
{SCALAR}_t[::1] y) except -1
124+
cdef void _residual(self,
125+
{SCALAR}_t[::1] x,
126+
{SCALAR}_t[::1] rhs,
127+
{SCALAR}_t[::1] result,
128+
BOOL_t simpleResidual=*)

base/PyNucleus_base/LinearOperator_{SCALAR}.pxi

+78-15
Original file line numberDiff line numberDiff line change
@@ -25,19 +25,32 @@ cdef class {SCALAR_label}LinearOperator:
2525
{SCALAR}_t[:, ::1] y) except -1:
2626
return -1
2727

28+
cdef INDEX_t matvecTrans(self,
29+
{SCALAR}_t[::1] x,
30+
{SCALAR}_t[::1] y) except -1:
31+
return -1
32+
33+
cdef INDEX_t matvecTrans_no_overwrite(self,
34+
{SCALAR}_t[::1] x,
35+
{SCALAR}_t[::1] y) except -1:
36+
return -1
37+
2838
def __call__(self,
2939
{SCALAR}_t[::1] x,
3040
{SCALAR}_t[::1] y,
31-
BOOL_t no_overwrite=False):
32-
if no_overwrite:
33-
self.matvec_no_overwrite(x, y)
41+
BOOL_t no_overwrite=False,
42+
BOOL_t trans=False):
43+
if not trans:
44+
if no_overwrite:
45+
self.matvec_no_overwrite(x, y)
46+
else:
47+
self.matvec(x, y)
3448
else:
35-
self.matvec(x, y)
49+
self.matvecTrans(x, y)
3650

3751
def dot(self, {SCALAR}_t[::1] x):
3852
cdef:
39-
np.ndarray[{SCALAR}_t, ndim=1] y = np.zeros(self.num_rows,
40-
dtype={SCALAR})
53+
np.ndarray[{SCALAR}_t, ndim=1] y = np.zeros(self.num_rows, dtype={SCALAR})
4154
self(x, y)
4255
return y
4356

@@ -274,6 +287,10 @@ cdef class {SCALAR_label}LinearOperator:
274287
def getMemorySize(self):
275288
return -1
276289

290+
@property
291+
def T(self):
292+
return {SCALAR_label}Transpose_Linear_Operator(self)
293+
277294

278295
cdef class {SCALAR_label}TimeStepperLinearOperator({SCALAR_label}LinearOperator):
279296
def __init__(self,
@@ -283,8 +300,7 @@ cdef class {SCALAR_label}TimeStepperLinearOperator({SCALAR_label}LinearOperator)
283300
{SCALAR}_t facM=1.0):
284301
assert M.num_columns == S.num_columns
285302
assert M.num_rows == S.num_rows
286-
super({SCALAR_label}TimeStepperLinearOperator, self).__init__(M.num_rows,
287-
M.num_columns)
303+
super({SCALAR_label}TimeStepperLinearOperator, self).__init__(M.num_rows, M.num_columns)
288304
self.M = M
289305
self.S = S
290306
self.facM = facM
@@ -306,8 +322,8 @@ cdef class {SCALAR_label}TimeStepperLinearOperator({SCALAR_label}LinearOperator)
306322
return 0
307323

308324
cdef INDEX_t matvec_no_overwrite(self,
309-
{SCALAR}_t[::1] x,
310-
{SCALAR}_t[::1] y) except -1:
325+
{SCALAR}_t[::1] x,
326+
{SCALAR}_t[::1] y) except -1:
311327
if self.facS == 1.0:
312328
self.S.matvec_no_overwrite(x, y)
313329
elif self.facS != 0.:
@@ -374,8 +390,7 @@ cdef class {SCALAR_label}Multiply_Linear_Operator({SCALAR_label}LinearOperator):
374390
def __init__(self,
375391
{SCALAR_label}LinearOperator A,
376392
{SCALAR}_t factor):
377-
super({SCALAR_label}Multiply_Linear_Operator, self).__init__(A.num_rows,
378-
A.num_columns)
393+
super({SCALAR_label}Multiply_Linear_Operator, self).__init__(A.num_rows, A.num_columns)
379394
self.A = A
380395
self.factor = factor
381396

@@ -420,7 +435,7 @@ cdef class {SCALAR_label}Multiply_Linear_Operator({SCALAR_label}LinearOperator):
420435
elif isinstance(x, {SCALAR_label}Multiply_Linear_Operator) and isinstance(self, ({SCALAR}, float)):
421436
return {SCALAR_label}Multiply_Linear_Operator(x.A, x.factor*self)
422437
elif isinstance(x, COMPLEX):
423-
return ComplexMultiply_Linear_Operator(wrapRealToComplex(self.A), self.factor*x)
438+
return ComplexMultiply_Linear_Operator(wrapRealToComplex(self.A), self.factor*x)
424439
else:
425440
return super({SCALAR_label}Multiply_Linear_Operator, self).__mul__(x)
426441

@@ -442,8 +457,7 @@ cdef class {SCALAR_label}Product_Linear_Operator({SCALAR_label}LinearOperator):
442457
{SCALAR_label}LinearOperator B,
443458
{SCALAR}_t[::1] temporaryMemory=None):
444459
assert A.num_columns == B.num_rows, '{} and {} are not compatible'.format(A.num_columns, B.num_rows)
445-
super({SCALAR_label}Product_Linear_Operator, self).__init__(A.num_rows,
446-
B.num_columns)
460+
super({SCALAR_label}Product_Linear_Operator, self).__init__(A.num_rows, B.num_columns)
447461
self.A = A
448462
self.B = B
449463
if temporaryMemory is not None:
@@ -626,3 +640,52 @@ cdef class {SCALAR_label}VectorLinearOperator:
626640

627641
def getEntry_py(self, INDEX_t I, INDEX_t J, {SCALAR}_t[::1] val):
628642
return self.getEntry(I, J, val)
643+
644+
645+
cdef class {SCALAR_label}Transpose_Linear_Operator({SCALAR_label}LinearOperator):
646+
def __init__(self,
647+
{SCALAR_label}LinearOperator A):
648+
super({SCALAR_label}Transpose_Linear_Operator, self).__init__(A.num_columns, A.num_rows)
649+
self.A = A
650+
651+
cdef INDEX_t matvec(self,
652+
{SCALAR}_t[::1] x,
653+
{SCALAR}_t[::1] y) except -1:
654+
self.A.matvecTrans(x, y)
655+
return 0
656+
657+
cdef INDEX_t matvec_no_overwrite(self,
658+
{SCALAR}_t[::1] x,
659+
{SCALAR}_t[::1] y) except -1:
660+
self.A.matvecTrans_no_overwrite(x, y)
661+
return 0
662+
663+
def isSparse(self):
664+
return self.A.isSparse()
665+
666+
def to_csr(self):
667+
return self.A.to_csr().T
668+
669+
def to_csr_linear_operator(self):
670+
if isinstance(self.A, {SCALAR_label}Dense_LinearOperator):
671+
return {SCALAR_label}Dense_LinearOperator(self.A.transpose())
672+
else:
673+
B = self.A.transpose()
674+
Bcsr = {SCALAR_label}CSR_LinearOperator(B.indices, B.indptr, B.data)
675+
Bcsr.num_rows = B.shape[0]
676+
Bcsr.num_columns = B.shape[1]
677+
return Bcsr
678+
679+
def toarray(self):
680+
return self.A.transpose().toarray()
681+
682+
def get_diagonal(self):
683+
return np.array(self.A.diagonal, copy=False)
684+
685+
diagonal = property(fget=get_diagonal)
686+
687+
def __repr__(self):
688+
return 'transpose({})'.format(self.A)
689+
690+
def getMemorySize(self):
691+
return self.A.getMemorySize()

0 commit comments

Comments
 (0)