Skip to content

Commit

Permalink
chore(linalg3): settle deprecation warnings
Browse files Browse the repository at this point in the history
  • Loading branch information
nauaneed committed Sep 22, 2024
1 parent d349325 commit 6e5be66
Show file tree
Hide file tree
Showing 2 changed files with 26 additions and 26 deletions.
14 changes: 7 additions & 7 deletions pysph/base/linalg3.pxd
Original file line number Diff line number Diff line change
Expand Up @@ -3,11 +3,11 @@

"""Routines for eigen decomposition of symmetric 3x3 matrices.
"""
cdef double det(double a[3][3]) noexcept nogil
cdef void get_eigenvalues(double a[3][3], double *result) noexcept nogil
cdef void eigen_decomposition(double A[3][3], double V[3][3], double *d) noexcept nogil
cdef void transform(double A[3][3], double P[3][3], double res[3][3]) noexcept nogil
cdef void transform_diag(double *A, double P[3][3], double res[3][3]) noexcept nogil
cdef void transform_diag_inv(double *A, double P[3][3], double res[3][3]) nogil
cdef double det(double [3][3]a) noexcept nogil
cdef void get_eigenvalues(double [3][3]a, double *result) noexcept nogil
cdef void eigen_decomposition(double [3][3]A, double [3][3]V, double *d) noexcept nogil
cdef void transform(double [3][3]A, double [3][3]P, double [3][3]res) noexcept nogil
cdef void transform_diag(double *A, double [3][3]P, double [3][3]res) noexcept nogil
cdef void transform_diag_inv(double *A, double [3][3]P, double [3][3]res) nogil

cdef void get_eigenvalvec(double A[3][3], double *R, double *e)
cdef void get_eigenvalvec(double [3][3]A, double *R, double *e)
38 changes: 19 additions & 19 deletions pysph/base/linalg3.pyx
Original file line number Diff line number Diff line change
Expand Up @@ -30,7 +30,7 @@ cdef inline double hypot2(double x, double y) noexcept nogil:
return sqrt(x*x+y*y)


cdef double det(double a[3][3]) noexcept nogil:
cdef double det(double [3][3]a) noexcept nogil:
'''Determinant of symmetrix matrix
'''
return (a[0][0]*a[1][1]*a[2][2] + 2*a[1][2]*a[0][2]*a[0][1] -
Expand All @@ -46,12 +46,12 @@ cpdef double py_det(double[:,:] m):
# d:11,22,33; s:23,13,12
# d:00,11,22; s:12,02,01

cdef void get_eigenvalues(double a[3][3], double *result) noexcept nogil:
cdef void get_eigenvalues(double [3][3]a, double *result) noexcept nogil:
'''Compute the eigenvalues of symmetric matrix a and return in
result array.
'''
cdef double m = (a[0][0]+ a[1][1]+a[2][2])/3.0
cdef double K[3][3]
cdef double [3][3]K
memcpy(&K[0][0], &a[0][0], sizeof(double)*9)
K[0][0], K[1][1], K[2][2] = a[0][0] - m, a[1][1] - m, a[2][2] - m
cdef double q = det(K)*0.5
Expand Down Expand Up @@ -89,7 +89,7 @@ cpdef py_get_eigenvalues(double[:,:] m):


##############################################################################
cdef void get_eigenvector_np(double A[n][n], double r, double *res):
cdef void get_eigenvector_np(double [n][n]A, double r, double *res):
''' eigenvector of symmetric matrix for given eigenvalue `r` using numpy
'''
cdef numpy.ndarray[ndim=2, dtype=numpy.float64_t] mat=numpy.empty((3,3)), evec
Expand All @@ -109,7 +109,7 @@ cdef void get_eigenvector_np(double A[n][n], double r, double *res):
for i in range(3):
res[i] = evec[idx,i]

cdef void get_eigenvector(double A[n][n], double r, double *res):
cdef void get_eigenvector(double [n][n]A, double r, double *res):
''' get eigenvector of symmetric 3x3 matrix for given eigenvalue `r`
uses a fast method to get eigenvectors with a fallback to using numpy
Expand All @@ -134,22 +134,22 @@ cpdef py_get_eigenvector(double[:,:] A, double r):
get_eigenvector(<double(*)[n]>&A[0,0], r, &_d[0])
return d

cdef void get_eigenvec_from_val(double A[n][n], double *R, double *e):
cdef void get_eigenvec_from_val(double [n][n]A, double *R, double *e):
cdef int i, j
cdef double res[3]
cdef double [3]res
for i in range(3):
get_eigenvector(A, e[i], &res[0])
for j in range(3):
R[j*3+i] = res[j]


cdef bint _nearly_diagonal(double A[n][n]) noexcept nogil:
cdef bint _nearly_diagonal(double [n][n]A) noexcept nogil:
return (
(SQR(A[0][0]) + SQR(A[1][1]) + SQR(A[2][2])) >
1e8*(SQR(A[0][1]) + SQR(A[0][2]) + SQR(A[1][2]))
)

cdef void get_eigenvalvec(double A[n][n], double *R, double *e):
cdef void get_eigenvalvec(double [n][n]A, double *R, double *e):
'''Get the eigenvalues and eigenvectors of symmetric 3x3 matrix.
A is the input 3x3 matrix.
Expand Down Expand Up @@ -194,7 +194,7 @@ def py_get_eigenvalvec(double[:,:] A):

##############################################################################

cdef void transform(double A[3][3], double P[3][3], double res[3][3]) noexcept nogil:
cdef void transform(double [3][3]A, double [3][3]P, double [3][3]res) noexcept nogil:
'''Compute the transformation P.T*A*P and add it into res.
'''
cdef int i, j, k, l
Expand All @@ -204,8 +204,8 @@ cdef void transform(double A[3][3], double P[3][3], double res[3][3]) noexcept n
for l in range(3):
res[i][j] += P[k][i]*A[k][l]*P[l][j] # P.T*A*P

cdef void transform_diag(double *A, double P[3][3],
double res[3][3]) noexcept nogil:
cdef void transform_diag(double *A, double [3][3]P,
double [3][3]res) noexcept nogil:
'''Compute the transformation P.T*A*P and add it into res.
A is diagonal and contains the diagonal entries alone.
Expand All @@ -216,8 +216,8 @@ cdef void transform_diag(double *A, double P[3][3],
for k in range(3):
res[i][j] += P[k][i]*A[k]*P[k][j] # P.T*A*P

cdef void transform_diag_inv(double *A, double P[3][3],
double res[3][3]) noexcept nogil:
cdef void transform_diag_inv(double *A, double [3][3]P,
double [3][3]res) noexcept nogil:
'''Compute the transformation P*A*P.T and set it into res.
A is diagonal and contains just the diagonal entries.
'''
Expand Down Expand Up @@ -257,7 +257,7 @@ def py_transform_diag_inv(double[:] A, double[:,:] P):
return res


cdef double * tred2(double V[n][n], double *d, double *e) noexcept nogil:
cdef double * tred2(double [n][n]V, double *d, double *e) noexcept nogil:
'''Symmetric Householder reduction to tridiagonal form
This is derived from the Algol procedures tred2 by
Expand Down Expand Up @@ -376,7 +376,7 @@ cdef double * tred2(double V[n][n], double *d, double *e) noexcept nogil:
return d


cdef void tql2(double V[n][n], double *d, double *e) noexcept nogil:
cdef void tql2(double [n][n]V, double *d, double *e) noexcept nogil:
'''Symmetric tridiagonal QL algo for eigendecomposition
This is derived from the Algol procedures tql2, by
Expand Down Expand Up @@ -492,18 +492,18 @@ cdef void tql2(double V[n][n], double *d, double *e) noexcept nogil:
V[j][k] = p


cdef void zero_matrix_case(double V[n][n], double *d) noexcept nogil:
cdef void zero_matrix_case(double [n][n]V, double *d) noexcept nogil:
cdef int i, j
for i in range(3):
d[i] = 0.0
for j in range(3):
V[i][j] = (i==j)

cdef void eigen_decomposition(double A[n][n], double V[n][n], double *d) noexcept nogil:
cdef void eigen_decomposition(double [n][n]A, double [n][n]V, double *d) noexcept nogil:
'''Get eigenvalues and eigenvectors of matrix A.
V is output eigenvectors and d are the eigenvalues.
'''
cdef double e[n]
cdef double [n]e
cdef int i, j
# Scale the matrix, as if the matrix is tiny, floating point errors
# creep up leading to zero division errors in tql2. This is
Expand Down

0 comments on commit 6e5be66

Please sign in to comment.