NAG CL Interface
f02jcc (real_​gen_​quad)

Settings help

CL Name Style:


1 Purpose

f02jcc solves the quadratic eigenvalue problem
(λ2A+λB+C) x = 0 ,  
where A, B and C are real n×n matrices.
The function returns the 2n eigenvalues, λj, for j=1,2,,2n, and can optionally return the corresponding right eigenvectors, xj and/or left eigenvectors, yj as well as estimates of the condition numbers of the computed eigenvalues and backward errors of the computed right and left eigenvectors. A left eigenvector satisfies the equation
yH (λ2A+λB+C) = 0 ,  
where yH is the complex conjugate transpose of y.
λ is represented as the pair (α,β) , such that λ = α/β. Note that the computation of α/β may overflow and indeed β may be zero.

2 Specification

#include <nag.h>
void  f02jcc (Nag_ScaleType scal, Nag_LeftVecsType jobvl, Nag_RightVecsType jobvr, Nag_CondErrType sense, double tol, Integer n, double a[], Integer pda, double b[], Integer pdb, double c[], Integer pdc, double alphar[], double alphai[], double beta[], double vl[], Integer pdvl, double vr[], Integer pdvr, double s[], double bevl[], double bevr[], Integer *iwarn, NagError *fail)
The function may be called by the names: f02jcc or nag_eigen_real_gen_quad.

3 Description

The quadratic eigenvalue problem is solved by linearizing the problem and solving the resulting 2n×2n generalized eigenvalue problem. The linearization is chosen to have favourable conditioning and backward stability properties. An initial preprocessing step is performed that reveals and deflates the zero and infinite eigenvalues contributed by singular leading and trailing matrices.
The algorithm is backward stable for problems that are not too heavily damped, that is B10A·C.
Further details on the algorithm are given in Hammarling et al. (2013).

4 References

Fan H -Y, Lin W.-W and Van Dooren P. (2004) Normwise scaling of second order polynomial matrices SIAM J. Matrix Anal. Appl. 26, 1 252–256
Gaubert S and Sharify M (2009) Tropical scaling of polynomial matrices Lecture Notes in Control and Information Sciences Series 389 291–303 Springer–Verlag
Hammarling S, Munro C J and Tisseur F (2013) An algorithm for the complete solution of quadratic eigenvalue problems ACM Trans. Math. Software. 39(3):18:1–18:119 http://eprints.maths.manchester.ac.uk/2061/

5 Arguments

1: scal Nag_ScaleType Input
On entry: determines the form of scaling to be performed on A, B and C.
scal=Nag_NoScale
No scaling.
scal=Nag_CondFanLinVanDooren (the recommended value)
Fan, Lin and Van Dooren scaling if B A×C <10 and no scaling otherwise where Z is the Frobenius norm of Z.
scal=Nag_FanLinVanDooren
Fan, Lin and Van Dooren scaling.
scal=Nag_TropicalLargest
Tropical scaling with largest root.
scal=Nag_TropicalSmallest
Tropical scaling with smallest root.
Constraint: scal=Nag_NoScale, Nag_CondFanLinVanDooren, Nag_FanLinVanDooren, Nag_TropicalLargest or Nag_TropicalSmallest.
2: jobvl Nag_LeftVecsType Input
On entry: if jobvl=Nag_NotLeftVecs, do not compute left eigenvectors.
If jobvl=Nag_LeftVecs, compute the left eigenvectors.
If sense=Nag_CondOnly, Nag_BackErrLeft, Nag_BackErrBoth, Nag_CondBackErrLeft, Nag_CondBackErrRight or Nag_CondBackErrBoth, jobvl must be set to Nag_LeftVecs.
Constraint: jobvl=Nag_NotLeftVecs or Nag_LeftVecs.
3: jobvr Nag_RightVecsType Input
On entry: if jobvr=Nag_NotRightVecs, do not compute right eigenvectors.
If jobvr=Nag_RightVecs, compute the right eigenvectors.
If sense=Nag_CondOnly, Nag_BackErrRight, Nag_BackErrBoth, Nag_CondBackErrLeft, Nag_CondBackErrRight or Nag_CondBackErrBoth, jobvr must be set to Nag_RightVecs.
Constraint: jobvr=Nag_NotRightVecs or Nag_RightVecs.
4: sense Nag_CondErrType Input
On entry: determines whether, or not, condition numbers and backward errors are computed.
sense=Nag_NoCondBackErr
Do not compute condition numbers, or backward errors.
sense=Nag_CondOnly
Just compute condition numbers for the eigenvalues.
sense=Nag_BackErrLeft
Just compute backward errors for the left eigenpairs.
sense=Nag_BackErrRight
Just compute backward errors for the right eigenpairs.
sense=Nag_BackErrBoth
Compute backward errors for the left and right eigenpairs.
sense=Nag_CondBackErrLeft
Compute condition numbers for the eigenvalues and backward errors for the left eigenpairs.
sense=Nag_CondBackErrRight
Compute condition numbers for the eigenvalues and backward errors for the right eigenpairs.
sense=Nag_CondBackErrBoth
Compute condition numbers for the eigenvalues and backward errors for the left and right eigenpairs.
Constraint: sense=Nag_NoCondBackErr, Nag_CondOnly, Nag_BackErrLeft, Nag_BackErrRight, Nag_BackErrBoth, Nag_CondBackErrLeft, Nag_CondBackErrRight or Nag_CondBackErrBoth.
5: tol double Input
On entry: tol is used as the tolerance for making decisions on rank in the deflation procedure. If tol is zero on entry then n × machine precision is used in place of tol, where machine precision is as returned by function X02AJC. A diagonal element of a triangular matrix, R, is regarded as zero if |rjj| tol×size(X), or n×machine precision×size(X) when tol is zero, where size(X) is based on the size of the absolute values of the elements of the matrix X containing the matrix R. See Hammarling et al. (2013) for the motivation. If tol is -1.0 on entry then no deflation is attempted. The recommended value for tol is zero.
6: n Integer Input
On entry: the order of the matrices A, B and C.
Constraint: n0.
7: a[dim] double Input/Output
Note: the dimension, dim, of the array a must be at least pda×n.
The (i,j)th element of the matrix A is stored in a[(j-1)×pda+i-1].
On entry: the n×n matrix A.
On exit: a is used as internal workspace, but if jobvl=Nag_LeftVecs or jobvr=Nag_RightVecs, a is restored on exit.
8: pda Integer Input
On entry: the stride separating matrix row elements in the array a.
Constraint: pdan.
9: b[dim] double Input/Output
Note: the dimension, dim, of the array b must be at least pdb×n.
The (i,j)th element of the matrix B is stored in b[(j-1)×pdb+i-1].
On entry: the n×n matrix B.
On exit: b is used as internal workspace, but is restored on exit.
10: pdb Integer Input
On entry: the stride separating matrix row elements in the array b.
Constraint: pdbn.
11: c[dim] double Input/Output
Note: the dimension, dim, of the array c must be at least pdc×n.
The (i,j)th element of the matrix C is stored in c[(j-1)×pdc+i-1].
On entry: the n×n matrix C.
On exit: c is used as internal workspace, but if jobvl=Nag_LeftVecs or jobvr=Nag_RightVecs, c is restored on exit.
12: pdc Integer Input
On entry: the stride separating matrix row elements in the array c.
Constraint: pdcn.
13: alphar[2×n] double Output
On exit: alphar[j-1], for j=1,2,,2n, contains the real part of αj for the jth eigenvalue pair (αj,βj) of the quadratic eigenvalue problem.
14: alphai[2×n] double Output
On exit: alphai[j-1], for j=1,2,,2n, contains the imaginary part of αj for the jth eigenvalue pair (αj,βj) of the quadratic eigenvalue problem. If alphai[j-1] is zero then the jth eigenvalue is real; if alphai[j-1] is positive then the jth and (j+1)th eigenvalues are a complex conjugate pair, with alphai[j] negative.
15: beta[2×n] double Output
On exit: beta[j-1], for j=1,2,,2n, contains the second part of the jth eigenvalue pair (αj,βj) of the quadratic eigenvalue problem, with βj0. Infinite eigenvalues have βj set to zero.
16: vl[dim] double Output
Note: the dimension, dim, of the array vl must be at least 2×n when jobvl=Nag_LeftVecs.
where VL(i,j) appears in this document, it refers to the array element vl[(j-1)×pdvl+i-1].
On exit: if jobvl=Nag_LeftVecs, the left eigenvectors yj are stored one after another in vl, in the same order as the corresponding eigenvalues. If the jth eigenvalue is real, then yj=VL(:,j), the jth column of VL. If the jth and (j+1)th eigenvalues form a complex conjugate pair, then yj= VL(:,j) + i× VL(:,j+1) and yj+1 = VL(:,j) -i × VL(:,j+1) . Each eigenvector will be normalized with length unity and with the element of largest modulus real and positive.
If jobvl=Nag_NotLeftVecs, vl is not referenced and may be NULL.
17: pdvl Integer Input
On entry: the stride separating matrix row elements in the array vl.
Constraint: pdvln.
18: vr[dim] double Output
Note: the dimension, dim, of the array vr must be at least 2×n when jobvr=Nag_RightVecs.
where VR(i,j) appears in this document, it refers to the array element vr[(j-1)×pdvr+i-1].
On exit: if jobvr=Nag_RightVecs, the right eigenvectors xj are stored one after another in vr, in the same order as the corresponding eigenvalues. If the jth eigenvalue is real, then xj = VR(:,j) , the jth column of VR. If the jth and (j+1)th eigenvalues form a complex conjugate pair, then xj = VR(:,j) +i× VR(:,j+1) and xj+1= VR(:,j) -i× VR(:,j+1) . Each eigenvector will be normalized with length unity and with the element of largest modulus real and positive.
If jobvr=Nag_NotRightVecs, vr is not referenced and may be NULL.
19: pdvr Integer Input
On entry: the stride separating matrix row elements in the array vr.
Constraint: pdvrn.
20: s[dim] double Output
Note: the dimension, dim, of the array s must be at least
  • 2×n when sense=Nag_CondOnly, Nag_CondBackErrLeft, Nag_CondBackErrRight or Nag_CondBackErrBoth.
also: computing the condition numbers of the eigenvalues requires that both the left and right eigenvectors be computed.
On exit: if sense=Nag_CondOnly, Nag_CondBackErrLeft, Nag_CondBackErrRight or Nag_CondBackErrBoth, s[j-1] contains the condition number estimate for the jth eigenvalue (large condition numbers imply that the problem is near one with multiple eigenvalues). Infinite condition numbers are returned as the largest model double number (X02ALC).
If sense=Nag_NoCondBackErr, Nag_BackErrLeft, Nag_BackErrRight or Nag_BackErrBoth, s is not referenced and may be NULL.
21: bevl[dim] double Output
Note: the dimension, dim, of the array bevl must be at least
  • 2×n when sense=Nag_BackErrLeft, Nag_BackErrBoth, Nag_CondBackErrLeft or Nag_CondBackErrBoth.
On exit: if sense=Nag_BackErrLeft, Nag_BackErrBoth, Nag_CondBackErrLeft or Nag_CondBackErrBoth, bevl[j-1] contains the backward error estimate for the computed left eigenpair (λj,yj) .
If sense=Nag_NoCondBackErr, Nag_CondOnly, Nag_BackErrRight or Nag_CondBackErrRight, bevl is not referenced and may be NULL.
22: bevr[dim] double Output
Note: the dimension, dim, of the array bevr must be at least
  • 2×n when sense=Nag_BackErrRight, Nag_BackErrBoth, Nag_CondBackErrRight or Nag_CondBackErrBoth.
On exit: if sense=Nag_BackErrRight, Nag_BackErrBoth, Nag_CondBackErrRight or Nag_CondBackErrBoth, bevr[j-1] contains the backward error estimate for the computed right eigenpair (λj,xj) .
If sense=Nag_NoCondBackErr, Nag_CondOnly, Nag_BackErrLeft or Nag_CondBackErrLeft, bevr is not referenced and may be NULL.
23: iwarn Integer * Output
On exit: iwarn will be positive if there are warnings, otherwise iwarn will be 0.
If fail.code= NE_NOERROR then:
  • if iwarn=1 then one, or both, of the matrices A and C is zero. In this case no scaling is performed, even if scal=Nag_CondFanLinVanDooren;
  • if iwarn=2 then the matrices A and C are singular, or nearly singular, so the problem is potentially ill-posed;
  • if iwarn=3 then both the conditions for iwarn=1 and iwarn=2 above, apply. If iwarn=4, b10A·C and backward stability cannot be guaranteed.
If fail.code= NE_ITERATIONS_QZ, f08xcc has flagged that iwarn eigenvalues are invalid.
If fail.code= NE_ITERATIONS_QZ, f08wcc has flagged that iwarn eigenvalues are invalid.
24: fail NagError * Input/Output
The NAG error argument (see Section 7 in the Introduction to the NAG Library CL Interface).

6 Error Indicators and Warnings

NE_ALLOC_FAIL
Dynamic memory allocation failed.
See Section 3.1.2 in the Introduction to the NAG Library CL Interface for further information.
NE_ARRAY_SIZE
On entry, pda=value and n=value.
Constraint: pdan.
On entry, pdb=value and n=value.
Constraint: pdbn.
On entry, pdc=value and n=value.
Constraint: pdcn.
On entry, pdvl=value, n=value and jobvl=value.
Constraint: when jobvl=Nag_LeftVecs, pdvln.
On entry, pdvr=value, n=value and jobvr=value.
Constraint: when jobvr=Nag_RightVecs, pdvrn.
NE_BAD_PARAM
On entry, argument value had an illegal value.
NE_INT
On entry, n=value.
Constraint: n0.
NE_INTERNAL_ERROR
An internal error has occurred in this function. Check the function call and any array sizes. If the call is correct then please contact NAG for assistance.
See Section 7.5 in the Introduction to the NAG Library CL Interface for further information.
NE_INVALID_VALUE
On entry, sense=value and jobvl=value.
Constraint: when jobvl=Nag_NotLeftVecs, sense=Nag_NoCondBackErr or Nag_BackErrRight,
when jobvl=Nag_LeftVecs, sense=Nag_CondOnly, Nag_BackErrLeft, Nag_BackErrBoth, Nag_CondBackErrLeft, Nag_CondBackErrRight or Nag_CondBackErrBoth.
On entry, sense=value and jobvr=value.
Constraint: when jobvr=Nag_NotRightVecs, sense=Nag_NoCondBackErr or Nag_BackErrLeft,
when jobvr=Nag_RightVecs, sense=Nag_CondOnly, Nag_BackErrRight, Nag_BackErrBoth, Nag_CondBackErrLeft, Nag_CondBackErrRight or Nag_CondBackErrBoth.
NE_ITERATIONS_QZ
The QZ iteration failed in f08wcc.
iwarn returns the value of info returned by f08wcc. This failure is unlikely to happen, but if it does, please contact NAG.
The QZ iteration failed in f08xcc.
iwarn returns the value of info returned by f08xcc. This failure is unlikely to happen, but if it does, please contact NAG.
NE_NO_LICENCE
Your licence key may have expired or may not have been installed correctly.
See Section 8 in the Introduction to the NAG Library CL Interface for further information.
NE_SINGULAR
The quadratic matrix polynomial is nonregular (singular).

7 Accuracy

The algorithm is backward stable for problems that are not too heavily damped, that is B10 A·C .

8 Parallelism and Performance

f02jcc is threaded by NAG for parallel execution in multithreaded implementations of the NAG Library.
f02jcc makes calls to BLAS and/or LAPACK routines, which may be threaded within the vendor library used by this implementation. Consult the documentation for the vendor library for further information.
Please consult the X06 Chapter Introduction for information on how to control and interrogate the OpenMP environment used within this function. Please also consult the Users' Note for your implementation for any additional implementation-specific information.

9 Further Comments

None.

10 Example

To solve the quadratic eigenvalue problem
(λ2A+λB+C) x = 0  
where
A = ( 1 2 2 3 1 1 3 2 1 ) ,  B = ( 3 2 1 2 1 3 1 3 2 )  and  C = ( 1 1 2 2 3 1 3 1 2 ) .  
The example also returns the left eigenvectors, condition numbers for the computed eigenvalues and backward errors of the computed right and left eigenpairs.

10.1 Program Text

Program Text (f02jcce.c)

10.2 Program Data

Program Data (f02jcce.d)

10.3 Program Results

Program Results (f02jcce.r)