scipy used for direct balancing method
AntonioWR opened this issue · 1 comments
Describe the bug
In librom.py under balreal_direct_py, Scipy Schur() decomposition, solve_discrete_lyapunov() and svd() do not accept csr_matrix type for state matrix A.
My temporal way of dealing with it is to transfer it to a dense matrix by the todense() function shown below.
if Schur:
# decompose A
# A is a sparse matrix in csr_matrix format, can not be directly passed into scipy _solver.py so is transformed here first.
A = A.todense()
Atri, U = scalg.schur(A)
# solve Lyapunov
BBtri = np.dot(U.T, np.dot(B, np.dot(B.T, U)))
CCtri = np.dot(U.T, np.dot(C.T, np.dot(C, U)))
Wctri = sollyap(Atri, BBtri)
Wotri = sollyap(Atri.T, CCtri)
# reconstruct Wo,Wc
Wc = np.dot(U, np.dot(Wctri, U.T))
Wo = np.dot(U, np.dot(Wotri, U.T))
else:
# A is a sparse matrix in csr_matrix format, can not be directly passed into scipy _solver.py so is transformed here first.
A = A.todense()
Wc = sollyap(A, np.dot(B, B.T))
Wo = sollyap(A.T, np.dot(C.T, C))
To Reproduce
I am using the Goland wing example, but I change the ROM method to direct balancing method.
If model related, please try and reproduce with one of the provided models (found in /cases/
). Else, please provided the input files (a .fem.h5
, aero.h5
and case.sharpy
).
Expected behavior
A clear and concise description of what you expected to happen.
Screenshots
If applicable, add screenshots to help explain your problem.
System Info (please complete the following information):
- OS: Ubuntu 20.04
- SHARPy Version 1.3
Additional context
Add any other context about the problem here.
Hi Antonio,
Many thanks for raising this, and apologies for not getting back to you earlier.
I think your approach works well, because I haven't found a scipy
alternative to solve the Lyapunov equation for sparse matrices.
If you want you can open a pull request to fix this little issue with your solution. We would greatly appreciate that! (maybe you can add a warning message if the matrices are sparse saying that they need to be dense and that you are transforming them, or something along those lines)
Thanks,
Norberto