linalg  1.4.3
A linear algebra library that provides a user-friendly interface to several BLAS and LAPACK routines.
linalg_eigen Module Reference

linalg_eigen More...

Data Types

interface  eigen
 Computes the eigenvalues, and optionally the eigenvectors, of a matrix. More...
 

Functions/Subroutines

subroutine eigen_symm (vecs, a, vals, work, olwork, err)
 Computes the eigenvalues, and optionally the eigenvectors of a real, symmetric matrix. More...
 
subroutine eigen_asymm (a, vals, vecs, work, olwork, err)
 Computes the eigenvalues, and optionally the right eigenvectors of a square matrix. More...
 
subroutine eigen_gen (a, b, alpha, beta, vecs, work, olwork, err)
 Computes the eigenvalues, and optionally the right eigenvectors of a square matrix assuming the structure of the eigenvalue problem is A*X = lambda*B*X. More...
 

Detailed Description

linalg_eigen

Purpose
Provides routines for computing the eigenvalues and eigenvectors of matrices.

Function/Subroutine Documentation

◆ eigen_asymm()

subroutine linalg_eigen::eigen_asymm ( real(real64), dimension(:,:), intent(inout)  a,
complex(real64), dimension(:), intent(out)  vals,
complex(real64), dimension(:,:), intent(out), optional  vecs,
real(real64), dimension(:), intent(out), optional, pointer  work,
integer(int32), intent(out), optional  olwork,
class(errors), intent(inout), optional, target  err 
)
private

Computes the eigenvalues, and optionally the right eigenvectors of a square matrix.

Parameters
[in,out]aOn input, the N-by-N matrix on which to operate. On output, the contents of this matrix are overwritten.
[out]valsAn N-element array containing the eigenvalues of the matrix. The eigenvalues are not sorted.
[out]vecsAn optional N-by-N matrix, that if supplied, signals to compute the right eigenvectors (one per column). If not provided, only the eigenvalues will be computed.
[out]workAn optional input, that if provided, prevents any local memory allocation. If not provided, the memory required is allocated within. If provided, the length of the array must be at least olwork.
[out]olworkAn optional output used to determine workspace size. If supplied, the routine determines the optimal size for work, and returns without performing any actual calculations.
[out]errAn optional errors-based object that if provided can be used to retrieve information relating to any errors encountered during execution. If not provided, a default implementation of the errors class is used internally to provide error handling. Possible errors and warning messages that may be encountered are as follows.
  • LA_ARRAY_SIZE_ERROR: Occurs if any of the input arrays are not sized appropriately.
  • LA_OUT_OF_MEMORY_ERROR: Occurs if local memory must be allocated, and there is insufficient memory available.
  • LA_CONVERGENCE_ERROR: Occurs if the algorithm failed to converge.
Notes
This routine utilizes the LAPACK routine DGEEV.

Definition at line 186 of file linalg_eigen.f90.

◆ eigen_gen()

subroutine linalg_eigen::eigen_gen ( real(real64), dimension(:,:), intent(inout)  a,
real(real64), dimension(:,:), intent(inout)  b,
complex(real64), dimension(:), intent(out)  alpha,
real(real64), dimension(:), intent(out), optional  beta,
complex(real64), dimension(:,:), intent(out), optional  vecs,
real(real64), dimension(:), intent(out), optional, pointer  work,
integer(int32), intent(out), optional  olwork,
class(errors), intent(inout), optional, target  err 
)
private

Computes the eigenvalues, and optionally the right eigenvectors of a square matrix assuming the structure of the eigenvalue problem is A*X = lambda*B*X.

Parameters
[in,out]aOn input, the N-by-N matrix A. On output, the contents of this matrix are overwritten.
[in,out]bOn input, the N-by-N matrix B. On output, the contents of this matrix are overwritten.
[out]alphaAn N-element array that, if beta is not supplied, contains the eigenvalues. If beta is supplied however, the eigenvalues must be computed as ALPHA / BETA. This however, is not as trivial as it seems as it is entirely possible, and likely, that ALPHA / BETA can overflow or underflow. With that said, the values in ALPHA will always be less than and usually comparable with the NORM(A).
[out]betaAn optional N-element array that if provided forces alpha to return the numerator, and this array contains the denominator used to determine the eigenvalues as ALPHA / BETA. If used, the values in this array will always be less than and usually comparable with the NORM(B).
[out]vecsAn optional N-by-N matrix, that if supplied, signals to compute the right eigenvectors (one per column). If not provided, only the eigenvalues will be computed.
[out]workAn optional input, that if provided, prevents any local memory allocation. If not provided, the memory required is allocated within. If provided, the length of the array must be at least olwork.
[out]olworkAn optional output used to determine workspace size. If supplied, the routine determines the optimal size for work, and returns without performing any actual calculations.
[out]errAn optional errors-based object that if provided can be used to retrieve information relating to any errors encountered during execution. If not provided, a default implementation of the errors class is used internally to provide error handling. Possible errors and warning messages that may be encountered are as follows.
  • LA_ARRAY_SIZE_ERROR: Occurs if any of the input arrays are not sized appropriately.
  • LA_OUT_OF_MEMORY_ERROR: Occurs if local memory must be allocated, and there is insufficient memory available.
  • LA_CONVERGENCE_ERROR: Occurs if the algorithm failed to converge.
Usage
As an example, consider the eigenvalue problem arising from a mechanical system of masses and springs such that the masses are described by a mass matrix M, and the arrangement of springs are described by a stiffness matrix K.
! This is an example illustrating the use of the eigenvalue and eigenvector
! routines to solve a free vibration problem of 3 masses connected by springs.
!
! k1 k2 k3 k4
! |-\/\/\-| m1 |-\/\/\-| m2 |-\/\/\-| m3 |-\/\/\-|
!
! As illustrated above, the system consists of 3 masses connected by springs.
! Spring k1 and spring k4 connect the end masses to ground. The equations of
! motion for this system are as follows.
!
! | m1 0 0 | |x1"| | k1+k2 -k2 0 | |x1| |0|
! | 0 m2 0 | |x2"| + | -k2 k2+k3 -k3 | |x2| = |0|
! | 0 0 m3| |x3"| | 0 -k3 k3+k4| |x3| |0|
!
! Notice: x1" = the second time derivative of x1.
program example
use linalg_constants, only : dp, int32
implicit none
! Define the model parameters
real(real64), parameter :: pi = 3.14159265359d0
real(real64), parameter :: m1 = 0.5d0
real(real64), parameter :: m2 = 2.5d0
real(real64), parameter :: m3 = 0.75d0
real(real64), parameter :: k1 = 5.0d6
real(real64), parameter :: k2 = 10.0d6
real(real64), parameter :: k3 = 10.0d6
real(real64), parameter :: k4 = 5.0d6
! Local Variables
integer(int32) :: i, j
real(real64) :: m(3,3), k(3,3), natFreq(3)
complex(real64) :: vals(3), modeShapes(3,3)
! Define the mass matrix
m = reshape([m1, 0.0d0, 0.0d0, 0.0d0, m2, 0.0d0, 0.0d0, 0.0d0, m3], [3, 3])
! Define the stiffness matrix
k = reshape([k1 + k2, -k2, 0.0d0, -k2, k2 + k3, -k3, 0.0d0, -k3, k3 + k4], &
[3, 3])
! Compute the eigenvalues and eigenvectors.
call eigen(k, m, vals, vecs = modeshapes)
! Compute the natural frequency values, and return them with units of Hz.
! Notice, all eigenvalues and eigenvectors are real for this example.
natfreq = sqrt(real(vals)) / (2.0d0 * pi)
! Display the natural frequency and mode shape values. Notice, the eigen
! routine does not necessarily sort the values.
print '(A)', "Modal Information (Not Sorted):"
do i = 1, size(natfreq)
print '(AI0AF8.4A)', "Mode ", i, ": (", natfreq(i), " Hz)"
print '(F10.3)', (real(modeShapes(j,i)), j = 1, size(natfreq))
end do
end program
Notes
This routine utilizes the LAPACK routine DGGEV.

Definition at line 466 of file linalg_eigen.f90.

◆ eigen_symm()

subroutine linalg_eigen::eigen_symm ( logical, intent(in)  vecs,
real(real64), dimension(:,:), intent(inout)  a,
real(real64), dimension(:), intent(out)  vals,
real(real64), dimension(:), intent(out), optional, pointer  work,
integer(int32), intent(out), optional  olwork,
class(errors), intent(inout), optional, target  err 
)
private

Computes the eigenvalues, and optionally the eigenvectors of a real, symmetric matrix.

Parameters
[in]vecsSet to true to compute the eigenvectors as well as the eigenvalues; else, set to false to just compute the eigenvalues.
[in,out]aOn input, the N-by-N symmetric matrix on which to operate. On output, and if vecs is set to true, the matrix will contain the eigenvectors (one per column) corresponding to each eigenvalue in vals. If vecs is set to false, the lower triangular portion of the matrix is overwritten.
[out]valsAn N-element array that will contain the eigenvalues sorted into ascending order.
[out]workAn optional input, that if provided, prevents any local memory allocation. If not provided, the memory required is allocated within. If provided, the length of the array must be at least olwork.
[out]olworkAn optional output used to determine workspace size. If supplied, the routine determines the optimal size for work, and returns without performing any actual calculations.
[out]errAn optional errors-based object that if provided can be used to retrieve information relating to any errors encountered during execution. If not provided, a default implementation of the errors class is used internally to provide error handling. Possible errors and warning messages that may be encountered are as follows.
  • LA_ARRAY_SIZE_ERROR: Occurs if any of the input arrays are not sized appropriately.
  • LA_OUT_OF_MEMORY_ERROR: Occurs if local memory must be allocated, and there is insufficient memory available.
  • LA_CONVERGENCE_ERROR: Occurs if the algorithm failed to converge.
Notes
This routine utilizes the LAPACK routine DSYEV.

Definition at line 69 of file linalg_eigen.f90.