scipy used for direct balancing method

[AntonioWR]

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.

[ngoiz]

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