NAG Toolbox: nag_lapack_dggevx (f08wb)

Purpose

Syntax

Description

References

Parameters

Compulsory Input Parameters

1:     $\mathrm{balanc}$ – string (length ≥ 1)

Specifies the balance option to be performed.

$\mathbf{balanc}=\text{'N'}$

Do not diagonally scale or permute.

$\mathbf{balanc}=\text{'P'}$

Permute only.

$\mathbf{balanc}=\text{'S'}$

Scale only.

$\mathbf{balanc}=\text{'B'}$

Both permute and scale.

Computed reciprocal condition numbers will be for the matrices after permuting and/or balancing. Permuting does not change condition numbers (in exact arithmetic), but balancing does. In the absence of other information, $\mathbf{balanc}=\text{'B'}$ is recommended.

Constraint: $\mathbf{balanc}=\text{'N'}$, $\text{'P'}$, $\text{'S'}$ or $\text{'B'}$.

2:     $\mathrm{jobvl}$ – string (length ≥ 1)

If $\mathbf{jobvl}=\text{'N'}$, do not compute the left generalized eigenvectors.

If $\mathbf{jobvl}=\text{'V'}$, compute the left generalized eigenvectors.

Constraint: $\mathbf{jobvl}=\text{'N'}$ or $\text{'V'}$.

3:     $\mathrm{jobvr}$ – string (length ≥ 1)

If $\mathbf{jobvr}=\text{'N'}$, do not compute the right generalized eigenvectors.

If $\mathbf{jobvr}=\text{'V'}$, compute the right generalized eigenvectors.

Constraint: $\mathbf{jobvr}=\text{'N'}$ or $\text{'V'}$.

4:     $\mathrm{sense}$ – string (length ≥ 1)

Determines which reciprocal condition numbers are computed.

$\mathbf{sense}=\text{'N'}$

None are computed.

$\mathbf{sense}=\text{'E'}$

Computed for eigenvalues only.

$\mathbf{sense}=\text{'V'}$

Computed for eigenvectors only.

$\mathbf{sense}=\text{'B'}$

Computed for eigenvalues and eigenvectors.

Constraint: $\mathbf{sense}=\text{'N'}$, $\text{'E'}$, $\text{'V'}$ or $\text{'B'}$.

5:     $\mathrm{a}\left(\mathit{lda},:\right)$ – double array

The first dimension of the array a must be at least $\max \left(1,\mathbf{n}\right)$.

The second dimension of the array a must be at least $\max \left(1,\mathbf{n}\right)$.

The matrix $A$ in the pair $\left(A,B\right)$.

6:     $\mathrm{b}\left(\mathit{ldb},:\right)$ – double array

The first dimension of the array b must be at least $\max \left(1,\mathbf{n}\right)$.

The second dimension of the array b must be at least $\max \left(1,\mathbf{n}\right)$.

The matrix $B$ in the pair $\left(A,B\right)$.

Optional Input Parameters

1:     $\mathrm{n}$ – int64int32nag_int scalar

Default: the first dimension of the arrays a, b and the second dimension of the arrays a, b. (An error is raised if these dimensions are not equal.)

$n$, the order of the matrices $A$ and $B$.

Constraint: $\mathbf{n}\ge 0$.

Output Parameters

1:     $\mathrm{a}\left(\mathit{lda},:\right)$ – double array

The first dimension of the array a will be $\max \left(1,\mathbf{n}\right)$.

The second dimension of the array a will be $\max \left(1,\mathbf{n}\right)$.

a has been overwritten. If $\mathbf{jobvl}=\text{'V'}$ or $\mathbf{jobvr}=\text{'V'}$ or both, then $A$ contains the first part of the real Schur form of the ‘balanced’ versions of the input $A$ and $B$.

2:     $\mathrm{b}\left(\mathit{ldb},:\right)$ – double array

The first dimension of the array b will be $\max \left(1,\mathbf{n}\right)$.

The second dimension of the array b will be $\max \left(1,\mathbf{n}\right)$.

b has been overwritten.

3:     $\mathrm{alphar}\left(\mathbf{n}\right)$ – double array

The element $\mathbf{alphar}\left(j\right)$ contains the real part of ${\alpha }_j$.

4:     $\mathrm{alphai}\left(\mathbf{n}\right)$ – double array

The element $\mathbf{alphai}\left(j\right)$ contains the imaginary part of ${\alpha }_j$.

5:     $\mathrm{beta}\left(\mathbf{n}\right)$ – double array

$\left(\mathbf{alphar}\left(\mathit{j}\right)+\mathbf{alphai}\left(\mathit{j}\right)\times i\right)/\mathbf{beta}\left(\mathit{j}\right)$, for $\mathit{j}=1,2,\dots ,\mathbf{n}$, will be the generalized eigenvalues.

If $\mathbf{alphai}\left(j\right)$ is zero, then the $j$th eigenvalue is real; if positive, then the $j$th and $\left(j+1\right)$st eigenvalues are a complex conjugate pair, with $\mathbf{alphai}\left(j+1\right)$ negative.

Note:  the quotients $\mathbf{alphar}\left(j\right)/\mathbf{beta}\left(j\right)$ and $\mathbf{alphai}\left(j\right)/\mathbf{beta}\left(j\right)$ may easily overflow or underflow, and $\mathbf{beta}\left(j\right)$ may even be zero. Thus, you should avoid naively computing the ratio ${\alpha }_j/{\beta }_j$. However, $\max \left|{\alpha }_j\right|$ will always be less than and usually comparable with ${‖\mathbf{a}‖}_2$ in magnitude, and $\max \left|{\beta }_j\right|$ will always be less than and usually comparable with ${‖\mathbf{b}‖}_2$.

6:     $\mathrm{vl}\left(\mathit{ldvl},:\right)$ – double array

The first dimension, $\mathit{ldvl}$, of the array vl will be

• if $\mathbf{jobvl}=\text{'V'}$, $\mathit{ldvl}=\max \left(1,\mathbf{n}\right)$;
• otherwise $\mathit{ldvl}=1$.

The second dimension of the array vl will be $\max \left(1,\mathbf{n}\right)$ if $\mathbf{jobvl}=\text{'V'}$ and $1$ otherwise.

If $\mathbf{jobvl}=\text{'V'}$, the left generalized eigenvectors $u_j$ are stored one after another in the columns of vl, in the same order as the corresponding eigenvalues. Each eigenvector will be scaled so the largest component will have $\left|\text{real part}\right|+\left|\text{imag. part}\right|=1$.

If $\mathbf{jobvl}=\text{'N'}$, vl is not referenced.

7:     $\mathrm{vr}\left(\mathit{ldvr},:\right)$ – double array

The first dimension, $\mathit{ldvr}$, of the array vr will be

• if $\mathbf{jobvr}=\text{'V'}$, $\mathit{ldvr}=\max \left(1,\mathbf{n}\right)$;
• otherwise $\mathit{ldvr}=1$.

The second dimension of the array vr will be $\max \left(1,\mathbf{n}\right)$ if $\mathbf{jobvr}=\text{'V'}$ and $1$ otherwise.

If $\mathbf{jobvr}=\text{'V'}$, the right generalized eigenvectors $v_j$ are stored one after another in the columns of vr, in the same order as the corresponding eigenvalues. Each eigenvector will be scaled so the largest component will have $\left|\text{real part}\right|+\left|\text{imag. part}\right|=1$.

If $\mathbf{jobvr}=\text{'N'}$, vr is not referenced.

8:     $\mathrm{ilo}$ – int64int32nag_int scalar

9:     $\mathrm{ihi}$ – int64int32nag_int scalar

ilo and ihi are integer values such that $\mathbf{a}\left(i,j\right)=0$ and $\mathbf{b}\left(i,j\right)=0$ if $i>j$ and $j=1,2,\dots ,\mathbf{ilo}-1$ or $i=\mathbf{ihi}+1,\dots ,\mathbf{n}$.

If $\mathbf{balanc}=\text{'N'}$ or $\text{'S'}$, $\mathbf{ilo}=1$ and $\mathbf{ihi}=\mathbf{n}$.

10:   $\mathrm{lscale}\left(\mathbf{n}\right)$ – double array

Details of the permutations and scaling factors applied to the left side of $A$ and $B$.

If ${\mathit{pl}}_j$ is the index of the row interchanged with row $j$, and ${\mathit{dl}}_j$ is the scaling factor applied to row $j$, then:

• $\mathbf{lscale}\left(\mathit{j}\right)={\mathit{pl}}_{\mathit{j}}$, for $\mathit{j}=1,2,\dots ,\mathbf{ilo}-1$;
• $\mathbf{lscale}={\mathit{dl}}_{\mathit{j}}$, for $\mathit{j}=\mathbf{ilo},\dots ,\mathbf{ihi}$;
• $\mathbf{lscale}={\mathit{pl}}_{\mathit{j}}$, for $\mathit{j}=\mathbf{ihi}+1,\dots ,\mathbf{n}$.

The order in which the interchanges are made is n to $\mathbf{ihi}+1$, then $1$ to $\mathbf{ilo}-1$.

11:   $\mathrm{rscale}\left(\mathbf{n}\right)$ – double array

Details of the permutations and scaling factors applied to the right side of $A$ and $B$.

If ${\mathit{pr}}_j$ is the index of the column interchanged with column $j$, and ${\mathit{dr}}_j$ is the scaling factor applied to column $j$, then:

• $\mathbf{rscale}\left(\mathit{j}\right)={\mathit{pr}}_{\mathit{j}}$, for $\mathit{j}=1,2,\dots ,\mathbf{ilo}-1$;
• if $\mathbf{rscale}={\mathit{dr}}_{\mathit{j}}$, for $\mathit{j}=\mathbf{ilo},\dots ,\mathbf{ihi}$;
• if $\mathbf{rscale}={\mathit{pr}}_{\mathit{j}}$, for $\mathit{j}=\mathbf{ihi}+1,\dots ,\mathbf{n}$.

The order in which the interchanges are made is n to $\mathbf{ihi}+1$, then $1$ to $\mathbf{ilo}-1$.

12:   $\mathrm{abnrm}$ – double scalar

The $1$-norm of the balanced matrix $A$.

13:   $\mathrm{bbnrm}$ – double scalar

The $1$-norm of the balanced matrix $B$.

14:   $\mathrm{rconde}\left(:\right)$ – double array

The dimension of the array rconde will be $\max \left(1,\mathbf{n}\right)$

If $\mathbf{sense}=\text{'E'}$ or $\text{'B'}$, the reciprocal condition numbers of the eigenvalues, stored in consecutive elements of the array. For a complex conjugate pair of eigenvalues two consecutive elements of rconde are set to the same value. Thus $\mathbf{rconde}\left(j\right)$, $\mathbf{rcondv}\left(j\right)$, and the $j$th columns of vl and vr all correspond to the $j$th eigenpair.

If $\mathbf{sense}=\text{'V'}$, rconde is not referenced.

15:   $\mathrm{rcondv}\left(:\right)$ – double array

The dimension of the array rcondv will be $\max \left(1,\mathbf{n}\right)$

If $\mathbf{sense}=\text{'V'}$ or $\text{'B'}$, the estimated reciprocal condition numbers of the eigenvectors, stored in consecutive elements of the array. For a complex eigenvector two consecutive elements of rcondv are set to the same value.

If $\mathbf{sense}=\text{'E'}$, rcondv is not referenced.

16:   $\mathrm{info}$ – int64int32nag_int scalar

$\mathbf{info}=0$ unless the function detects an error (see Error Indicators and Warnings).

Error Indicators and Warnings

Cases prefixed with W are classified as warnings and do not generate an error of type NAG:error_n. See nag_issue_warnings.

   $\mathbf{info}=-i$

If $\mathbf{info}=-i$, parameter $i$ had an illegal value on entry. The parameters are numbered as follows:

1: balanc, 2: jobvl, 3: jobvr, 4: sense, 5: n, 6: a, 7: lda, 8: b, 9: ldb, 10: alphar, 11: alphai, 12: beta, 13: vl, 14: ldvl, 15: vr, 16: ldvr, 17: ilo, 18: ihi, 19: lscale, 20: rscale, 21: abnrm, 22: bbnrm, 23: rconde, 24: rcondv, 25: work, 26: lwork, 27: iwork, 28: bwork, 29: info.

It is possible that info refers to a parameter that is omitted from the MATLAB interface. This usually indicates that an error in one of the other input parameters has caused an incorrect value to be inferred.

W  $\mathbf{info}=1\hspace{0.17em}\text{to}\hspace{0.17em}\mathbf{n}$

The $QZ$ iteration failed. No eigenvectors have been calculated, but $\mathbf{alphar}\left(j\right)$, $\mathbf{alphai}\left(j\right)$, and $\mathbf{beta}\left(j\right)$ should be correct for $j=\mathbf{info}+1,\dots ,\mathbf{n}$.

   $\mathbf{info}=\mathbf{n}+1$

Unexpected error returned from nag_lapack_dhgeqz (f08xe).

   $\mathbf{info}=\mathbf{n}+2$

Error returned from nag_lapack_dtgevc (f08yk).

Accuracy

Further Comments

Example

Open in the MATLAB editor: f08wb_example

function f08wb_example


fprintf('f08wb example results\n\n');

% Generalized eigenvalues and right eigenvectors of (A, B):
n = 4;
a = [3.9, 12.5, -34.5, -0.5;
     4.3, 21.5, -47.5,  7.5;
     4.3, 21.5, -43.5,  3.5;
     4.4, 26,   -46,    6];
b = [1,    2,    -3,    1;
     1,    3,    -5,    4;
     1,    3,    -4,    3;
     1,    3,    -4,    4];

balanc = 'Balance';
jobvl = 'No vectors (left)';
jobvr = 'Vectors (right)';
sense = 'Both reciprocal condition numbers';
[~, ~, alphar, alphai, beta, ~, VR, ~, ~, ~, ~, abnrm, bbnrm, ...
	     rconde, rcondv, info] = ...
  f08wb( ...
	 balanc, jobvl, jobvr, sense, a, b);

epsilon = x02aj;
small   = x02am;
absnrm  = sqrt(abnrm^2+bbnrm^2);
tol     = epsilon*absnrm;

for j=1:n

  % display information on the jth eigenvalue
  if (abs(alphar(j)) + abs(alphai(j)))*small >= abs(beta(j))
    fprintf('\n%4d: Eigenvalue is numerically infinite or undetermined\n',j);
    fprintf('%4d: alphar = %11.4e, alphai = %11.4e, beta = %11.4e\n', ...
	    j, alphar(j), alphai(j), beta(j));
  else
    fprintf('\n%4d: Eigenvalue =  ', j);
    if alphai(j)==0
      fprintf('%10.4e\n',alphar(j)/beta(j));
    elseif alphai(j)>=0
      fprintf('%10.4e + %10.4ei\n',alphar(j)/beta(j), alphai(j)/beta(j));
    else
      fprintf('%10.4e - %10.4ei\n',alphar(j)/beta(j), -alphai(j)/beta(j));
    end
  end

  if rconde(j) > 0
    fprintf('\n      Condition number            = %8.1e\n', 1/rconde(j));
    fprintf('      Error bound                 = %8.1e\n', tol/rconde(j));
  else
    fprintf('\n      Reciprocal condition number = %8.1e\n', rconde(j));
    fprintf('      Error bound is infinite\n');
  end

  fprintf('\n      Eigenvector:\n');
  if alphai(j) == 0
    fprintf('%30.4e\n',VR(:, j));
  else
    if alphai(j) > 0
      k = j;
    else   
      k = j-1;
    end
    for l = 1:n
      if VR(l,k+1)>=0
        fprintf('%30.4e + %10.4ei\n',VR(l, k), VR(l, k+1));
      else
        fprintf('%30.4e - %10.4ei\n',VR(l, k), -VR(l, k+1));
      end
    end
  end

  if rcondv(j) > 0
    fprintf('\n      Condition number            = %8.1e\n', 1/rcondv(j));
    fprintf('      Error bound                 = %8.1e\n', tol/rcondv(j));
  else
    fprintf('\n      Reciprocal condition number = %8.1e\n', rcondv(j));
    fprintf('      Error bound is infinite\n');
  end
end

f08wb example results


   1: Eigenvalue =  2.0000e+00

      Condition number            =  1.1e+01
      Error bound                 =  2.5e-14

      Eigenvector:
                   -1.0000e+00
                   -5.7143e-03
                   -6.2857e-02
                   -6.2857e-02

      Condition number            =  8.0e+00
      Error bound                 =  1.9e-14

   2: Eigenvalue =  3.0000e+00 + 4.0000e+00i

      Condition number            =  6.1e+00
      Error bound                 =  1.4e-14

      Eigenvector:
                   -4.2550e-01 - 5.7450e-01i
                   -8.5099e-02 - 1.1490e-01i
                   -1.4298e-01 - 8.6125e-04i
                   -1.4298e-01 - 8.6125e-04i

      Condition number            =  2.6e+01
      Error bound                 =  6.2e-14

   3: Eigenvalue =  3.0000e+00 - 4.0000e+00i

      Condition number            =  6.1e+00
      Error bound                 =  1.4e-14

      Eigenvector:
                   -4.2550e-01 - 5.7450e-01i
                   -8.5099e-02 - 1.1490e-01i
                   -1.4298e-01 - 8.6125e-04i
                   -1.4298e-01 - 8.6125e-04i

      Condition number            =  2.6e+01
      Error bound                 =  6.2e-14

   4: Eigenvalue =  4.0000e+00

      Condition number            =  1.9e+00
      Error bound                 =  4.6e-15

      Eigenvector:
                   -1.0000e+00
                   -1.1111e-02
                    3.3333e-02
                   -1.5556e-01

      Condition number            =  1.4e+01
      Error bound                 =  3.3e-14