F07 Chapter Contents
F07 Chapter Introduction
NAG Library Manual

# NAG Library Routine DocumentF07MDF (DSYTRF)

Note:  before using this routine, please read the Users' Note for your implementation to check the interpretation of bold italicised terms and other implementation-dependent details.

## 1  Purpose

F07MDF (DSYTRF) computes the Bunch–Kaufman factorization of a real symmetric indefinite matrix.

## 2  Specification

 SUBROUTINE F07MDF ( UPLO, N, A, LDA, IPIV, WORK, LWORK, INFO)
 INTEGER N, LDA, IPIV(*), LWORK, INFO REAL (KIND=nag_wp) A(LDA,*), WORK(max(1,LWORK)) CHARACTER(1) UPLO
The routine may be called by its LAPACK name dsytrf.

## 3  Description

F07MDF (DSYTRF) factorizes a real symmetric matrix $A$, using the Bunch–Kaufman diagonal pivoting method. $A$ is factorized as either $A=PUD{U}^{\mathrm{T}}{P}^{\mathrm{T}}$ if ${\mathbf{UPLO}}=\text{'U'}$ or $A=PLD{L}^{\mathrm{T}}{P}^{\mathrm{T}}$ if ${\mathbf{UPLO}}=\text{'L'}$, where $P$ is a permutation matrix, $U$ (or $L$) is a unit upper (or lower) triangular matrix and $D$ is a symmetric block diagonal matrix with $1$ by $1$ and $2$ by $2$ diagonal blocks; $U$ (or $L$) has $2$ by $2$ unit diagonal blocks corresponding to the $2$ by $2$ blocks of $D$. Row and column interchanges are performed to ensure numerical stability while preserving symmetry.
This method is suitable for symmetric matrices which are not known to be positive definite. If $A$ is in fact positive definite, no interchanges are performed and no $2$ by $2$ blocks occur in $D$.

## 4  References

Golub G H and Van Loan C F (1996) Matrix Computations (3rd Edition) Johns Hopkins University Press, Baltimore

## 5  Parameters

1:     UPLO – CHARACTER(1)Input
On entry: specifies whether the upper or lower triangular part of $A$ is stored and how $A$ is to be factorized.
${\mathbf{UPLO}}=\text{'U'}$
The upper triangular part of $A$ is stored and $A$ is factorized as $PUD{U}^{\mathrm{T}}{P}^{\mathrm{T}}$, where $U$ is upper triangular.
${\mathbf{UPLO}}=\text{'L'}$
The lower triangular part of $A$ is stored and $A$ is factorized as $PLD{L}^{\mathrm{T}}{P}^{\mathrm{T}}$, where $L$ is lower triangular.
Constraint: ${\mathbf{UPLO}}=\text{'U'}$ or $\text{'L'}$.
2:     N – INTEGERInput
On entry: $n$, the order of the matrix $A$.
Constraint: ${\mathbf{N}}\ge 0$.
3:     A(LDA,$*$) – REAL (KIND=nag_wp) arrayInput/Output
Note: the second dimension of the array A must be at least $\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left(1,{\mathbf{N}}\right)$.
On entry: the $n$ by $n$ symmetric indefinite matrix $A$.
• If ${\mathbf{UPLO}}=\text{'U'}$, the upper triangular part of $A$ must be stored and the elements of the array below the diagonal are not referenced.
• If ${\mathbf{UPLO}}=\text{'L'}$, the lower triangular part of $A$ must be stored and the elements of the array above the diagonal are not referenced.
On exit: the upper or lower triangle of $A$ is overwritten by details of the block diagonal matrix $D$ and the multipliers used to obtain the factor $U$ or $L$ as specified by UPLO.
4:     LDA – INTEGERInput
On entry: the first dimension of the array A as declared in the (sub)program from which F07MDF (DSYTRF) is called.
Constraint: ${\mathbf{LDA}}\ge \mathrm{max}\phantom{\rule{0.125em}{0ex}}\left(1,{\mathbf{N}}\right)$.
5:     IPIV($*$) – INTEGER arrayOutput
Note: the dimension of the array IPIV must be at least $\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left(1,{\mathbf{N}}\right)$.
On exit: details of the interchanges and the block structure of $D$. More precisely,
• if ${\mathbf{IPIV}}\left(i\right)=k>0$, ${d}_{ii}$ is a $1$ by $1$ pivot block and the $i$th row and column of $A$ were interchanged with the $k$th row and column;
• if ${\mathbf{UPLO}}=\text{'U'}$ and ${\mathbf{IPIV}}\left(i-1\right)={\mathbf{IPIV}}\left(i\right)=-l<0$, $\left(\begin{array}{cc}{d}_{i-1,i-1}& {\stackrel{-}{d}}_{i,i-1}\\ {\stackrel{-}{d}}_{i,i-1}& {d}_{ii}\end{array}\right)$ is a $2$ by $2$ pivot block and the $\left(i-1\right)$th row and column of $A$ were interchanged with the $l$th row and column;
• if ${\mathbf{UPLO}}=\text{'L'}$ and ${\mathbf{IPIV}}\left(i\right)={\mathbf{IPIV}}\left(i+1\right)=-m<0$, $\left(\begin{array}{cc}{d}_{ii}& {d}_{i+1,i}\\ {d}_{i+1,i}& {d}_{i+1,i+1}\end{array}\right)$ is a $2$ by $2$ pivot block and the $\left(i+1\right)$th row and column of $A$ were interchanged with the $m$th row and column.
6:     WORK($\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left(1,{\mathbf{LWORK}}\right)$) – REAL (KIND=nag_wp) arrayWorkspace
On exit: if ${\mathbf{INFO}}=0$, ${\mathbf{WORK}}\left(1\right)$ contains the minimum value of LWORK required for optimum performance.
7:     LWORK – INTEGERInput
On entry: the dimension of the array WORK as declared in the (sub)program from which F07MDF (DSYTRF) is called, unless ${\mathbf{LWORK}}=-1$, in which case a workspace query is assumed and the routine only calculates the optimal dimension of WORK (using the formula given below).
Suggested value: for optimum performance LWORK should be at least ${\mathbf{N}}×\mathit{nb}$, where $\mathit{nb}$ is the block size.
Constraint: ${\mathbf{LWORK}}\ge 1$ or ${\mathbf{LWORK}}=-1$.
8:     INFO – INTEGEROutput
On exit: ${\mathbf{INFO}}=0$ unless the routine detects an error (see Section 6).

## 6  Error Indicators and Warnings

Errors or warnings detected by the routine:
${\mathbf{INFO}}<0$
If ${\mathbf{INFO}}=-i$, the $i$th parameter had an illegal value. An explanatory message is output, and execution of the program is terminated.
${\mathbf{INFO}}>0$
If ${\mathbf{INFO}}=i$, $d\left(i,i\right)$ is exactly zero. The factorization has been completed, but the block diagonal matrix $D$ is exactly singular, and division by zero will occur if it is used to solve a system of equations.

## 7  Accuracy

If ${\mathbf{UPLO}}=\text{'U'}$, the computed factors $U$ and $D$ are the exact factors of a perturbed matrix $A+E$, where
 $E≤cnεPUDUTPT ,$
$c\left(n\right)$ is a modest linear function of $n$, and $\epsilon$ is the machine precision.
If ${\mathbf{UPLO}}=\text{'L'}$, a similar statement holds for the computed factors $L$ and $D$.

## 8  Further Comments

The elements of $D$ overwrite the corresponding elements of $A$; if $D$ has $2$ by $2$ blocks, only the upper or lower triangle is stored, as specified by UPLO.
The unit diagonal elements of $U$ or $L$ and the $2$ by $2$ unit diagonal blocks are not stored. The remaining elements of $U$ or $L$ are stored in the corresponding columns of the array A, but additional row interchanges must be applied to recover $U$ or $L$ explicitly (this is seldom necessary). If ${\mathbf{IPIV}}\left(\mathit{i}\right)=\mathit{i}$, for $\mathit{i}=1,2,\dots ,n$ (as is the case when $A$ is positive definite), then $U$ or $L$ is stored explicitly (except for its unit diagonal elements which are equal to $1$).
The total number of floating point operations is approximately $\frac{1}{3}{n}^{3}$.
A call to F07MDF (DSYTRF) may be followed by calls to the routines:
• F07MEF (DSYTRS) to solve $AX=B$;
• F07MGF (DSYCON) to estimate the condition number of $A$;
• F07MJF (DSYTRI) to compute the inverse of $A$.
The complex analogues of this routine are F07MRF (ZHETRF) for Hermitian matrices and F07NRF (ZSYTRF) for symmetric matrices.

## 9  Example

This example computes the Bunch–Kaufman factorization of the matrix $A$, where
 $A= 2.07 3.87 4.20 -1.15 3.87 -0.21 1.87 0.63 4.20 1.87 1.15 2.06 -1.15 0.63 2.06 -1.81 .$

### 9.1  Program Text

Program Text (f07mdfe.f90)

### 9.2  Program Data

Program Data (f07mdfe.d)

### 9.3  Program Results

Program Results (f07mdfe.r)