Linear buckling analysis (structsolve.linear_buckling)#

structsolve.lb(K, KG, tol=0, sparse_solver=True, silent=False, num_eigvalues=25, num_eigvalues_print=5, skip_null_cols=False, max_dense_size=2000, check_rtol=0.001, symmetric=None)#

Linear buckling analysis

Calculates the eigenvalues \(\lambda\) and eigenvectors \(\{u\}\) of the eigenvalue problem:

\[([K] + \lambda [K_G])\{u\} = \{0\}\]

where \(\lambda\) is the load multiplier of the load state that generated \([K_G]\). It can also be used for more general eigenvalue analyses, e.g. if K is the tangent stiffness matrix of a given load state. The null rows and columns of K are removed from both matrices before solving the eigenvalue problem, unless skip_null_cols=True.

Parameters:
Ksparse_matrix

Stiffness matrix. Should include all constant terms of the initial stress stiffness matrix, aerodynamic matrix and so forth when applicable.

KGsparse_matrix

Initial stress stiffness matrix that multiplies the load multiplier \(\lambda\) of the eigenvalue problem.

tolfloat, optional

A float tolerance passed to scipy.sparse.linalg.eigsh() and scipy.sparse.linalg.eigs().

sparse_solverbool, optional

With True, the default, the sparse solvers always run first, whatever the size of the problem, and a dense solver is only used as a last resort, see the Notes. num_eigvalues load multipliers are returned. With False, the full problem is explicitly solved with a dense solver, without trying the sparse solvers: scipy.linalg.eigh(), which requires a positive definite K, or scipy.linalg.eig() for unsymmetric matrices, returning all the load multipliers. It is a reference for validations and for the callers that need all of them.

silentbool, optional

A boolean to tell whether the log messages should be printed.

num_eigvaluesint, optional

Number of load multipliers calculated with the sparse solver, and number of returned eigenvectors.

num_eigvalues_printint, optional

Number of eigenvalues to print.

skip_null_colsbool, optional

If True, skip the removal of null columns from the matrices. Use only when K is known to be non-singular.

max_dense_sizeint, optional

Maximum size of the dense problem solved as a last resort, when all the sparse solvers failed, i.e. the number of non-null rows of KG, or of non-null rows or columns for unsymmetric matrices, see the Notes. Larger problems raise a RuntimeError instead, and 0 disables the dense fallback. Until version 0.5.3, it was the size below which the dense solver was tried first.

check_rtolfloat, optional

Maximum relative residual of the eigenpairs of the sparse solver, see the Notes.

symmetricbool or None, optional

Whether K and KG are symmetric. None checks it with is_symmetric(), True skips the check and uses the solvers of symmetric matrices, False forces the solvers of unsymmetric matrices, required e.g. when KG includes the load stiffness of follower loads, see the Notes.

Returns:
eigvalsndarray

The load multipliers \(\lambda\), calculated as -1/eigval from the eigenvalues eigval of KG u = eigval K u. The positive load multipliers come first in increasing order, followed by the negative ones, for both solvers. For unsymmetric matrices the array is complex when some eigenvalues are complex, which come last.

eigvecsndarray

The \(i^{th}\) eigenvector is eigvecs[:, i], with the size of the original matrices. Only num_eigvalues eigenvectors are returned.

Raises:
RuntimeError

When no eigenvalue solver passed the verification, with the list of the attempted solvers and their failures.

Warns:
DenseFallbackWarning

When the sparse solvers failed and the dense solver was used.

RuntimeWarning

For unsymmetric matrices, when the dense fallback found no positive real load multiplier.

Notes

Solver policy. With sparse_solver=True the sparse solvers run first for every problem size, in the order given below, and each result must pass the verification. Only when all of them failed (verification, ARPACK error or no convergence, singular factorization, no positive real load multiplier…) a dense solver is used, if the dense problem has at most max_dense_size dofs, with a DenseFallbackWarning that gives the reasons of the failures. Its result is verified as well. A RuntimeError listing every attempted solver is raised when the problem is too large for the dense fallback, or when the dense solution also failed. The dense fallback keeps structsolve robust where ARPACK is unreliable, normal runs never pay for a dense solution.

Symmetric matrices. The sparse solvers are:

  1. scipy.sparse.linalg.eigsh() in Cayley mode, with the shift estimated from the matrices by _estimate_sigma();

  2. the same with a 10 times larger shift.

The dense fallback condenses out the dofs where KG has null rows, typically the in-plane dofs, with a sparse factorization of K, and solves the condensed problem with scipy.linalg.eigh(), see _eigh_condensed(). The eigenpairs are verified as follows:

  • the relative residual \(||K u + \lambda K_G u|| / (||K u|| + |\lambda| ||K_G u||)\) of each eigenpair must not be larger than check_rtol;

  • when K is positive definite, no load multiplier below the highest positive one returned may be missing: the number of negative pivots of \([K] + s [K_G]\), for \(s\) slightly below that multiplier, is the number of load multipliers in \((0, s)\) (Sylvester’s law of inertia). The factorization is done with SuperLU, without row interchanges. When a row interchange was needed and the count is unknown, \([K] + s [K_G]\) must be positive definite for \(s\) slightly below the lowest positive load multiplier. See _check_eigenpairs().

Unsymmetric matrices, e.g. KG with the load stiffness of a follower pressure that is not conservative, may have complex eigenvalues. The sparse solver is a search of shifts along the positive real axis with scipy.sparse.linalg.eigs(), see _eigs_shift_search(), which finds the lowest positive real load multiplier even when many negative or complex multipliers have a smaller modulus, verifies the residuals, and checks the parity of the number of real multipliers below the critical one. The returned eigenpairs are the num_eigvalues load multipliers nearest the last shift of the search, which is zero unless the multipliers of smallest modulus did not include a positive real one, not necessarily the load multipliers of smallest modulus. The dense fallback, _eig_condensed(), keeps the dofs where either the row or the column of KG is not null, which for follower loads are almost all of them, and returns the num_eigvalues first load multipliers in the order described above; when none of them is positive and real, i.e. no divergence load exists, e.g. for a system that only loses stability by flutter, they are returned with a warning. The inertia check requires symmetric matrices and is skipped. Eigenvalues whose imaginary part is below 1e-8 of their modulus are returned as real, with real eigenvectors. The arrays are real when the returned eigenpairs are real, even if other eigenvalues of the solver are complex.

Intel MKL. ARPACK, used by scipy.sparse.linalg.eigsh(), returns wrong eigenpairs at random, without any error, or raises ArpackError -8 when linked against Intel MKL 2024.2.0 to 2025.0.0, e.g. in some Anaconda builds of SciPy: the dsteqr routine of these MKL versions returns wrong eigenvectors for matrices larger than 32 x 32, which ARPACK uses with the size ncv. When such a LAPACK is detected, ncv is capped at 32 and more than 15 eigenpairs are computed in several passes, combined with a Rayleigh-Ritz projection, see structsolve.arpackutils, where the cap can be overridden, e.g. with the environment variable STRUCTSOLVE_ARPACK_MAX_NCV=0 to lift it with MKL 2025.0.1 or newer. The non-symmetric drivers of scipy.sparse.linalg.eigs() are not affected.

With Intel MKL 2025.0.0, the real LAPACK drivers dgeev and dggev of scipy.linalg.eig() crash the Python process for some matrices when computing the eigenvectors, with any number of threads. The dense solvers of unsymmetric matrices therefore call scipy.linalg.eig() with complex matrices, i.e. zgeev and zggev, which are not affected.

structsolve.linear_buckling.is_symmetric(A, rtol=1e-10)#

Tell whether ||A - A^T|| <= rtol ||A||, in the Frobenius norm

A cheap check, of the order of the number of non-zero entries of A, used by lb() to select the solvers of symmetric or unsymmetric matrices.

class structsolve.DenseFallbackWarning#

Issued by lb() and freq() when the sparse eigenvalue solvers failed and the dense solver was used as a last resort

See also Eigenvalue solvers for the solver policy.