3.2.9. eigen Command

This command is used to perform the eigenvalue analysis.

eigen <type> <$solver> $numEigenvalues

Argument	Type	Description
$numEigenvalues	integer	number of eigenvalues required.
$solver	string	optional string detailing type of solver: -genBandArpack (default), -fullGenLapack, -symmBandLapack. Python users may also use the PythonSparse eigenvalue solver.
$type	string	optional string indicating type of eigenvalue problem to solve: ‘general’ (default) or ‘standard’

Returns

A list containing the eigenvalues

Note

The eigenvectors are stored at the nodes and can be printed out using a Node Recorder, the nodeEigenvector command, or the Print command.
The default eigensolver is able to solve only for N-1 eigenvalues, where N is the number of inertial DOFs. When running into this limitation the -fullGenLapack solver can be used instead of the default Arpack solver.
The -symmBandLapack option works only standard eigenvalue analysis of the stiffness matrix, i.e., K*x = lam*x
Python users: A PythonSparse eigenvalue solver is available as well (OpenSeesPy only).

Warning

The default eigen solver utilizes ARPACK, an iterative eigenvalue solver. Like other iterative methods, ARPACK can struggle to accurately compute eigenvectors when the model contains repeated (degenerate) eigenvalues or clusters of very closely spaced eigenvalues. The issue arises from difficulties in forming an orthonormal basis for the eigenspace numerically.

In such cases, ARPACK may fail to find all of the requested eigenvalues, or it may return correct eigenvalues but incorrect eigenvectors. When modal damping is used, inaccurate eigenvectors may cause issues with the damping forces calculated. If you are unsure if this is an issue with your model, switch the solver used and check the reslts again – this check works as the eigenvectors will vary greatly with the chosen solver when this issue exists.

Fortunately, small adjustments to the model can prevent this numerical issue. In actuality, as material properties (e.g., Young’s modulus) or nodal masses are rarely identical and by introducing slight variations in your model, ARPACK is more likely to form the basis correctly.

As an illustration of the problem, consider the following script (provided by mhscott):

model Basic -ndm 1 -ndf 1
node 0 0; fix 0 1
node 1 0; mass 1 1.0
node 2 0; mass 2 1.0
node 3 0; mass 3 1.0
node 4 0; mass 4 1.0

set k 610
uniaxialMaterial Elastic 1 $k

element zeroLength 1 0 1 -mat 1 -dir 1
element zeroLength 2 1 2 -mat 1 -dir 1
element zeroLength 3 0 3 -mat 1 -dir 1
element zeroLength 4 3 4 -mat 1 -dir 1

systetm ProfileSPD
numberer RCM

set nModes 2
set eigV [eigen $nModes]
puts "Eigenvalues: $eigV"

# print eigenvectors
for {set i 1} {$i <= 4} {incr i 1}  {
   puts "$i [nodeEigenvector $i 1] [nodeEigenvector $i 2]"
}

will produce:

Eigenvalues:  232.99926686256409880116   232.99926686256412722287
  0.50235827201282656773    0.15501409223134476889
  0.81283275864641846287    0.25081806996552691302
  0.15501409223134465787   -0.50235827201282667875
  0.25081806996552685751   -0.81283275864641857389

wheras if I change the mass at node 4 using:

node 4 0; mass 4 1.0000000001

I would get:

Eigenvalues:  232.99926679816474006657   232.99926686256424090971
 -0.00000004574590262152    0.52573111211913170493
 -0.00000007401842538890    0.85065080835203654708
  0.52573111211146139610    0.00000004574590273254
  0.85065080819431726500    0.00000007401842544441

For this in an academic setting dealing with small problems, the -fullGenLapack option will provide correct eigenvalues and eigenvectors to the original problem:

will produce:

Eigenvalues:   232.99926684570411339337   232.99926686256412722287
  0.00000000000000000000   -0.52573111211913370333
  0.00000000000000000000   -0.85065080835203987775
  0.52573111209361389484    0.00000000000000000000
  0.85065080832527950605   -0.00000000000000000000

3.2.9.1. Theory

A generalized eigenvalue problem for two symmetric matrices $K$ and $M$ of size $n \times n$ is given by:

\[\left (K - \lambda M \right ) \Phi = 0\]

where:

$K$ is the stiffness matrix

$M$ is the mass matrix

$\lambda$ is the eigenvalue

and $\Phi$ is the associated eigenvector

Example

The following example shows how to use the eigen command to obtain a list of eigenvalues.

Tcl Code

# obtain 10 eigenvalues using the default solver (-genBandArpack)
set eigenvalues [eigen 10]

# or, obtain 10 eigenvalues explicitly specifying the solver
set eigenvalues [eigen -fullGenLapack 10]

# obtain 10 eigenvalues of the stiffness matrix
set eigenvalues [eigen standard -symmBandLapack 10]

Python Code

# obtain 10 eigenvalues using the default solver (-genBandArpack)
eigenvalues = eigen(10)

# or, obtain 10 eigenvalues explicitly specifying the solver
eigenvalues = eigen('-fullGenLapack', 10)

# obtain 10 eigenvalues of the stiffness matrix
eigenvalues = eigen('standard','-symmBandLapack',10)

# PythonSparse: use a custom Python solver (e.g., SciPy eigsh)
eigenvalues = eigen('PythonSparse', 10, {'solver': solver, 'scheme': 'CSR'})

3.2.9.2. PythonSparse Eigenvalue Solver

The PythonSparse eigen solver delegates the generalized eigenvalue problem $\mathbf{K}\phi = \lambda \mathbf{M}\phi$ to user-defined Python objects. Sparse stiffness (K) and mass (M) matrix buffers are exposed via zero-copy memoryviews, enabling the use of SciPy, CuPy, or custom eigensolvers without recompiling OpenSees. See example below for a SciPy eigsh implementation.

eigen PythonSparse $numEigenvalues <dict>?: Perform eigenvalue analysis using a user-defined Python solver. The dict is the solver configuration; keys below.

Argument	Type	Description
$numEigenvalues	int	Number of eigenvalues (and eigenvectors) to compute.
dict	dict	Solver configuration; see dict options table below.

Key	Type	Description
solver	Python object	A Python object with a `solve(**kwargs)` method. Required.
scheme	string	Sparse format: `'CSR'`, `'CSC'`, `'COO'`. Default: `'CSR'`.

Note

The PythonSparse eigen interface is available only when using the Python interpreter (OpenSeesPy).

The solver object must implement a solve(**kwargs) method. OpenSees passes the matrix data as keyword arguments: memoryviews (buffer-like objects) for arrays, and scalars for metadata. Use numpy.frombuffer to interpret memoryviews as NumPy arrays without copying. The structure format depends on scheme:

CSR/CSC format — compressed storage uses index_ptr and indices:

Keyword	Type	Description
index_ptr	memoryview	Row pointers (CSR) or column pointers (CSC), int32.
indices	memoryview	Column indices (CSR) or row indices (CSC), int32.
k_values	memoryview	Stiffness matrix K coefficients (float64).
m_values	memoryview	Mass matrix M coefficients (float64).
eigenvalues	memoryview	Output buffer for eigenvalues (float64, writable). Write in place.
eigenvectors	memoryview	Output buffer for eigenvectors (float64, writable). Write in place, row-major: mode 0, mode 1, …
num_eqn	int	Number of equations.
nnz	int	Number of nonzeros (K and M share sparsity pattern).
matrix_status	string	`'UNCHANGED'`, `'COEFFICIENTS_CHANGED'`, or `'STRUCTURE_CHANGED'`.
storage_scheme	string	`'CSR'`, `'CSC'`, or `'COO'`.
num_modes	int	Number of modes requested.
generalized	bool	True for generalized (K, M); False for standard (K only).
find_smallest	bool	True: smallest eigenvalues (e.g., frequency); False: largest.

COO format — coordinate storage uses row_indices and col_indices instead of index_ptr/indices:

Keyword	Type	Description
row_indices	memoryview	Row indices for each nonzero, int32 (COO only).
col_indices	memoryview	Column indices for each nonzero, int32 (COO only).
k_values	memoryview	Stiffness matrix K coefficients (float64).
m_values	memoryview	Mass matrix M coefficients (float64).
eigenvalues	memoryview	Output buffer for eigenvalues (float64, writable). Write in place.
eigenvectors	memoryview	Output buffer for eigenvectors (float64, writable). Write in place.
num_eqn	int	Number of equations.
nnz	int	Number of nonzeros.
matrix_status	string	`'UNCHANGED'`, `'COEFFICIENTS_CHANGED'`, or `'STRUCTURE_CHANGED'`.
storage_scheme	string	`'COO'`.
num_modes	int	Number of modes requested.
generalized	bool	True for generalized (K, M); False for standard (K only).
find_smallest	bool	True: smallest eigenvalues; False: largest.

The solve method should return None on success, or raise an exception on failure.

Warning

In-place output required

You must write eigenvalues and eigenvectors directly into the eigenvalues and eigenvectors buffers. OpenSees reads from these buffers; it does not use return values.

Do: Use in-place assignment, e.g. np.frombuffer(eigenvalues, ...)[:] = eigvals and np.frombuffer(eigenvectors, ...)[:] = eigvecs.T.flatten()

Do not: Return the results or assign to new arrays without copying back into the buffers. OpenSees will use uninitialized or stale data.

Example — SciPy Generalized Eigenvalue Solver

This example uses SciPy and its eigsh (ARPACK symmetric eigensolver) to solve the generalized eigenvalue problem:

import numpy as np
import scipy.sparse as sp
import scipy.sparse.linalg as sp_linalg
import openseespy.opensees as ops

class SciPyGeneralizedEigenSolver:
    def __init__(self, maxiter=None, tol=0.0):
        self.maxiter = maxiter
        self.tol = tol
        self._k_matrix = None
        self._m_matrix = None

    def solve(self, **kwargs):
        index_ptr = kwargs['index_ptr']
        indices = kwargs['indices']
        k_values = kwargs['k_values']
        m_values = kwargs['m_values']
        eigenvalues = kwargs['eigenvalues']
        eigenvectors = kwargs['eigenvectors']
        num_eqn = kwargs['num_eqn']
        nnz = kwargs['nnz']
        matrix_status = kwargs['matrix_status']
        num_modes = kwargs['num_modes']
        find_smallest = kwargs['find_smallest']

        indptr = np.frombuffer(index_ptr, dtype=np.int32, count=num_eqn + 1)
        idx = np.frombuffer(indices, dtype=np.int32, count=nnz)
        k_data = np.frombuffer(k_values, dtype=np.float64, count=nnz)
        m_data = np.frombuffer(m_values, dtype=np.float64, count=nnz)

        if matrix_status == 'STRUCTURE_CHANGED' or self._k_matrix is None:
            self._k_matrix = sp.csr_matrix((k_data, idx, indptr),
                                          shape=(num_eqn, num_eqn))
            self._m_matrix = sp.csr_matrix((m_data, idx, indptr),
                                          shape=(num_eqn, num_eqn))
        elif matrix_status == 'COEFFICIENTS_CHANGED':
            self._k_matrix.data[:] = k_data
            self._m_matrix.data[:] = m_data

        eigsh_kwargs = {
            'k': num_modes,
            'M': self._m_matrix,
            'which': 'LM',
        }
        if find_smallest:
            eigsh_kwargs['sigma'] = 0.0
        if self.maxiter is not None:
            eigsh_kwargs['maxiter'] = self.maxiter
        if self.tol > 0.0:
            eigsh_kwargs['tol'] = self.tol

        eigvals, eigvecs = sp_linalg.eigsh(self._k_matrix, **eigsh_kwargs)

        np.frombuffer(eigenvalues, dtype=np.float64,
                      count=num_modes)[:] = eigvals[:num_modes]
        np.frombuffer(eigenvectors, dtype=np.float64,
                      count=num_modes * num_eqn)[:] = eigvecs.T.flatten()

        return None

solver = SciPyGeneralizedEigenSolver()
eigenvalues = ops.eigen('PythonSparse', 10, {'solver': solver, 'scheme': 'CSR'})