# NAG FL InterfaceReplacement Calls Advice

The following list gives the names of routines that are suitable replacements for routines that have either been withdrawn or deprecated since Mark 20.
The list indicates the minimum change necessary, but many of the replacement routines have additional flexibility and you may wish to take advantage of new features. It is strongly recommended that you consult the routine documents.

## C05 – Roots of One or More Transcendental Equations

Withdrawn at Mark 25.
Replaced by c05ayf.

#### Old Code

```FUNCTION F(XX)
...
END FUNCTION F
...

#### New Code

```FUNCTION F(XX,IUSER,RUSER)
...
INTEGER,            INTENT(INOUT) :: IUSER(*)
REAL (KIND=nag_wp), INTENT(INOUT) :: RUSER(*)
...
END FUNCTION F
...
INTEGER            :: IUSER(1)
REAL (KIND=nag_wp) :: RUSER(1)
...
CALL c05ayf(A,B,EPS,ETA,F,X,IUSER,RUSER,IFAIL)```

### c05agf

Withdrawn at Mark 25.
Replaced by c05auf.

#### Old Code

```FUNCTION F(XX)
...
END FUNCTION F
...
CALL c05agf(X,H,EPS,ETA,F,A,B,IFAIL)```

#### New Code

```FUNCTION F(XX,IUSER,RUSER)
...
INTEGER,            INTENT(INOUT) :: IUSER(*)
REAL (KIND=nag_wp), INTENT(INOUT) :: RUSER(*)
...
END FUNCTION F
...
INTEGER            :: IUSER(1)
REAL (KIND=nag_wp) :: RUSER(1)
...
CALL c05auf(X,H,EPS,ETA,F,A,B,IUSER,RUSER,IFAIL)```

### c05ajf

Withdrawn at Mark 25.
Replaced by c05awf.

#### Old Code

```FUNCTION F(XX)
...
END FUNCTION F
...
CALL c05ajf(X,EPS,ETA,F,NFMAX,IFAIL)```

#### New Code

```FUNCTION F(XX,IUSER,RUSER)
...
INTEGER,            INTENT(INOUT) :: IUSER(*)
REAL (KIND=nag_wp), INTENT(INOUT) :: RUSER(*)
...
END FUNCTION F
...
INTEGER            :: IUSER(1)
REAL (KIND=nag_wp) :: RUSER(1)
...
CALL c05awf(X,EPS,ETA,F,NFMAX,IUSER,RUSER,IFAIL)```

### c05nbf

Withdrawn at Mark 25.
Replaced by c05qbf.

#### Old Code

```SUBROUTINE FCN(N,X,FVEC,IFLAG)
...
END SUBROUTINE FCN
...
CALL c05nbf(FCN,N,X,FVEC,XTOL,WA,LWA,IFAIL)```

#### New Code

```SUBROUTINE FCN(N,X,FVEC,IUSER,RUSER,IFLAG)
...
INTEGER,            INTENT(INOUT) :: IUSER(*)
REAL (KIND=nag_wp), INTENT(INOUT) :: RUSER(*)
...
END FUNCTION FCN
...
INTEGER            :: IUSER(1)
REAL (KIND=nag_wp) :: RUSER(1)
...
CALL c05qbf(FCN,N,X,FVEC,XTOL,IUSER,RUSER,IFAIL)```

### c05ncf

Withdrawn at Mark 25.
Replaced by c05qcf.

#### Old Code

```SUBROUTINE FCN(N,X,FVEC,IFLAG)
...
END SUBROUTINE FCN
...
REAL (KIND=nag_wp) :: FJAC(LDFJAC,N)
...
CALL c05ncf(FCN,N,X,FVEC,XTOL,MAXFEV,ML,MU,EPSFCN,DIAG,MODE,FACTOR, &
NPRINT,NFEV,FJAC,LDFJAC,R,LR,QTF,W,IFAIL)```

#### New Code

```SUBROUTINE FCN(N,X,FVEC,IUSER,RUSER,IFLAG)
...
INTEGER,            INTENT(INOUT) :: IUSER(*)
REAL (KIND=nag_wp), INTENT(INOUT) :: RUSER(*)
...
END FUNCTION FCN
...
INTEGER            :: IUSER(1)
REAL (KIND=nag_wp) :: FJAC(N,N), RUSER(1)
...
CALL c05qcf(FCN,N,X,FVEC,XTOL,MAXFEV,ML,MU,EPSFCN,MODE,DIAG,FACTOR, &
NPRINT,NFEV,FJAC,R,QTF,IUSER,RUSER,IFAIL)```

### c05ndf

Withdrawn at Mark 25.
Replaced by c05qdf.

#### Old Code

```REAL (KIND=nag_wp) :: FJAC(LDFJAC,N)
...
CALL c05ndf(IREVCM,N,X,FVEC,XTOL,ML,MU,EPSFCN,DIAG,MODE,FACTOR, &
FJAC,LDFJAC,R,LR,QTF,W,IFAIL)```

#### New Code

```REAL (KIND=nag_wp) :: FJAC(N,N), RWSAV(4*N+20)
INTEGER            :: IWSAV(17)
...
CALL c05qdf(IREVCM,N,X,FVEC,XTOL,ML,MU,EPSFCN,MODE,DIAG,FACTOR, &
FJAC,R,QTF,IWSAV,RWSAV,IFAIL)```

### c05pbf/​c05pba

Withdrawn at Mark 25.
Replaced by c05rbf.

#### Old Code

```SUBROUTINE FCN_C05PBF(N,X,FVEC,FJAC,LDFJAC,IFLAG)
...
END SUBROUTINE FCN_C05PBF
...
REAL (KIND=nag_wp) :: FJAC(LDFJAC,N)
...
CALL C05PBF(FCN_C05PBF,N,X,FVEC,FJAC,LDFJAC,XTOL,WA,LWA,IFAIL)
or
SUBROUTINE FCN_C05PBA(N,X,FVEC,FJAC,LDFJAC,IFLAG,IUSER,RUSER)
...
END SUBROUTINE FCN_C05PBA
...
REAL (KIND=nag_wp) :: FJAC(LDFJAC,N)
...
CALL C05PBA(FCN_C05PBA,N,X,FVEC,FJAC,LDFJAC,XTOL,WA,LWA,IUSER,RUSER,IFAIL)```

#### New Code

```SUBROUTINE FCN(N,X,FVEC,FJAC,IUSER,RUSER,IFLAG)
...
END SUBROUTINE FCN
...
REAL (KIND=nag_wp) :: FJAC(N,N)
...
CALL c05rbf(FCN,N,X,FVEC,FJAC,XTOL,IUSER,RUSER,IFAIL)```

### c05pcf/​c05pca

Withdrawn at Mark 25.
Replaced by c05rcf.

#### Old Code

```SUBROUTINE FCN_C05PCF(N,X,FVEC,FJAC,LDFJAC,IFLAG)
...
END SUBROUTINE FCN_C05PCF
...
REAL (KIND=nag_wp) :: FJAC(LDFJAC,N)
...
CALL C05PCF(FCN_C05PCF,N,X,FVEC,FJAC,LDFJAC,XTOL,MAXFEV,DIAG,MODE,FACTOR, &
NPRINT,NFEV,NJEV,R,LR,QTF,W,IFAIL)
or
SUBROUTINE FCN_C05PCA(N,X,FVEC,FJAC,LDFJAC,IFLAG,IUSER,RUSER)
...
END SUBROUTINE FCN_C05PCA
...
REAL (KIND=nag_wp) :: FJAC(LDFJAC,N)
...
CALL C05PCA(FCN_C05PCA,N,X,FVEC,FJAC,LDFJAC,XTOL,MAXFEV,DIAG,MODE,FACTOR, &
NPRINT,NFEV,NJEV,R,LR,QTF,W,IUSER,RUSER,IFAIL)```

#### New Code

```SUBROUTINE FCN(N,X,FVEC,FJAC,IUSER,RUSER,IFLAG)
...
INTEGER,            INTENT(INOUT) :: IUSER(*)
REAL (KIND=nag_wp), INTENT(INOUT) :: RUSER(*)
...
END FUNCTION FCN
...
REAL (KIND=nag_wp) :: FJAC(N,N)
...
CALL c05rcf(FCN,N,X,FVEC,FJAC,XTOL,MAXFEV,MODE,DIAG,FACTOR, &
NPRINT,NFEV,NJEV,R,QTF,IUSER,RUSER,IFAIL)```

### c05pdf/​c05pda

Withdrawn at Mark 25.
Replaced by c05rdf.

#### Old Code

```REAL (KIND=nag_wp) :: FJAC(LDFJAC,N), RWSAV(10)
INTEGER            :: IWSAV(15)
...
CALL C05PDF(IREVCM,N,X,FVEC,FJAC,LDFJAC,XTOL,DIAG,MODE,FACTOR, &
R,LR,QTF,W,IFAIL)
or
CALL C05PDA(IREVCM,N,X,FVEC,FJAC,LDFJAC,XTOL,DIAG,MODE,FACTOR, &
R,LR,QTF,W,LWSAV,IWSAV,RWSAV,IFAIL)```

#### New Code

```REAL (KIND=nag_wp) :: FJAC(N,N), RWSAV(4*N+10)
INTEGER            :: IWSAV(17)
...
CALL c05rdf(IREVCM,N,X,FVEC,FJAC,XTOL,MODE,DIAG,FACTOR,        &
R,QTF,IWSAV,RWSAV,IFAIL)```

### c05zaf

Withdrawn at Mark 25.
Replaced by c05zdf.

#### Old Code

`CALL c05zaf(M,N,X,FVEC,FJAC,LDFJAC,XP,FVECP,MODE,ERR)`

#### New Code

```IFAIL = 0
CALL c05zdf(MODE,M,N,X,FVEC,FJAC,LDFJAC,XP,FVECP,ERR,IFAIL)```
The array xp must now have dimension n regardless of the value of mode, and likewise err must now have dimension m regardless. The argument ifail is the standard NAG argument for error trapping. If you are unfamiliar with this argument you should refer to Section 4 in the Introduction to the NAG Library FL Interface for details.

## C06 – Summation of Series

### c06dbf

Withdrawn at Mark 25.
Replaced by c06dcf.

#### Old Code

```DO I = 1, LX
RES(I) = c06dbf(X(I),C,N,S)
END DO```

#### New Code

```XMIN = -1.0D0
XMAX = 1.0D0

SELECT CASE (S)
CASE (1,2,3)
S_USE = S
CASE DEFAULT
S_USE = 2
END SELECT

IFAIL = 0
CALL c06dcf(X,LX,XMIN,XMAX,C,N,S_USE,RES,IFAIL)```
The old routine c06dbf returns a single sum at a time, whereas the new routine c06dcf returns a vector of lx values at once. The values supplied in x to c06dcf are un-normalized original variable values in the range $\left[{\mathbf{xmin}},{\mathbf{xmax}}\right]$. The argument ifail is the standard NAG argument for error trapping. If you are unfamiliar with this argument you should refer to Section 4 in the Introduction to the NAG Library FL Interface for details.

### c06eaf

Withdrawn at Mark 26.
Replaced by c06paf.
c06paf removes restrictions on sequence length and combines transform directions.

#### Old Code

`CALL c06eaf(X,N,IFAIL)`

#### New Code

`CALL c06paf('F',X,N,WORK,IFAIL)`
where work is a real array of length $3×{\mathbf{n}}+100$ and the dimension of the array x has been extended from the original n to ${\mathbf{n}}+2$. The output values x are stored in a different order with real and imaginary parts stored contiguously. The mapping of output elements is as follows:
• ${\mathbf{x}}\left(2×\mathit{i}\right)←{\mathbf{x}}\left(\mathit{i}\right)$, for $\mathit{i}=0,1,\dots ,{\mathbf{n}}/2$ and
${\mathbf{x}}\left(2×\mathit{i}+1\right)←{\mathbf{x}}\left({\mathbf{n}}-\mathit{i}\right)$, for $\mathit{i}=1,2,\dots ,\left({\mathbf{n}}+1\right)/2$.

### c06ebf

Withdrawn at Mark 26.
Replaced by c06paf.
c06paf removes restrictions on sequence length and combines transform directions.

#### Old Code

`CALL c06ebf(X,N,IFAIL)`

#### New Code

`CALL c06paf('B',X,N,WORK,IFAIL)`
where work is a real array of length $3×{\mathbf{n}}+100$ and the dimension of the array x has been extended from the original n to ${\mathbf{n}}+2$. The input values of x are stored in a different order with real and imaginary parts stored contiguously. Also c06paf performs the inverse transform without the need to first conjugate. If prior conjugation of original array x is assumed then the mapping of input elements is:
• ${\mathbf{x}}\left(2×\mathit{i}\right)←{\mathbf{x}}\left(\mathit{i}\right)$, for $\mathit{i}=0,1,\dots ,{\mathbf{n}}/2$ and
${\mathbf{x}}\left(2×\mathit{i}+1\right)←{\mathbf{x}}\left({\mathbf{n}}-\mathit{i}\right)$, for $\mathit{i}=1,2,\dots ,\left({\mathbf{n}}-1\right)/2$.

### c06ecf

Withdrawn at Mark 26.
Replaced by c06pcf.
c06pcf removes restrictions on sequence length, combines transform directions and uses complex types.

#### Old Code

`CALL c06ecf(X,Y,N,IFAIL)`

#### New Code

`CALL c06pcf('F',Z,N,WORK,IFAIL)`
where work is a complex array of length $2×{\mathbf{n}}+15$ and Z is a complex array of length n such that $\text{Z}\left(\mathit{i}\right)=\mathrm{CMPLX}\left(\text{X}\left(\mathit{i}\right),\text{Y}\left(\mathit{i}\right)\right)$, for $\mathit{i}=0,1,\dots ,{\mathbf{n}}-1$ on input and output.

### c06ekf

Withdrawn at Mark 26.
Replaced by c06fkf.
c06fkf removes restrictions on sequence length.

#### Old Code

`CALL c06ekf(IJOB,X,Y,N,IFAIL)`

#### New Code

`CALL c06fkf(IJOB,X,Y,N,WORK,IFAIL)`
where work is a real array of length n.

### c06fpf

Scheduled for withdrawal at Mark 28.
Replaced by c06ppf or c06pqf.
c06pqf provides a simpler interface for both forward and backward transforms. c06ppf retains original input ordering at the expense of efficiency.

#### Old Code

`CALL c06fpf(M,N,X,INIT,TRIG,WORK,IFAIL)`

#### New Code

`CALL c06pqf('F',N,M,X,WORK,IFAIL)`
where the dimension of work has been extended from ${\mathbf{m}}×{\mathbf{n}}$ to ${\mathbf{m}}×{\mathbf{n}}+2×{\mathbf{n}}$ (to include TRIG) and the dimension of the array x has been extended from the original ${\mathbf{n}}×{\mathbf{m}}$ to $\left({\mathbf{n}}+2\right)×{\mathbf{m}}$.
The input values are stored slightly differently to allow for two extra storage spaces at the end of each sequence.
The mapping of input elements is as follows:
• (Here x begins at element zero ${\mathbf{x}}\left(0\right)$.)
for $\mathit{j}=0,1,\dots ,{\mathbf{m}}-1$
• $J=j×\left({\mathbf{n}}+2\right)$;
• ${\mathbf{x}}\left(J+\mathit{i}\right)←{\mathbf{x}}\left(\mathit{i}×{\mathbf{m}}+j\right)$, for $\mathit{i}=0,1,\dots ,{\mathbf{n}}-1$;
• ${\mathbf{x}}\left(J+{\mathbf{n}}\right)$ and ${\mathbf{x}}\left(J+{\mathbf{n}}+1\right)$ need not be set.
The output values x are stored in a different order with real and imaginary parts of each Hermitian sequence stored contiguously.
The mapping of output elements is as follows:
• (Here x begins at element zero ${\mathbf{x}}\left(0\right)$.)
for $\mathit{j}=0,1,\dots ,{\mathbf{m}}-1$
• $J=j×\left({\mathbf{n}}+2\right)$;
• ${\mathbf{x}}\left(J+2×\mathit{i}\right)←{\mathbf{x}}\left(\mathit{i}×{\mathbf{m}}+j\right)$, for $\mathit{i}=0,1,\dots ,{\mathbf{n}}/2$ [real parts];
• ${\mathbf{x}}\left(J+2×\mathit{i}+1\right)←{\mathbf{x}}\left(\left({\mathbf{n}}-\mathit{i}\right)×{\mathbf{m}}+j\right)$, for $\mathit{i}=1,2,\dots ,\left({\mathbf{n}}-1\right)/2$ [imaginary parts];
• ${\mathbf{x}}\left(J+1\right)$ is set to zero;
• ${\mathbf{x}}\left(J+{\mathbf{n}}+1\right)$ is set to zero when n is even.

### c06fqf

Scheduled for withdrawal at Mark 28.
Replaced by c06ppf or c06pqf.
c06pqf provides a simpler interface for both forward and backward transforms. c06ppf returns real sequences in row order at the expense of efficiency.

#### Old Code

`CALL c06fqf(M,N,X,INIT,TRIG,WORK,IFAIL)`

#### New Code

`CALL c06pqf('B',N,M,X,WORK,IFAIL)`
where the dimension of work has been extended from ${\mathbf{m}}×{\mathbf{n}}$ to ${\mathbf{m}}×{\mathbf{n}}+2×{\mathbf{n}}$ (to include TRIG) and the dimension of the array x has been extended from the original ${\mathbf{n}}×{\mathbf{m}}$ to $\left({\mathbf{n}}+2\right)×{\mathbf{m}}$.
The input values x are stored in a different order with real and imaginary parts of each Hermitian sequence stored contiguously.
The mapping of input elements is as follows:
• (Here x begins at element zero ${\mathbf{x}}\left(0\right)$.)
for $\mathit{j}=0,1,\dots ,{\mathbf{m}}-1$
• $J=j×\left({\mathbf{n}}+2\right)$;
• ${\mathbf{x}}\left(J+2×\mathit{i}\right)←{\mathbf{x}}\left(\mathit{i}×{\mathbf{m}}+j\right)$, for $\mathit{i}=0,1,\dots ,{\mathbf{n}}/2$ [real parts];
• ${\mathbf{x}}\left(J+2×\mathit{i}+1\right)←{\mathbf{x}}\left(\left({\mathbf{n}}-\mathit{i}\right)×{\mathbf{m}}+j\right)$, for $\mathit{i}=1,2,\dots ,\left({\mathbf{n}}-1\right)/2$ [imaginary parts];
• ${\mathbf{x}}\left(J+1\right)$ must be zero;
• ${\mathbf{x}}\left(J+{\mathbf{n}}+1\right)$ must zero when n is even.
The output values are stored slightly differently to allow for two extra storage spaces at the end of each sequence.
The mapping of output elements is as follows:
• (Here x begins at element zero ${\mathbf{x}}\left(0\right)$.)
for $\mathit{j}=0,1,\dots ,{\mathbf{m}}-1$
• $J=j×\left({\mathbf{n}}+2\right)$;
• ${\mathbf{x}}\left(J+\mathit{i}\right)←{\mathbf{x}}\left(\mathit{i}×{\mathbf{m}}+j\right)$, for $\mathit{i}=0,1,\dots ,{\mathbf{n}}-1$;
• ${\mathbf{x}}\left(J+{\mathbf{n}}\right)$ and ${\mathbf{x}}\left(J+{\mathbf{n}}+1\right)$ will be set to zero.

### c06frf

Withdrawn at Mark 26.
Replaced by c06psf.
c06psf provides a simpler interface for both forward and backward transforms.

#### Old Code

`call c06frf(M,N,X,Y,INIT,TRIG,WORK,IFAIL)`

#### New Code

```Do j = 1, m*n
cx(j) = cmplx(x(j),y(j),kind=nag_wp)
End Do
Call c06psf('F',M,N,CX,CWORK,IFAIL)
x(1:m*n) = real(cx(1:m*n))
y(1:m*n) = aimag(cx(1:m*n))```
where $\mathrm{cx}$ and $\mathrm{cwork}$ are complex array of length $\mathrm{m}×\mathrm{n}$ and $\mathrm{n}×\mathrm{m}+2×\mathrm{n}+15$ respectively.

### c06fuf

Withdrawn at Mark 26.
Replaced by c06puf.
c06puf provides a simpler interface for both forward and backward transforms.

#### Old Code

`Call c06fuf(M,N,X,Y,INIT,TRIGM,TRIGN,WORK,IFAIL)`

#### New Code

```Do j = 1, m*n
cx(j) = cmplx(x(j),y(j),kind=nag_wp)
End Do
Call c06puf('F',M,N,CX,CWORK,IFAIL)
x(1:m*n) = real(cx(1:m*n))
y(1:m*n) = aimag(cx(1:m*n))```
where $\mathrm{cx}$ and $\mathrm{cwork}$ are complex arrays of lengths $\mathrm{m}×\mathrm{n}$ and $\mathrm{n}×\mathrm{m}+2×\mathrm{n}+2×\mathrm{m}+30$ respectively.

### c06gbf

Withdrawn at Mark 26.
There is no replacement for this routine.

### c06gcf

Withdrawn at Mark 26.
There is no replacement for this routine.

### c06gqf

Withdrawn at Mark 26.
There is no replacement for this routine.

### c06gsf

Withdrawn at Mark 26.
There is no replacement for this routine.

### c06haf

Withdrawn at Mark 26.
Replaced by c06ref.
c06ref has a simpler interface, storing sequences by column.

#### Old Code

`Call c06haf(M,N,X,INIT,TRIG,WORK,IFAIL)`

#### New Code

`Call c06ref(M,N,Y,IFAIL)`
where $\mathit{y}\left(1:n-1,m\right)$ is a two-dimensional real array such that $\mathit{y}\left(1:n-1,j\right)=\mathit{x}\left(j:m×\left(n-1\right),m\right)$.

### c06hbf

Withdrawn at Mark 26.
Replaced by c06rff.
c06rff has a simpler interface, storing sequences by column.

#### Old Code

`Call c06hbf(M,N,X,INIT,TRIG,WORK,IFAIL)`

#### New Code

`Call c06rff(M,N,Y,IFAIL)`
where $\mathit{y}\left(0:n,m\right)$ is a two-dimensional real array such that $\mathit{y}\left(0:n,j\right)=\mathit{x}\left(j:m×\left(n+1\right),m\right)$.

### c06hcf

Withdrawn at Mark 26.
Replaced by c06rgf.
c06rgf has a simpler interface, storing sequences by column.

#### Old Code

`Call c06hcf(DIRECT,M,N,X,INIT,TRIG,WORK,IFAIL)`

#### New Code

`Call c06rgf(IDIR,M,N,Y,IFAIL)`
where $\mathit{y}\left(1:n,m\right)$ is a two-dimensional real array such that $\mathit{y}\left(1:n,j\right)=\mathit{x}\left(j:m×n,m\right)$; ${\mathbf{idir}}=1$ or $-1$ for forward and inverse transforms respectively.

### c06hdf

Withdrawn at Mark 26.
Replaced by c06rhf.
c06rhf has a simpler interface, storing sequences by column.

#### Old Code

`Call c06hdf(DIRECT,M,N,X,INIT,TRIG,WORK,IFAIL)`

#### New Code

`Call c06rhf(IDIR,M,N,Y,IFAIL)`
where $\mathit{y}\left(0:n-1,m\right)$ is a two-dimensional real array such that $\mathit{y}\left(0:n-1,j\right)=\mathit{x}\left(j:m×n,m\right)$; ${\mathbf{idir}}=1$ or $-1$ for forward and inverse transforms respectively.

### d01baf

Withdrawn at Mark 26.
Replaced by d01uaf.
Withdrawn to provide thread safety in passing of data to user supplied function and a simpler interface to select the quadrature rule.

#### Old Code

```FUNCTION FUN(x)

...
real(kind=nag_wp) :: FUN
real(kind=nag_wp), intent(in) :: X
FUN = ...
END FUNCTION

DINEST = d01baf(D01XXX,A,B,N,FUN,IFAIL)```

#### New Code

```SUBROUTINE F(X,NX,FV,IFLAG,IUSER,RUSER)
...
! see example below
...
END SUBROUTINE F
...
integer :: key
integer, allocatable :: iuser(:)
real(kind=nag_wp), allocatable :: ruser(:)

! set KEY according to quadrature formula
! KEY =  0 : (D01XXX=D01BAZ)
! KEY = -3 : (D01XXX=D01BAY)
! KEY = -4 : (D01XXX=D01BAW)
! KEY = -5 : (D01XXX=D01BAX)
! KEY = ABS(KEY) for normal weights
KEY = 0

allocate(iuser(liuser), ruser(lruser))

CALL d01uaf(KEY,A,B,N,F,DINEST,IUSER,RUSER,IFAIL)```
iuser and ruser are arrays available to allow you to pass information to the user-supplied subroutine f.
iflag is an integer which you may use to force an immediate exit from d01uaf in case of an error in the user-supplied subroutine f.
f may be used to call the original fun as follows, although it may be more efficient to recode the integrand.
```     SUBROUTINE F(X,NX,FV,IFLAG,IUSER,RUSER)
...
integer, intent(in) :: NX
integer, intent(inout) :: iflag
real(kind=nag_wp), intent(in) :: X(NX)
real(kind=nag_wp), intent(out) :: fv(nx)
real(kind=wp), intent(inout) :: ruser(*)
integer, intent(inout) :: iuser(*)
integer :: j
external FUN

do j=1,nx
FV(j)  = FUN(x(j))
enddo

END SUBROUTINE F```

### d01bbf

Withdrawn at Mark 26.
Replaced by d01tbf.
Withdrawn to provide thread safety in passing of data to the user-supplied routine and a simpler interface to select the quadrature rule.

#### Old Code

`CALL d01bbf(D01XXX,A,B,ITYPE,N,WEIGHT,ABSCIS,IFAIL)`

#### New Code

```Integer :: key
CALL d01tbf(KEY,A,B,N,WEIGHT,ABSICS,IFAIL)```
The supplied subroutines D01XXX and the argument itype have been combined into a single argument key. ${\mathbf{key}}<0$ is equivalent to $\text{ITYPE}=1$ (adjusted weights). ${\mathbf{key}}>0$ is equivalent to $\mathbf{itype}=0$ (normal weights). $\left|{\mathbf{key}}\right|$ indicates the quadrature rule.
• $\left|{\mathbf{key}}\right|=0$ : Gauss–Legendre ($\mathbf{D01XXX}=\mathbf{D01BAZ}$)
• $\left|{\mathbf{key}}\right|=3$ : Gauss–Laguerre ($\mathbf{D01XXX}=\mathbf{D01BAX}$)
• $\left|{\mathbf{key}}\right|=4$ : Gauss–Hermite ($\mathbf{D01XXX}=\mathbf{D01BAW}$)
• $\left|{\mathbf{key}}\right|=5$ : Rational Gauss ($\mathbf{D01XXX}=\mathbf{D01BAY}$)

### d01rbf

Scheduled for withdrawal at Mark 28.
There is no replacement for this routine.
Withdrawn as a separate diagnostic routine is not required. The details of the computation, as stored in the parameters icom and com, are specified in Section 10.1 in d01raf.
See Section 10 in d01raf for further details.

## D02 – Ordinary Differential Equations

### d02pcf

Withdrawn at Mark 26.
Replaced by d02pef and associated D02P routines.

#### Old Code

```CALL d02pvf(N,TSTART,YINIT,TEND,TOL,THRESH,METHOD,'U',ERRASS, &
HSTART,W,LW,IFAIL)
...
CALL d02pcf(F,TWANT,T,Y,YP,YMAX,W,IFAIL)```

#### New Code

```IF (.Not. ERRASS) METHOD = -METHOD
CALL d02pqf(N,TSTART,TEND,YINIT,TOL,THRESH,METHOD,HSTART,IWSAV, &
RWSAV,IFAIL)
...
CALL d02pef(F2,N,TWANT,T,Y,YP,YMAX,IUSER,RUSER,IWSAV,RWSAV,IFAIL)```
iwsav is an integer array of length $130$ and rwsav is a real array of length $350+32×{\mathbf{n}}$.
iuser and ruser are arrays available to allow you to pass information to the user defined routine F2 (see f in d02pef).
The definition of F2 (see f in d02pef) can use the original routine f as follows:
```SUBROUTINE F2(T,N,Y,YP,IUSER,RUSER)
!     .. Scalar Arguments ..
Real (Kind=wp), Intent (In)      :: t
Integer, Intent (In)             :: n
!     .. Array Arguments ..
Real (Kind=wp), Intent (Inout)   :: ruser(1)
Real (Kind=wp), Intent (In)      :: y(n)
Real (Kind=wp), Intent (Out)     :: yp(n)
Integer, Intent (Inout)          :: iuser(1)
!     .. Procedure Arguments ..
External                         :: f
!     .. Executable Statements ..
Continue

Call f(t,y,yp)

Return
End Subroutine F2```

### d02pdf

Withdrawn at Mark 26.
Replaced by d02pff or d02pgf and associated D02P routines.
These replacements were made primarily for reasons of threadsafety. d02pgf also offers a reverse communication approach.

#### Old Code

```CALL d02pvf(N,TSTART,YINIT,TEND,TOL,THRESH,METHOD,'U',ERRASS, &
HSTART,W,LW,IFAIL)
...
CALL d02pdf(F,T,Y,YP,WORK,IFAIL)```

#### New Code

```IF (.Not. ERRASS) METHOD = -METHOD
CALL d02pqf(N,TSTART,TEND,YINIT,TOL,THRESH,METHOD,HSTART,IWSAV, &
RWSAV,IFAIL)
...
CALL d02pff(F2,N,T,Y,YP,IUSER,RUSER,IWSAV,RWSAV,IFAIL)```
iwsav is an integer array of length $130$ and rwsav is a real array of length $350+32×{\mathbf{n}}$.
iuser and ruser are arrays available to allow you to pass information to the user defined routine F2 (see f in d02pef).
The definition of F2 (see f in d02pef) can use the original routine f as follows:
```SUBROUTINE F2(T,N,Y,YP,IUSER,RUSER)
!     .. Scalar Arguments ..
Real (Kind=wp), Intent (In)      :: t
Integer, Intent (In)             :: n
!     .. Array Arguments ..
Real (Kind=wp), Intent (Inout)   :: ruser(1)
Real (Kind=wp), Intent (In)      :: y(n)
Real (Kind=wp), Intent (Out)     :: yp(n)
Integer, Intent (Inout)          :: iuser(1)
!     .. Procedure Arguments ..
External                         :: f
!     .. Executable Statements ..
Continue

Call f(t,y,yp)

Return
End Subroutine F2```

### d02pvf

Withdrawn at Mark 26.
Replaced by d02pqf.
See d02pcf and d02pdf for further information.

### d02pwf

Withdrawn at Mark 26.
Replaced by d02prf.

#### Old Code

`CALL d02pwf(TENDNU,IFAIL)`

#### New Code

`CALL d02prf(TENDNU,IWSAV,RWSAV,IFAIL)`
iwsav is an integer array of length $130$ and rwsav is a real array of length $350$.

### d02pxf

Withdrawn at Mark 26.
Replaced by d02psf.

#### Old Code

```CALL d02pxf(TWANT,REQEST,NWANT,YWANT,YPWANT,F,WORK,WRKINT, &
LENINT,IFAIL)```

#### New Code

```If (REQEST=='S' .or. REQEST=='s') Then
IDERIV = 0
Else if (REQEST=='D' .or. REQEST=='d') Then
IDERIV = 1
Else
IDERIV = 2
End If
CALL d02psf(TWANT,IDERIV,NWANT,YWANT,YPWANT,F2,WORKINT, &
LENINT,IUSER,RUSER,IWSAV,RWSAV,IFAIL)```
iwsav is an integer array of length $130$ and rwsav is a real array of length $350+32×{\mathbf{n}}$.
iuser and ruser are arrays available to allow you to pass information to the user defined routine F2 (see f in d02psf).
wcomm is a real array of length lwcomm. See the routine document for d02psf for further information.
The definition of F2 (see f in d02psf) can use the original routine f as follows:
```SUBROUTINE F2(T,N,Y,YP,IUSER,RUSER)
!     .. Scalar Arguments ..
Real (Kind=wp), Intent (In)      :: t
Integer, Intent (In)             :: n
!     .. Array Arguments ..
Real (Kind=wp), Intent (Inout)   :: ruser(1)
Real (Kind=wp), Intent (In)      :: y(n)
Real (Kind=wp), Intent (Out)     :: yp(n)
Integer, Intent (Inout)          :: iuser(1)
!     .. Procedure Arguments ..
External                         :: f
!     .. Executable Statements ..
Continue

Call f(t,y,yp)

Return
End Subroutine F2```

### d02pyf

Withdrawn at Mark 26.
Replaced by d02ptf.

#### Old Code

`Call d02pyf(TOTFCN,STPCST,WASTE,STPSOK,HNEXT,IFAIL)`

#### New Code

```Call d02ptf(TOTFCN,STPCST,WASTE,STPSOK,HNEXT,IWSAV, &
RWSAV,IFAIL)```

### d02pzf

Withdrawn at Mark 26.
Replaced by d02puf.

#### Old Code

`Call d02pzf(RMSERR,ERRMAX,TERRMX,WORK,IFAIL)`

#### New Code

`Call d02puf(N,RMSERR,ERRMAX,TERRMX,IWSAV,RWSAV,IFAIL)`
n must be unchanged from that passed to d02pqf.
iwsav is an integer array of length $130$ and rwsav is a real array of length $350+32×{\mathbf{n}}$.

### d02tkf

Withdrawn at Mark 27.
Replaced by d02tlf.

#### Old Code

`Call d02tkf(FFUN,FJAC,GAFUN,GBFUN,GAJAC,GBJAC,GUESS,RCOMM,ICOMM,IFAIL)`

#### New Code

```Call d02tlf(FFUN,FJAC,GAFUN,GBFUN,GAJAC,GBJAC,GUESS,RCOMM,ICOMM,IUSER, &
RUSER,IFAIL)```
The arrays iuser and ruser are also supplied as an additional two arguments to the seven user-supplied routines. These arrays are free to use to supply information to the seven routine arguments.

## D03 – Partial Differential Equations

### d03ryf

Withdrawn at Mark 27.
There is no replacement for this routine.

## E01 – Interpolation

### e01sef

Withdrawn at Mark 20.
Replaced by e01sgf.

#### Old Code

`CALL e01sef(M,X,Y,F,RNW,RNQ,NW,NQ,FNODES,MINNQ,WRK,IFAIL)`

#### New Code

`CALL e01sgf(M,X,Y,F,NW,NQ,IQ,LIQ,RQ,LRQ,IFAIL)`
e01sef has been superseded by e01sgf which gives improved accuracy, facilities for obtaining gradient values and a consistent interface with e01tgf for interpolation of scattered data in three dimensions.
The interpolant generated by the two routines will not be identical, but similar results may be obtained by using the same values of nw and nq. Details of the interpolant are passed to the evaluator through the arrays iq and rq rather than fnodes and rnw.

### e01sff

Withdrawn at Mark 20.
Replaced by e01shf.

#### Old Code

`CALL e01sff(M,X,Y,F,RNW,FNODES,PX,PY,PF,IFAIL)`

#### New Code

`CALL e01shf(M,X,Y,F,IQ,LIQ,RQ,LRQ,1,PX,PY,PF,QX,QY,IFAIL)`
The two calls will not produce identical results due to differences in the generation routines e01sef and e01sgf. Details of the interpolant are passed from e01sgf through the arrays iq and rq rather than fnodes and rnw.
e01shf also returns gradient values in qx and qy and allows evaluation at arrays of points rather than just single points.

## E02 – Curve and Surface Fitting

### e02acf

Withdrawn at Mark 27.
Replaced by e02alf.

#### Old Code

`CALL e02acf(X, Y, N, A, M1, REF)`

#### New Code

`CALL e02alf(N, X, Y, M1, A, REF, IFAIL)`

## E04 – Minimizing or Maximizing a Function

### e04ccf/​e04cca

Withdrawn at Mark 24.
Replaced by e04cbf.

#### Old Code

```CALL e04ccf(N,X,F,TOL,IW,W1,W2,W3,W4,W5,W6,FUNCT,MONIT,MAXCAL,   &
IFAIL)
or
CALL E04CCA(N,X,F,TOL,IW,W1,W2,W3,W4,W5,W6,FUNCT2,MONIT2,MAXCAL, &
IUSER,RUSER,IFAIL)```

#### New Code

```CALL e04cbf(N,X,F,TOLF,TOLX,FUNCT2,MONIT3,MAXCAL,IUSER,RUSER,    &
IFAIL)```
``` SUBROUTINE         MONIT3(FMIN,FMAX,SIM,N,NCALL,SERROR,VRATIO,IUSER,
RUSER)
INTEGER            N, NCALL, IUSER(*)
REAL (KIND=nag_wp) FMIN, FMAX, SIM(N+1,N), SERROR, VRATIO, RUSER(*)

CALL MONIT2(FMIN,FMAX,SIM,N,N+1,NCALL,IUSER,RUSER)
!      Add code here to monitor the values of SERROR and VRATIO, if necessary
RETURN
END```

### e04dgf/​e04dga

Deprecated.
Replaced by e04kff.
The new solver e04kff is part of the NAG Optimization suite (see Section 3.1 in the E04 Chapter Introduction), therefore the definition of the objective function values and gradients need to be split into two separate subroutines. e04kff offers a significant improvement in performance over e04dgf/​e04dga as well as additional functionality, such as the addition of variable bounds and monitoring.
Callbacks

#### Old Code

```      Subroutine objfun(mode,n,x,objf,objgrd,nstate,iuser,ruser)

!       .. Implicit None Statement ..
Implicit None
!       .. Scalar Arguments ..
Real (Kind=nag_wp), Intent (Out) :: objf
Integer, Intent (Inout)        :: mode
Integer, Intent (In)           :: n, nstate
!       .. Array Arguments ..
Real (Kind=nag_wp), Intent (Out) :: objgrd(n)
Real (Kind=nag_wp), Intent (Inout) :: ruser(*)
Real (Kind=nag_wp), Intent (In) :: x(n)
Integer, Intent (Inout)        :: iuser(*)
!       .. Executable Statements ..

!       Compute objective at point x
objf = ...

If (mode==2) Then
!         Compute objective gradient at point x
objgrd(1:n) = (/.../)
End If

Return
End Subroutine objfun```

#### New Code

```      Subroutine objfun(nvar,x,fx,inform,iuser,ruser,cpuser)

!       .. Use Statements ..
Use, Intrinsic                 :: iso_c_binding, Only: c_ptr
!       .. Implicit None Statement ..
Implicit None
!       .. Scalar Arguments ..
Type (c_ptr), Intent (In)      :: cpuser
Real (Kind=nag_wp), Intent (Out) :: fx
Integer, Intent (Inout)        :: inform
Integer, Intent (In)           :: nvar
!       .. Array Arguments ..
Real (Kind=nag_wp), Intent (Inout) :: ruser(*)
Real (Kind=nag_wp), Intent (In) :: x(nvar)
Integer, Intent (Inout)        :: iuser(*)
!       .. Executable Statements ..

!       Compute objective at point x
fx = ...

Return
End Subroutine objfun

Subroutine objgrd(nvar,x,nnzfd,fdx,inform,iuser,ruser,cpuser)

!       .. Use Statements ..
Use, Intrinsic                 :: iso_c_binding, Only: c_ptr
!       .. Implicit None Statement ..
Implicit None
!       .. Scalar Arguments ..
Type (c_ptr), Intent (In)      :: cpuser
Integer, Intent (Inout)        :: inform
Integer, Intent (In)           :: nnzfd, nvar
!       .. Array Arguments ..
Real (Kind=nag_wp), Intent (Inout) :: fdx(nvar), ruser(*)
Real (Kind=nag_wp), Intent (In) :: x(nvar)
Integer, Intent (Inout)        :: iuser(*)
!       .. Executable Statements ..

!       Compute objective gradient at point x
fdx(1:nnzfd) = (/.../)

Return
End Subroutine objgrd```
Main Call

#### Old Code

```      ifail = -1
Call e04dgf(n,objfun,iter,objf,objgrd,x,iwork,work,iuser,ruser,ifail)```

#### New Code

```!     .. Use Statements ..
Use, Intrinsic                   :: iso_c_binding, Only: c_null_ptr, c_ptr

!     ...

!     Initialize problem with n variables
ifail = 0
Call e04raf(handle,n,ifail)

!     (dependent on all variables)
ifail = 0
Call e04rgf(handle,n,(/(j,j=1,n)/),ifail)

!     Solve the problem
cpuser = c_null_ptr
ifail = -1
Call e04kff(handle,objfun,objgrd,e04kfu,n,x,rinfo,stats,iuser,ruser,     &
cpuser,ifail)

iter = stats(8)
objf = rinfo(1)

!     Free the handle memory
ifail = 0
Call e04rzf(handle,ifail)```

### e04jcf

Deprecated.
Replaced by e04jdf and e04jef.
e04jdf and e04jef are a part of the new NAG Optimization Suite (see Section 3.1 in the E04 Chapter Introduction) which allows you to define and solve various problems in a uniform manner. They also offer various algorithmic additions, such as a performance improvement on noisy problems, a possibility to progress towards the solution earlier than after $n$ initial function evaluations, heuristic stopping criteria, recovery from unavailable function evaluations, and various other algorithmic updates and tuning. In addition, e04jef offers a reverse communication interface (see Section 3.2 in the E04 Chapter Introduction) which might be useful in some environments or if you can parallelize your function evaluations.
Callbacks
The callback functions are almost identical but the new one has added cpuser argument and the order of arguments is slightly different.

#### Old Code

`      Subroutine objfun(n,x,f,iuser,ruser,inform)`

#### New Code

`      Subroutine objfun(n,x,f,inform,iuser,ruser,cpuser)`
Main Call

#### Old Code

```      ifail = -1
Call e04jcf(objfun,n,npt,x,bl,bu,rhobeg,rhoend,e04jcp,maxcal,f,nf,iuser, &
ruser,ifail)```

#### New Code

```!     .. Use Statements ..
Use, Intrinsic                   :: iso_c_binding, Only: c_null_ptr, c_ptr

!     ...

!     Initialize problem with n variables
ifail = 0
Call e04raf(handle,n,ifail)

!     Add nonlinear objective function which depends on all variables
ifail = 0
Call e04rgf(handle,n,(/(j,j=1,n)/),ifail)

!     Add bounds for the variables
Call e04rhf(handle,n,bl,bu,ifail)

!     Pass npt, rhobeg, rhoend, maxcal as options if different from defaults
ifail = 0
Call e04zmf(handle,'DFO Number Interp Points = <your npt>',ifail)
Call e04zmf(handle,'DFO Starting Trust Region = <your rhobeg>',ifail)
Call e04zmf(handle,'DFO Trust Region Tolerance = <your rhoend>',ifail)
Call e04zmf(handle,'DFO Max Objective Calls = <your maxcal>',ifail)

!     Solve the problem
cpuser = c_null_ptr
ifail = -1
Call e04jdf(handle,objfun,e04jdu,n,x,rinfo,stats,iuser,ruser,cpuser,ifail)

!     retrieve the objective value and the total number of calls made to the
!     objective function
f = rinfo(1)
nf = stats(1)

!     Free the handle memory
ifail = 0
Call e04rzf(handle,ifail)```

### e04mzf

Deprecated.
Replaced by e04mxf.

#### Old Code

```mpslst = .false.    ! or = .true.
Call e04mzf(infile,maxn,maxm,maxnnz,xbldef,xbudef,mpslst,n,m,nnz,iobj,
ncolh,a,irowa,iccola,bl,bu,start,pnames,nname,crname,xs,istate,ifail)
if (ifail==1 .or. ifail==2 .or. ifail==3) then
! not enough memory, allocate bigger arrays as given in M, N, NNZ
! and call e04mzf again
else if (ifail>=4 .and. ifail<=16) then
! MPS input file formating error
! stop
else if (ifail==17) then
! wrong arguments to e04mzf
! stop
else if (ifail==0) then
! data successfully read, call solver
end if```

#### New Code

```mpslst = 0      ! or = 1
Call e04mxf(infile,maxn,maxm,maxnnz,maxncolh,maxnnzh,maxlintvar,
mpslst,n,m,nnz,ncolh,nnzh,lintvar,iobj,a,irowa,iccola,bl,bu,pnames,
nname,crname,h,irowh,iccolh,minmax,intvar,ifail)

if (ifail==2) then
! not enough memory, allocate bigger arrays as given in M, N, NNZ, NCOLH,
! NNZH, LINTVAR and call e04mxf again
else if (ifail>=3 .and. ifail<=35) then
! MPS input file formatting error
! stop
else if (ifail==36) then
! wrong input argument
! stop
else if (ifail=-999) then
! internal memory allocation error
! stop
else if (ifail==0 .or. ifail==1) then
! data successfully read (with possible warning)
start = 'C'
do j=1,n
xs(j) = min(max(0.0_nag_wp,bl(j)),bu(j))
istate(j) = 0.0_nag_wp
end do
! call solver
end if```
e04mxf has extended the functionality of e04mzf and the interface has changed substantially. If there are integer variables, a quadratic part of the objective function or OBJSENSE section, e04mxf will read them and return them in the new arguments (lintvar, intvar, ncolh, nnzh, h, irowh, iccolh and minmax), with e04mzf these caused a file formatting error. The new routine e04mxf might also accept a slightly misformatted input file and return a warning ${\mathbf{ifail}}={\mathbf{1}}$.
The type of the argument mpslst has changed from logical to integer.
The parameters xbldef and xbudef of e04mzf were removed and fixed in e04mxf to their default values $0$ and ${10}^{20}$, respectively. Note that value ${10}^{20}$ within bounds is interpreted in our solvers as $+\infty$ (unconstrained).
Parameters start, xs and istate were present in e04mzf only for the convenience of calling the solver routine e04nkf/​e04nka and have been removed from e04mxf.
Should you need any assistance, please do not hesitate to contact NAG.

### e04unf

Withdrawn at Mark 22.
Replaced by e04usf/​e04usa.

#### Old Code

```CALL e04unf(M,N,NCLIN,NCNLN,LDA,LDCJ,LDFJ,     &
LDR,A,BL,BU,Y,CONFUN,OBJFUN,ITER,  &
ISTATE,C,CJAC,F,FJAC,CLAMDA,OBJF,  &
R,X,IWORK,LIWORK,WORK,LWORK,IUSER, &
RUSER,IFAIL)```

#### New Code

```CALL e04usf(M,N,NCLIN,NCNLN,LDA,LDCJ,LDFJ,     &
LDR,A,BL,BU,Y,CONFUN,OBJFUN,ITER,  &
ISTATE,C,CJAC,F,FJAC,CLAMDA,OBJF,  &
R,X,IWORK,LIWORK,WORK,LWORK,IUSER, &
RUSER,IFAIL)```
The specification of the subroutine objfun must also be changed as follows:

#### Old Code

```SUBROUTINE         OBJFUN(MODE,M,N,LDFJ,X,F,FJAC,NSTATE,IUSER,RUSER)
INTEGER            MODE,M,N,LDFJ,NSTATE,IUSER(*)
REAL (KIND=nag_wp) X(N),F(*),FJAC(LDFJ,*),RUSER(*)```

#### New Code

```SUBROUTINE         OBJFUN(MODE,M,N,LDFJ,NEEDFI,X,F,FJAC,NSTATE,      &
IUSER,RUSER)
INTEGER            MODE,M,N,LDFJ,NEEDFI,NSTATE,IUSER(*)
REAL (KIND=nag_wp) X(N),F(*),FJAC(LDFJ,*),RUSER(*)```
See the routine documents for further information.

### e04zcf/​e04zca

Withdrawn at Mark 24.
There is no replacement for this routine.

## F02 – Eigenvalues and Eigenvectors

### f02bjf

Withdrawn at Mark 23.
Replaced by f08waf.

#### Old Code

`CALL f02bjf(N,A,LDA,B,LDB,EPS1,ALFR,ALFI,BETA,MATV,V,LDV,ITER,IFAIL)`

#### New Code

```IF (MATV) THEN
JOBVR = 'V'
ELSE
JOBVR = 'N'
ENDIF
CALL dggev('N',JOBVR,N,A,LDA,B,LDB,ALFR,ALFI,BETA,VL,LDVL,           &
VR,LDVL,WORK,LWORK,INFO)
IF (INFO.EQ.0) THEN
...```

### f02eaf

Withdrawn at Mark 23.
Replaced by f08paf.

#### Old Code

`CALL f02eaf(JOB,N,A,LDA,WR,WI,Z,LDZ,WORK,LWORK,IFAIL)`

#### New Code

```LOGICAL SELECT
EXTERNAL SELECT
...
IF (JOB.EQ.'N') THEN
JOBVS = 'N'
ELSE
JOBVS = 'V'
END IF
CALL dgees(JOBVS,'N',SELECT,N,A,LDA,0,WR,WI,Z,LDZ,WORK, &
LWORK,BWORK,INFO)
IF (INFO.EQ.0) THEN
....
LOGICAL FUNCTION SELECT(AR,AI)
REAL (KIND=nag_wp) :: AR, AI
SELECT = .TRUE.
RETURN
END```

### f02ebf

Withdrawn at Mark 23.
Replaced by f08naf.

#### Old Code

```CALL f02ebf(JOB,N,A,LDA,WR,WI,VR,LDVR,VI,LDVI,WORK,LWORK, &
IFAIL)```

#### New Code

```IF (JOB.EQ.'N') THEN
JOBVR = 'N'
ELSE
JOBVR = 'V'
END IF
CALL dgeev('N',JOBVR,N,A,LDA,WR,WI,VL,LDVL,VR1,LDVR1,     &
WORK,LWORK,INFO)
IF (INFO.EQ.0) THEN
!       Eigenvector information is stored differently.
!       For complex conjugate pairs (that is, corresponding
!       to the j-th eigenvector such that WI(j) is nonzero,
!       and WI(j) = -WI(j+1)), the real and imaginary parts
!       of the first of the pair of eigenvectors are stored
!       as consecutive columns of VR1: VR1(:,j), VR1(:,j+1).
!       The second in the pair is just the conjugate of the
!       first, so can be constructed by negating the
!       elements in VR1(:,j+1).
!       If the j-th eigenvector is real (WI(j)=0), the
!       corresponding real eigenvector is stored in the
!       j-th column of VR1, VR1(1:N,j).```

### f02faf

Withdrawn at Mark 23.
Replaced by f08faf.

#### Old Code

`CALL f02faf(JOB,UPLO,N,A,LDA,W,WORK,LWORK,IFAIL)`

#### New Code

```CALL dsyev(JOB,UPLO,N,A,LDA,W,WORK,LWORK,INFO)
IF (INFO.EQ.0) THEN
...```
The minimum workspace requirement has not increased but the requirement for optimal performance might be different. The workspace query mechanism ($\mathbf{lwork}=-1$) should be used to determine the requirement for optimal performance.

### f02fcf

Withdrawn at Mark 23.
Replaced by f08fbf.

#### Old Code

```CALL f02fcf(JOB,RANGE,UPLO,N,A,LDA,WL,WU,IL,IU,MEST,M,   &
W,Z,LDZ,WORK,LWORK,IWORK,IFAIL)```

#### New Code

```CALL dsyevx(JOB,RANGE,UPLO,N,A,LDA,WL,WU,IL,IU,ABSTOL,M, &
W,Z,LDZ,WORK,LWORK,IWORK,JFAIL,INFO)
IF (INFO.EQ.0) THEN
...```
The minimum workspace requirement has not increased but the requirement for optimal performance might be different. The workspace query mechanism ($\mathbf{lwork}=-1$) should be used to determine the requirement for optimal performance.

### f02fdf

Withdrawn at Mark 23.
Replaced by f08saf.

#### Old Code

`CALL f02fdf(ITYPE,JOB,UPLO,N,A,LDA,B,LDB,W,WORK,LWORK,IFAIL)`

#### New Code

```CALL dsygv(ITYPE,JOB,UPLO,N,A,LDA,B,LDB,W,WORK,LWORK,INFO)
IF (INFO.EQ.0) THEN
...```
The minimum workspace requirement has not increased but the requirement for optimal performance might be different. The workspace query mechanism ($\mathbf{lwork}=-1$) should be used to determine the requirement for optimal performance.

### f02fhf

Withdrawn at Mark 23.
Replaced by f08uaf.

#### Old Code

` CALL f02fhf(N,MA,A,LDA,MB,B,LDB,D,WORK,LWORK,IFAIL)`

#### New Code

``` CALL dsbgv('N','U',N,MA,MB,A,LDA,B,LDB,D,Z,LDZ,WORK,INFO)
IF (INFO.EQ.0) THEN
...```
The order of eigenvalues in D changes from descending to ascending.
The minimum workspace requirement has changed to become $\mathtt{LWORK}=3×{\mathbf{n}}$

### f02gaf

Withdrawn at Mark 23.
Replaced by f08pnf.

#### Old Code

`CALL f02gaf(JOB,N,A,LDA,W,Z,LDZ,RWORK,WORK,LWORK,IFAIL)`

#### New Code

```LOGICAL BWORK(1)
LOGICAL SELECT
EXTERNAL SELECT
...
IF (JOB.EQ.'N') THEN
JOBVS = 'N'
ELSE
JOBVS = 'V'
END IF
CALL zgees(JOBVS,'N',SELECT,N,A,LDA,0,W,Z,LDZ,          &
WORK,LWORK,RWORK,BWORK,INFO)
IF (INFO.NE.0) THEN
...
LOGICAL FUNCTION SELECT(C)
COMPLEX*16 C
SELECT = .TRUE.
RETURN
END```
The minimum workspace requirement has not increased but the requirement for optimal performance might be different. The workspace query mechanism ($\mathbf{lwork}=-1$) should be used to determine the requirement for optimal performance.

### f02gbf

Withdrawn at Mark 23.
Replaced by f08nnf.

#### Old Code

`CALL f02gbf(JOB,N,A,LDA,W,V,LDV,RWORK,WORK,LWORK,IFAIL)`

#### New Code

```CALL zgeev('N',JOB,N,A,LDA,W,VL,LDVL,V,LDV,             &
WORK,LWORK,RWORK,INFO)
IF (INFO.EQ.0) THEN
...```

### f02gjf

Withdrawn at Mark 23.
Replaced by f08wnf.

#### Old Code

``` CALL f02gjf(N,AR,LDAR,AI,LDAI,BR,LDBR,BI,LDBI,EPS1,ALFR,  &
ALFI,BETA,MATV,VR,LDVR,VI,LDVI,ITER,IFAIL)```

#### New Code

``` IF (MATV) THEN
JOBVR = 'V'
ELSE
JOBVR = 'N'
END IF

!     Set A=AR + iAI and B = BR+iBI

CALL zggev('N',JOBVR,N,A,LDA,B,LDB,ALPHA,BETA1,VL,LDVL,   &
V,LDV,WORK,LWORK,RWORK,INFO)
IF (INFO.EQ.0) THEN
...```
Note that the separated real and imaginary parts of input and output data in f02gjf has been replaced by combined complex types in f08wnf.

### f02haf

Withdrawn at Mark 23.
Replaced by f08fnf.

#### Old Code

` CALL f02haf(JOB,UPLO,N,A,LDA,W,RWORK,WORK,LWORK,IFAIL)`

#### New Code

``` CALL zheev(JOB,UPLO,N,A,LDA,W,WORK,LWORK,RWORK,INFO)
IF (INFO.EQ.0) THEN
...```
The minimum workspace requirement has not increased but the requirement for optimal performance might be different. The workspace query mechanism ($\mathbf{lwork}=-1$) should be used to determine the requirement for optimal performance.

### f02hcf

Withdrawn at Mark 23.
Replaced by f08fpf.

#### Old Code

``` CALL f02hcf(JOB,RANGE,UPLO,N,A,LDA,WL,WU,IL,IU,MEST,M,   &
W,Z,LDZ,WORK,LWORK,RWORK,IWORK,IFAIL)```

#### New Code

``` CALL zheevx(JOB,RANGE,UPLO,N,A,LDA,WL,WU,IL,IU,ABSTOL,M, &
W,Z,LDZ,WORK,LWORK,RWORK,IWORK,JFAIL,INFO)
IF (INFO.EQ.0) THEN
...```
The minimum workspace requirement has not increased but the requirement for optimal performance might be different. The workspace query mechanism ($\mathbf{lwork}=-1$) should be used to determine the requirement for optimal performance.

### f02hdf

Withdrawn at Mark 23.
Replaced by f08snf.

#### Old Code

``` CALL f02hdf(ITYPE,JOB,UPLO,N,A,LDA,B,LDB,W,RWORK,WORK, &
LWORK,IFAIL)```

#### New Code

``` CALL zhegv(ITYPE,JOB,UPLO,N,A,LDA,B,LDB,W,WORK,LWORK,  &
RWORK,INFO)
IF (INFO.EQ.0) THEN
...```
The minimum workspace requirement has not increased but the requirement for optimal performance might be different. The workspace query mechanism ($\mathbf{lwork}=-1$) should be used to determine the requirement for optimal performance.

### f02sdf

Withdrawn at Mark 27.
Replaced by f12agf and f12fgf.
The replacement routines f12fgf (symmetric case) and f12agf (nonsymmetric case) are threaded for parallel execution in multithreaded implementations. These routines are based on the ARPACK package and make calls to BLAS/LAPACK routines. These may be threaded within the vendor library used by the implementation, which provides an additional opportunity for multithreaded performance.

#### Old Code

```CALL f02sdf(N,MA+1,MB+1,A,LDA,B,LDB,SYM,RELEP,RMU,VEC,D,IWORK,WORK, &
LWORK,IFAIL)```

#### New Code

```LICOMM = 140
LCOMM = 3*N + 3*NCV*NCV + 6*NCV + 60
ALLOCATE (COMM(LCOMM),DR(NCV),DI(NCV),RESID(N),V(N,NCV), &
ICOMM(LICOMM))
!   B is symmetric definite:
IF (B_symm_def) THEN
CALL f12aff(N,1,NCV,ICOMM,LICOMM,COMM,LCOMM,IFAIL)
CALL f12agf(KL,KU,A,LDA,B,LDB,RMU,0.0,NCONV,DR,DI,V,N,RESID, &
V,LDV,COMM,ICOMM,IFAIL)
VEC(1:N) = V(1:N,1)
ELSE
CALL f12aaf(N,NEV,NCV,ICOMM,LICOMM,COMM,LCOMM,IFAIL)
ALLOCATE(C(LDA,N),IPIV(N),X(N),MX(N))
C = A - RMU*B
CALL dgbtrf(N,N,KL,KU,C,LDA,IPIV,INFO)
IREVCM = 0
DO WHILE (IREVCM/=5)
CALL f12abf(IREVCM,RESID,V,LDV,X,MX,NSHIFT,COMM,ICOMM,IFAIL)
IF (IREVCM==-1 .OR. IREVCM==1) THEN
!     Perform  x <--- OP*x = inv[A-SIGMA*B]*Bx.
CALL dgbmv('N',N,N,KL,KU,ONE,B,LDB,X,1,ZERO,MX,1)
X(1:N) = MX(1:N)
CALL dgbtrs('N',N,KL,KU,1,C,LDA,IPIV,X,N,INFO)
END IF
END DO
! Post-process using f12acf to compute eigenvalue.
CALL F12ACF(NCONV,DR,DI,V,LDV,RMU,0.0,RESID,V,N,COMM,ICOMM,IFAIL)
LR = DR(1)/(DR(1)**2+DI(1)**2) + RMU
END IF```

### f02wdf

Withdrawn at Mark 27.
Replaced by f02wuf and f08aef.
This routine is replaced for multithreaded performance and the ability to benefit from vendor library performance (BLAS/LAPACK).
Note: only the multithreaded implementations of f02wdf were able to benefit from parallelism or vendor BLAS/LAPACK performance.
The Householder $QU$ factorization part of the functionality can be achieved with f08aef. The action ${Q}^{\mathrm{T}}b$ can be computed by a call to f08agf. The orthogonal matrix $Q$ can be explicitly constructed, in-place, by a subsequent call to f08aff.
If the singular value decomposition (SVD) of $U$ is required, the result of f08aef must be fed to f02wuf, remembering that the first orthogonal matrix of the SVD is called $Q$ in f02wuf and $R$ in f02wdf

#### Old Code

```IFAIL = 0
CALL f02wdf(M,N,A,LDA,WANTB,B,TOL,SVD,IRANK,Z,SV,WANTR,R, &
LDR,WANTPT,PT,LDPT,WORK,LWORK,IFAIL)```

#### New Code

```LWORK = -1
CALL dgeqrf(M,N,A,LDA,Z,WORK,LWORK,INFO)
LWORK = ANINT(WORK(1))
DEALLOCATE (WORK)
ALLOCATE (WORK(LWORK))
CALL dgeqrf(M,N,A,LDA,Z,WORK,LWORK,INFO)
NCOLB = 1
IF (WANTB) THEN
CALL dormqr('L','T',M,NCOLB,N,A,LDA,Z,B,M,WORK,LWORK,INFO)
END IF
IF (.NOT. SVD) THEN
! construct Q explicitly, overwrites A
CALL dorgqr(M,M,A,LDA,Z,WORK,LWORK,INFO)
ELSE
! SVD factorization, PT overwrites A
DEALLOCATE (WORK)
ALLOCATE (WORK(5*N))
CALL f02wuf(N,A,LDA,NCOLB,B,M,WANTR,R,LDR,SV,WANTPT,WORK,IFAIL)
! compute rank
IRANK = F06KLF(N,SV,1,TOL)
END IF```

### f02wef

Withdrawn at Mark 23.
Replaced by f08kbf.

#### Old Code

``` CALL f02wef(M,N,A,LDA,NCOLB,B,LDB,WANTQ,Q,LDQ,SV,WANTP,       &
PT,LDPT,WORK,IFAIL)```

#### New Code

``` IF (WANTQ) THEN
JOBU = 'A'
ELSE
JOBU = 'N'
END IF
IF (WANTP) THEN
JOBVT = 'A'
ELSE
JOBVT = 'N'
END IF
LWORK = -1
CALL dgesvd(JOBU,JOBVT,M,N,A,LDA,SV,Q,LDQ,PT,LDPT,WORK,LWORK,INFO)
LWORK = ANINT(WORK(1))
ALLOCATE (W(LWORK))

CALL dgesvd(JOBU,JOBVT,M,N,A,LDA,SV,Q,LDQ,PT,LDPT,W,LWORK,INFO)

DEALLOCATE (W)```
work must be a one-dimensional real array of length at least lwork given by: $\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left(1,3×\mathrm{min}\phantom{\rule{0.125em}{0ex}}\left({\mathbf{m}},{\mathbf{n}}\right)+\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left({\mathbf{m}},{\mathbf{n}}\right),5×\mathrm{min}\phantom{\rule{0.125em}{0ex}}\left({\mathbf{m}},{\mathbf{n}}\right)\right)$
Larger values of lwork, up to some optimal value, may improve performance.
Please note that the facility to return ${Q}^{\mathrm{T}}B$ is not provided so arguments $\mathbf{wantb}$ and $\mathbf{b}$ are not required. Instead, f08kbf has an option to return the entire ${\mathbf{m}}*{\mathbf{m}}$ orthogonal matrix $Q$, referred to as ${\mathbf{u}}$ in its documentation, through its 8th argument.

### f02xef

Withdrawn at Mark 23.
Replaced by f08kpf.

#### Old Code

``` CALL f02xef(M,N,A,LDA,NCOLB,B,LDB,WANTQ,Q,LDQ,SV,WANTP,        &
PH,LDPH,RWORK,CWORK,IFAIL)```

#### New Code

``` IF (WANTQ) THEN
JOBU = 'A'
ELSE
JOBU = 'N'
END IF
IF (WANTP) THEN
JOBVT = 'A'
ELSE
JOBVT = 'N'
END IF
LWORK = -1
CALL zgesvd(JOBU,JOBVT,M,N,A,LDA,SV,Q,LDQ,PT,LDPT,WORK,        &
LWORK,RWORK,INFO)
LWORK = ANINT(WORK(1))
ALLOCATE (W(LWORK))

CALL zgesvd(JOBU,JOBVT,M,N,A,LDA,SV,Q,LDQ,PT,LDPT,W,           &
LWORK,RWORK,INFO)

DEALLOCATE (W)```
work must be a one-dimensional complex array of length at least lwork given by $\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left(1,2×\mathrm{min}\phantom{\rule{0.125em}{0ex}}\left({\mathbf{m}},{\mathbf{n}}\right)+\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left({\mathbf{m}},{\mathbf{n}}\right)\right)$
work must be a one-dimensional real array of length $\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left(1,5×\mathrm{min}\phantom{\rule{0.125em}{0ex}}\left({\mathbf{m}},{\mathbf{n}}\right)\right)$.
Larger values of lwork, up to some optimal value, may improve performance.
Please note that the facility to return ${Q}^{\mathrm{H}}B$ is not provided so arguments wantb and $\mathbf{b}$ are not required. Instead, f08kpf has an option to return the entire ${\mathbf{m}}*{\mathbf{m}}$ unitary matrix $Q$, referred to as ${\mathbf{u}}$ in its documentation, through its 8th argument.

## F03 – Determinants

### f03aaf

Withdrawn at Mark 25.

#### Old Code

```IFAIL = 0
CALL f03aaf(A,LDA,N,DET,WKSPCE,IFAIL)```

#### New Code

```INTEGER IPIV(N)
...

CALL dgetrf(N,N,A,LDA,IPIV,INFO)
IFAIL = 0
CALL f03baf(N,A,LDA,IPIV,D,ID,IFAIL)
DET = D*2**ID```
Note: the real array wkspce has been replaced by the integer array ipiv for holding the pivots of the factorization.

### f03abf

Withdrawn at Mark 25.
Replaced by f07fdf and f03bff.

#### Old Code

```IFAIL = 0
CALL f03abf(A,LDA,N,DET,WKSPCE,IFAIL)```

#### New Code

```CALL dpotrf('U',N,A,LDA,INFO)
IFAIL = 0
CALL f03bff(N,A,LDA,D,ID,IFAIL)
DET = D*2**ID```
Note: the real array wkspce is no longer required. Also the upper triangular part of $A$, stored in a, has been replaced here by its Cholesky factorization; the lower triangular part of $A$ can be used and overwritten by replacing 'U' by 'L' in the call to dpotrf above.

### f03acf

Withdrawn at Mark 25.
Replaced by f07hdf and f03bhf.

#### Old Code

```IFAIL = 0
CALL f03acf(A,LDA,N,M,DET,RL,LDRL,M1,IFAIL)```

#### New Code

```CALL dpbtrf('L',N,M,AB,LDAB,INFO)
IFAIL = 0
CALL f03bhf('L',N,KD,AB,LDAB,D,ID,IFAIL)
DET = D*2**ID```
Note: the storage of $A$ in arrays a and ab is different. In fact ${\mathbf{ab}}\left(\mathit{i},\mathit{j}\right)=\mathbf{a}\left(\mathit{j},\mathit{i}\right)$, for $\mathit{i}=1,2,\dots ,m$ and $\mathit{j}=\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left(1,i-m\right),\dots ,i$ which conforms to the LAPACK banded storage scheme. The factorization is returned in ab rather than in a separate array (rl). The upper part of matrix $A$ can also be stored in ab on input to dpbtrf.

Withdrawn at Mark 25.
Replaced by f07arf and f03bnf.

#### Old Code

```IFAIL = 0

#### New Code

```INTEGER IPIV(N)
...
CALL zgetrf(N,N,A,LDA,IPIV,INFO)
IFAIL = 0
CALL f03bnf(N,A,LDA,IPIV,D,ID,IFAIL)
DETR = REAL(D)*2**ID(1)
DETI = AIMAG(D)*2**ID(2)```
Note: the real array wkspce has been replaced by the integer array ipiv for holding the pivots of the factorization. The real and imaginary parts of the determinant are independently scaled.

### f03aef

Withdrawn at Mark 25.
Replaced by f07fdf and f03bff.

#### Old Code

```IFAIL = 0
CALL f03aef(N,A,LDA,P,D1,ID,IFAIL)```

#### New Code

```CALL dpotrf('U',N,A,LDA,INFO)
IFAIL = 0
CALL f03bff(N,A,LDA,D1,ID,IFAIL)```
Note: the upper triangular part of $A$, stored in a, has been replaced here by its Cholesky factorization; the lower triangular part of $A$ can be used and overwritten by replacing ${\mathbf{uplo}}=\text{'U'}$ by ${\mathbf{uplo}}=\text{'L'}$ in the call to f07fdf above.

### f03aff

Withdrawn at Mark 25.

#### Old Code

```IFAIL = 0
CALL f03aff(N,EPS,A,LDA,D1,ID,P,IFAIL)```

#### New Code

```INTEGER IPIV(N)
...
CALL dgetrf(N,N,A,LDA,IPIV,INFO)
IFAIL = 0
CALL f03baf(N,A,LDA,IPIV,D1,ID,IFAIL)```
Note: real array p has been replaced by the integer array ipiv for holding the pivots of the factorization.

## F04 – Simultaneous Linear Equations

### f04aaf

Withdrawn at Mark 23.
Replaced by f07aaf.

#### Old Code

`CALL f04aaf(A,LDA,B,LDB,N,M,C,LDC,WKSPCE,IFAIL)`

#### New Code

```CALL dgesv(N,M,A,LDA,IPIV,B,LDB,INFO)
IF (INFO.EQ.0) THEN
...```

### f04abf

Scheduled for withdrawal at Mark 28.
Replaced by f07fbf.
f04abf and f04asf have been replaced by f07fbf for performance. The replacement routine is threaded by NAG and may also be threaded in the vendor library (BLAS/LAPACK).

#### Old Code

`CALL f04abf(A,LDA,B,LDB,N,M,C,LDC,WKSPCE,BB,LDBB,IFAIL)`

#### New Code

```CALL f04abf_wrap(a,lda,b,ldb,n,m,c,ldc,wkspce,bb,ldbb,ifail)

Subroutine f04abf_wrap(a,lda,b,ldb,n,m,c,ldc,wkspce,bb,ldbb,ifail)
!     .. Use Statements ..
Use nag_library, Only: dposvx, dsymm, nag_wp
!     .. Scalar Arguments ..
Integer, Intent (In)             :: lda, ldb, ldbb, ldc, m, n
Integer, Intent (Inout)          :: ifail
!     .. Array Arguments ..
Real (Kind=nag_wp), Intent (Inout) :: a(lda,*), b(ldb,*)
Real (Kind=nag_wp), Intent (Out)   :: bb(ldbb,m), c(ldc,m), wkspce(1)
!     .. Local Scalars ..
Real (Kind=nag_wp)               :: rcond, alpha, beta
Integer                          :: info, ldaf
Character (1)                    :: equed
!     .. Local Arrays ..
Real (Kind=nag_wp)               :: s(1)
Real (Kind=nag_wp), Allocatable  :: af(:,:), work(:), ferr(:), berr(:)
Integer, Allocatable             :: iwork(:)

ldaf = n
Allocate (af(ldaf,n),ferr(m),berr(m),work(3*n),iwork(n))
!     The NAG name equivalent of dposvx is f07fbf
Call dposvx('N','Upper',n,m,a,lda,af,ldaf,equed,s,b,ldb,  &
c,ldc,rcond,ferr,berr,work,iwork,info)

ifail = info
bb(1:n,1:m) = b(1:n,1:m)
alpha = -1.0_nag_wp
beta = 1.0_nag_wp
!     The NAG name equivalent of dgemm is f06yaf
Call dsymm('L','U',n,m,alpha,a,lda,c,ldc,beta,bb,ldbb)

End Subroutine f04abf_wrap```

### f04acf

Withdrawn at Mark 23.
Replaced by f07haf.

#### Old Code

` CALL f04acf(A,LDA,B,LDB,N,M,IR,C,LDC,RL,LDRL,M1,IFAIL)`

#### New Code

``` CALL dpbsv('U',N,M,IR,AB,LDAB,B,LDB,INFO)
IF (INFO.EQ.0) THEN
!        A and AB are stored differently.
!        AB may be regarded as the transpose of A, with the 'U' option.
!        Thus LDAB might be M+1
...```

Withdrawn at Mark 23.
Replaced by f07anf.

#### Old Code

` CALL f04adf(A,LDA,B,LDB,N,M,C,LDC,WKSPCE,IFAIL)`

#### New Code

``` CALL zgesv(N,M,A,LDA,IPIV,B,LDB,INFO)
IF (INFO.EQ.0) THEN
...```

### f04aef

Scheduled for withdrawal at Mark 28.
Replaced by f07abf.
f04aef and f04atf have been replaced by f07abf for performance. The replacement routine is threaded by NAG and may also be threaded in the vendor library (BLAS/LAPACK).

#### Old Code

`CALL f04aef(A,LDA,B,LDB,N,M,C,LDC,WKSPCE,AA,LDAA,BB,LDBB,IFAIL)`

#### New Code

```CALL f04aef_wrap(a,lda,b,ldb,n,m,c,ldc,wkspce,aa,ldaa,bb,ldbb,ifail)
Subroutine f04aef_wrap(a,lda,b,ldb,n,m,c,ldc,wkspce,aa,ldaa,bb,ldbb,ifail)
!     .. Use Statements ..
Use nag_library, Only: dgesvx, dgemm, nag_wp
!     .. Scalar Arguments ..
Integer, Intent (In)             :: lda, ldaa, ldb, ldbb, ldc, m, n
Integer, Intent (Inout)          :: ifail
!     .. Array Arguments ..
Real (Kind=nag_wp), Intent (Inout) :: a(lda,*), b(ldb,*)
Real (Kind=nag_wp), Intent (Out)   :: bb(ldbb,m), c(ldc,m), aa(ldaa,n), &
wkspce(1)
!     .. Local Scalars ..
Real (Kind=nag_wp)               :: rcond, alpha, beta
Integer                          :: info
Character (1)                    :: equed
!     .. Local Arrays ..
Real (Kind=nag_wp)               :: cscl(1), rscl(1)
Real (Kind=nag_wp), Allocatable  :: work(:), ferr(:), berr(:)
Integer, Allocatable             :: ipiv(:), iwork(:)

Allocate (berr(m),ferr(m),work(4*n),ipiv(n),iwork(n))
!     The NAG name equivalent of dgesvx is f07abf
Call dgesvx('N','N',n,m,a,lda,aa,ldaa,ipiv,equed,rscl,cscl,b,ldb, &
c,ldc,rcond,ferr,berr,work,iwork,info)

ifail = info
bb(1:n,1:m) = b(1:n,1:m)
alpha = -1.0_nag_wp
beta = 1.0_nag_wp
!     The NAG name equivalent of dgemm is f06yaf
Call dgemm('N','N',n,m,n,alpha,a,lda,c,ldc,beta,bb,ldbb)

End Subroutine f04aef_wrap```

### f04aff

Withdrawn at Mark 25.
There is no replacement for this routine.
The factorization and solution of a positive definite linear system can be handled by calls to routines from Chapter F07, e.g., f07fbf.
For example:

#### Old Code

```IFAIL = 0
CALL f03aef(N,A,LDA,P,D1,ID,IFAIL)
CALL f04aff(N,NRHS,A,LDA,P,B,LDB,EPS,X,LDX,BB,LDBB,K,IFAIL)```

#### New Code

```CALL dposvx('equil','upper',N,NRHS,A,LDA,AF,LDAF,'Yes',P,B, &
LDB,X,LDX,RCOND,FERR,BERR,WORK,IWORK,INFO)
IFAIL = 0
CALL f03bff(N,A,LDA,D1,ID,IFAIL)```

### f04agf

Withdrawn at Mark 25.
There is no replacement for this routine.
The factorization and solution of a positive definite linear system can be handled by calls to routines from Chapter F07, e.g., f07faf.
For example:

#### Old Code

```IFAIL = 0
CALL f03aef(N,A,LDA,P,D1,ID,IFAIL)
CALL f04agf(N,NRHS,A,LDA,P,B,LDB,X,LDX)```

#### New Code

```CALL dposv('upper',N,NRHS,A,LDA,B,LDB,INFO)
IFAIL = 0
CALL f03bff(N,A,LDA,D1,ID,IFAIL)```

### f04ahf

Withdrawn at Mark 25.
There is no replacement for this routine.
The factorization and solution of a real general linear system can be handled by calls to routines from the Chapter F07, e.g., f07abf.
For example:

#### Old Code

```IFAIL = 0
CALL f03aff(N,EPS,A,LDA,D1,ID,P,IFAIL)
CALL f04ahf(N,NRHS,A,LDA,AA,LDAA,P,B,LDB,EPS,X,LDX,BB,    &
LDBB,K,IFAIL)```

#### New Code

```CALL dgesvx('Equil','No trans',N,NRHS,A,LDA,AA,LDAA,IPIV, &
'Yes',R,C,B,LDB,X,LDX,RCOND,FERR,BERR,WORK,   &
IWORK,INFO)
IFAIL = 0
CALL f03baf(N,A,LDA,IPIV,D1,ID,IFAIL)```

### f04ajf

Withdrawn at Mark 25.
There is no replacement for this routine.
The factorization and solution of a real general linear system can be handled by calls to routines from Chapter F07, e.g., f07aaf.
For example:

#### Old Code

```IFAIL = 0
CALL f03aff(N,EPS,A,LDA,D1,ID,P,IFAIL)
CALL f04ajf(N,NRHS,A,LDA,P,B,LDB)```

#### New Code

```CALL dgesv(N,NRHS,A,LDA,IPIV,B,LDB,INFO)
IFAIL = 0
CALL f03baf(N,A,LDA,IPIV,D1,ID,IFAIL)```

### f04arf

Withdrawn at Mark 23.
Replaced by f07aaf.

#### Old Code

` CALL f04arf(A,LDA,B,N,C,WKSPCE,IFAIL)`

#### New Code

``` CALL dgesv(N,1,A,LDA,IPIV,B,N,INFO)
IF (INFO.EQ.0) THEN
...```

### f04asf

Scheduled for withdrawal at Mark 28.
Replaced by f07fbf.
f04abf and f04asf have been replaced by f07fbf for performance. The replacement routine is threaded by NAG and may also be threaded in the vendor library (BLAS/LAPACK).

#### Old Code

`CALL f04asf(A,LDA,B,N,C,WK1,WK2,IFAIL)`

#### New Code

```CALL f04asf_wrap(a,lda,b,n,c,wk1,wk2,ifail)
Subroutine f04asf_wrap(a,lda,b,n,c,wk1,wk2,ifail)
!     .. Use Statements ..
Use nag_library, Only: dposvx, nag_wp
!     .. Scalar Arguments ..
Integer, Intent (In)             :: lda, n
Integer, Intent (Inout)          :: ifail
!     .. Array Arguments ..
Real (Kind=nag_wp), Intent (Inout) :: a(lda,n), b(n)
Real (Kind=nag_wp), Intent (Out) :: c(n), wk1(1), wk2(1)
!     .. Local Scalars ..
Real (Kind=nag_wp)               :: rcond
Integer                          :: info, ldaf, m
Character (1)                    :: equed
!     .. Local Arrays ..
Real (Kind=nag_wp)               :: s(1), ferr(1), berr(1)
Real (Kind=nag_wp), Allocatable  :: af(:,:), work(:)
Integer, Allocatable             :: iwork(:)

ldaf = n
m = 1
Allocate (af(ldaf,n),work(3*n),iwork(n))
!     The NAG name equivalent of dposvx is f07fbf
Call dposvx('N','Upper',n,m,a,lda,af,ldaf,equed,s,b,n,  &
c,n,rcond,ferr,berr,work,iwork,info)
wk1(1) = rcond
wk2(1) = berr(1)

ifail = info
End Subroutine f04asf_wrap```

### f04atf

Scheduled for withdrawal at Mark 28.
Replaced by f07abf.
f04aef and f04atf have been replaced by f07abf for performance. The replacement routine is threaded by NAG and may also be threaded in the vendor library (BLAS/LAPACK).

#### Old Code

`Call f04atf(A,LDA,B,N,C,AA,LDAA,WKS1,WKS2,IFAIL)`

#### New Code

```CALL f04atf_wrap(a,lda,b,n,c,aa,ldaa,wks1,wks2,ifail)
Subroutine f04atf_wrap(a,lda,b,n,c,aa,ldaa,wks1,wks2,ifail)
!     .. Use Statements ..
Use nag_library, Only: dgesvx, nag_wp
!     .. Scalar Arguments ..
Integer, Intent (In)             :: lda, ldaa, n
Integer, Intent (Inout)          :: ifail
!     .. Array Arguments ..
Real (Kind=nag_wp), Intent (Inout) :: a(lda,*), b(n)
Real (Kind=nag_wp), Intent (Out)   :: c(n), aa(ldaa,n), wks1(1),
wks2(1)
!     .. Local Scalars ..
Real (Kind=nag_wp)               :: rcond
Integer                          :: info, m
Character (1)                    :: equed
!     .. Local Arrays ..
Real (Kind=nag_wp)               :: cscl(1), rscl(1), ferr(1),
berr(1)
Real (Kind=nag_wp), Allocatable  :: work(:)
Integer, Allocatable             :: ipiv(:), iwork(:)

m = 1
Allocate (work(4*n),ipiv(n),iwork(n))
!     The NAG name equivalent of dgesvx is f07abf
Call dgesvx('N','N',n,m,a,lda,aa,ldaa,ipiv,equed,rscl,cscl,b,n, &
c,n,rcond,ferr,berr,work,iwork,info)

ifail = info
wks1(1) = rcond
wks2(1) = berr(1)

End Subroutine f04atf_wrap```

### f04eaf

Withdrawn at Mark 23.
Replaced by f07caf.

#### Old Code

` CALL f04eaf(N,D,DU,DL,B,IFAIL)`

#### New Code

``` CALL dgtsv(N,1,DL(2),D,DU(2),B,N,INFO)
IF (INFO.EQ.0) THEN
...```

### f04faf

Withdrawn at Mark 23.
Replaced by f07jaf, or f07jdf and f07jef.

#### Old Code

` CALL f04faf(JOB,N,D,E,B,IFAIL)`

#### New Code

``` CALL dptsv(N,1,D,E(2),B,1,INFO)
...```

### f04jaf

Withdrawn at Mark 23.
Replaced by f08kaf.

#### Old Code

` CALL f04jaf(M,N,A,LDA,B,TOL,SIGMA,IRANK,WORK,LWORK,IFAIL)`

#### New Code

``` CALL dgelss(M,N,1,A,LDA,B,1,S,RCOND,IRANK,WORK,LWORK,INFO)
IF (INFO.EQ.0) THEN
!        Singular values now in S, not WORK.
!        The standard error is not computed
...```
The minimum workspace requirement has changed from $4×{\mathbf{n}}$ to $3×\mathrm{min}\phantom{\rule{0.125em}{0ex}}\left({\mathbf{n}},{\mathbf{m}}\right)+\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left(2×\mathrm{min}\phantom{\rule{0.125em}{0ex}}\left({\mathbf{n}},{\mathbf{m}}\right),\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left({\mathbf{m}},{\mathbf{n}}\right),1\right)$.

### f04jdf

Withdrawn at Mark 23.
Replaced by f08kaf.

#### Old Code

` CALL f04jdf(M,N,A,LDA,B,TOL,SIGMA,IRANK,WORK,LWORK,IFAIL)`

#### New Code

``` CALL dgelss(M,N,1,A,LDA,B,1,S,RCOND,IRANK,WORK,LWORK,INFO)
!     Note workspace requirements are different.
IF (INFO.EQ.0) THEN
!        Singular values now in S, not WORK.
!        The standard error is not computed
...```
The minimum workspace requirement has changed from $\mathbf{n}×\left(\mathbf{m}+4\right)$ to $3×\mathrm{min}\phantom{\rule{0.125em}{0ex}}\left({\mathbf{n}},{\mathbf{m}}\right)+\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left(2×\mathrm{min}\phantom{\rule{0.125em}{0ex}}\left({\mathbf{n}},{\mathbf{m}}\right),\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left({\mathbf{m}},{\mathbf{n}}\right),1\right)$.

### f04jlf

Withdrawn at Mark 23.
Replaced by f08zbf.

#### Old Code

` CALL f04jlf(M,N,P,A,LDA,B,LDB,D,X,Y,WORK,LWORK,IFAIL)`

#### New Code

``` CALL dggglm(M,N,P,A,LDA,B,LDB,D,X,Y,WORK,LWORK,INFO)
IF (INFO.EQ.0) THEN
...```
The minimum workspace requirement has not increased but the requirement for optimal performance might be different. The workspace query mechanism ($\mathbf{lwork}=-1$) should be used to determine the requirement for optimal performance.

### f04jmf

Withdrawn at Mark 23.
Replaced by f08zaf.

#### Old Code

` CALL f04jmf(M,N,P,A,LDA,B,LDB,C,D,X,WORK,LWORK,IFAIL)`

#### New Code

``` CALL dgglse(M,N,P,A,LDA,B,LDB,C,D,X,WORK,LWORK,INFO)
IF (INFO.EQ.0) THEN
...```
The minimum workspace requirement has not increased but the requirement for optimal performance might be different. The workspace query mechanism ($\mathbf{lwork}=-1$) should be used to determine the requirement for optimal performance.

### f04klf

Withdrawn at Mark 23.
Replaced by f08zpf.

#### Old Code

`CALL f04klf(M,N,P,A,LDA,B,LDB,D,X,Y,WORK,LWORK,IFAIL)`

#### New Code

```CALL zggglm(M,N,P,A,LDA,B,LDB,D,X,Y,WORK,LWORK,INFO)
IF (INFO.EQ.0) THEN
...```

### f04kmf

Withdrawn at Mark 23.
Replaced by f08znf.

#### Old Code

`CALL f04kmf(M,N,P,A,LDA,B,LDB,C,D,X,WORK,LWORK,IFAIL)`

#### New Code

```CALL zgglse(M,N,P,A,LDA,B,LDB,C,D,X,WORK,LWORK,INFO)
IF (INFO.EQ.0) THEN
...```

### f04ycf

Withdrawn at Mark 26.
Replaced by f04ydf.
f04ydf employs a better algorithm (see Higham and Tisseur (2000)).

#### Old Code

`CALL f04ycf(ICASE,N,X,ESTNRM,WORK,IWORK,IFAIL)`

#### New Code

`CALL f04ydf(IREVCM,M,N,X,LDX,Y,LDY,ESTNRM,T,SEED,WORK,IWORK,IFAIL)`
f04ydf returns an estimate of the $1$-norm of a rectangular $M×N$ matrix, whereas f04ycf only works with square matrices. The real array X, which was previously used to return matrix–vector products to f04ycf, has been replaced with two real arrays ${\mathbf{x}}\left({\mathbf{ldx}},*\right)$ and ${\mathbf{y}}\left({\mathbf{ldy}},*\right)$ which are used to return matrix-matrix products to f04ydf. Here, ${\mathbf{ldx}}\ge {\mathbf{n}}$, ${\mathbf{ldy}}\ge {\mathbf{m}}$ and the second dimensions of x and y are at least of size t, where you can choose argument t. The sizes of the workspace arrays work and iwork have been increased to ${\mathbf{m}}×{\mathbf{t}}$ and $2×{\mathbf{n}}+5×{\mathbf{t}}+20$ respectively. The integer seed provides a seed for the random number generator used by f04ydf. The integer icase has been replaced by irevcm, which can take the values $0$, $1$ or $2$. See the routine documentation for f04ydf further details about the reverse communication interface.

### f04zcf

Withdrawn at Mark 26.
Replaced by f04zdf.
f04zdf employs a better algorithm (see Higham and Tisseur (2000).

#### Old Code

`CALL f04zcf(ICASE,N,X,ESTNRM,WORK,IFAIL)`

#### New Code

`CALL f04zdf(IREVCM,M,N,X,LDX,Y,LDY,ESTNRM,T,SEED,WORK,RWORK,IWORK,IFAIL)`
f04zdf returns an estimate of the $1$-norm of a rectangular $M×N$ matrix, whereas f04zcf only works with square matrices. The complex array X, which was previously used to return matrix–vector products to f04zcf, has been replaced with two complex arrays ${\mathbf{x}}\left({\mathbf{ldx}},*\right)$ and ${\mathbf{y}}\left({\mathbf{ldy}},*\right)$ which are used to return matrix-matrix products to f04zdf. Here, ${\mathbf{ldx}}\ge {\mathbf{n}}$, ${\mathbf{ldy}}\ge {\mathbf{m}}$ and the second dimensions of x and y are at least of size t, where you can choose the argument t. The sizes of the workspace arrays work and iwork have been increased to ${\mathbf{m}}×{\mathbf{t}}$ and $2×{\mathbf{n}}+5×{\mathbf{t}}+20$ respectively and there is an additional real workspace array rwork of size $2×{\mathbf{n}}$. The integer seed provides a seed for the random number generator used by f04zdf. The integer icase has been replaced by irevcm, which can take the values $0$, $1$ or $2$. See the routine documentation for f04zdf for further details about the reverse communication interface.

## F08 – Least Squares and Eigenvalue Problems (LAPACK)

For each of the deprecated routines listed below, e.g., f08bef, the replacement routine has one additional argument lwork, used to supply the length of the array argument work or to perform a workspace query by setting ${\mathbf{lwork}}=-1$. It is safest to perform a workspace query first using a dummy workspace array of length $1$ and then to allocate the array of length equal to the optimal value returned in the dummy workspace, e.g.:
```LWORK = -1
CALL dgeqp3(...,DUMWRK,LWORK,...)
LWORK = NINT(DUMWRK(1))
ALLOCATE (WORK(LWORK))
CALL dgeqp3(...,WORK,LWORK,...)```

### f08bef

Deprecated.
Replaced by f08bff.

### f08bsf

Deprecated.
Replaced by f08btf.

### f08vaf

Deprecated.
Replaced by f08vcf.

### f08vef

Deprecated.
Replaced by f08vgf.

### f08vnf

Deprecated.
Replaced by f08vqf.

### f08vsf

Deprecated.
Replaced by f08vuf.

### f08waf

Deprecated.
Replaced by f08wcf.

### f08wef

Deprecated.
Replaced by f08wff.

### f08wnf

Deprecated.
Replaced by f08wqf.

### f08wsf

Deprecated.
Replaced by f08wtf.

### f08xaf

Deprecated.
Replaced by f08xcf.

### f08xnf

Deprecated.
Replaced by f08xqf.

## F11 – Large Scale Linear Systems

### f11baf

Withdrawn at Mark 21.
Replaced by f11bdf.

#### Old Code

```CALL f11baf(METHOD,PRECON,NORM,WEIGHT,ITERM,N,M,TOL,MAXITN, &
ANORM,SIGMAX,MONIT,LWREQ,IFAIL)```

#### New Code

```CALL f11bdf(METHOD,PRECON,NORM,WEIGHT,ITERM,N,M,TOL,MAXITN, &
ANORM,SIGMAX,MONIT,WORK,LWORK,LWREQ,IFAIL)```
f11bdf contains two additional arguments as follows:
See the routine document for further information.

### f11bbf

Withdrawn at Mark 21.
Replaced by f11bef.

#### Old Code

`CALL f11bbf(IREVCM,U,V,WORK,LWORK,IFAIL)`

#### New Code

`CALL f11bef(IREVCM,U,V,WGT,WORK,LWORK,IFAIL)`
wgt must be a one-dimensional real array of length at least $n$ (the order of the matrix) if weights are to be used in the termination criterion, and $1$ otherwise. Note that the call to f11bef requires the weights to be supplied in ${\mathbf{wgt}}\left(1:n\right)$ rather than ${\mathbf{work}}\left(1:n\right)$. The minimum value of the argument lwork may also need to be changed.

### f11bcf

Withdrawn at Mark 21.
Replaced by f11bff.

#### Old Code

`CALL f11bcf(ITN,STPLHS,STPRHS,ANORM,SIGMAX,IFAIL)`

#### New Code

`CALL f11bff(ITN,STPLHS,STPRHS,ANORM,SIGMAX,WORK,LWORK,IFAIL)`
f11bff contains two additional arguments as follows:
• ${\mathbf{work}}\left({\mathbf{lwork}}\right)$ – real array.
• lwork – integer.
See the routine document for further information.

### f11gaf

Withdrawn at Mark 22.
Replaced by f11gdf.

#### Old Code

```CALL f11gaf(METHOD,PRECON,SIGCMP,NORM,WEIGHT,ITERM,N,TOL,MAXITN,    &
ANORM,SIGMAX,SIGTOL,MAXITS,MONIT,LWREQ,IFAIL)```

#### New Code

```CALL f11gdf(METHOD,PRECON,SIGCMP,NORM,WEIGHT,ITERM,N,TOL,MAXITN,    &
ANORM,SIGMAX,SIGTOL,MAXITS,MONIT,LWREQ,WORK,LWORK,IFAIL)```
f11gdf contains two additional arguments as follows:
• ${\mathbf{work}}\left({\mathbf{lwork}}\right)$ – real array.
• lwork – integer.
See the routine document for further information.

### f11gbf

Withdrawn at Mark 22.
Replaced by f11gef.

#### Old Code

`CALL f11gbf(IREVCM,U,V,WORK,LWORK,IFAIL)`

#### New Code

`CALL f11gef(IREVCM,U,V,WGT,WORK,LWORK,IFAIL)`
wgt must be a one-dimensional real array of length at least $n$ (the order of the matrix) if weights are to be used in the termination criterion, and $1$ otherwise. Note that the call to f11gef requires the weights to be supplied in ${\mathbf{wgt}}\left(1:n\right)$ rather than ${\mathbf{work}}\left(1:n\right)$. The minimum value of the argument lwork may also need to be changed.

### f11gcf

Withdrawn at Mark 22.
Replaced by f11gff.

#### Old Code

`CALL f11gcf(ITN,STPLHS,STPRHS,ANORM,SIGMAX,ITS,SIGERR,IFAIL)`

#### New Code

```CALL f11gff(ITN,STPLHS,STPRHS,ANORM,SIGMAX,ITS,SIGERR,       &
WORK,LWORK,IFAIL)```
f11gff contains two additional arguments as follows:
• ${\mathbf{work}}\left({\mathbf{lwork}}\right)$ – real array.
• lwork – integer.
See the routine document for further information.

## G01 – Simple Calculations on Statistical Data

### g01aaf

Withdrawn at Mark 26.
Replaced by g01atf.
Withdrawn because on output, additional information was needed to allow large datasets to be processed in blocks and the results combined through a call to g01auf. This information is returned in rcomm.

#### Old Code

`CALL g01aaf(N,X,IWT,WT,XMEAN,S2,S3,S4,XMIN,XMAX,WTSUM,IFAIL)`

#### New Code

```PN = 0
CALL g01atf(N,X,IWT,WT,PN,XMEAN,S2,S3,S4,XMIN,XMAX,RCOMM,IFAIL)
IWT = PN
WTSUM = RCOMM(1)```

### g01agf

Withdrawn at Mark 27.
There is no replacement for this routine.

### g01ahf

Withdrawn at Mark 27.
There is no replacement for this routine.

### g01ajf

Withdrawn at Mark 27.
There is no replacement for this routine.

## G02 – Correlation and Regression Analysis

### g02jafg02jbfg02jcfg02jdfg02jef

Deprecated.
For each of the routines g02jaf, g02jbf, g02jcf, g02jdf and g02jef there is not a direct one-to-one replacement, rather they have been replaced with a new suite of routines. This new suite allows a linear mixed effects model to be specified using a modelling language; giving a more natural way of specifying the model, allowing interaction terms to be specified and means that it is no longer necessary to create dummy variables when the model contains categorical variables.
The new suite of routines consists of:
• g22ybf used to describe the dataset of interest. Calling this routine allows labels to be assign to variables, which can then be used when specifying the model.
• g22yaf multiple calls to this routine are used to specify the fixed and random part of the model. The model is specified using strings and a modelling language, for example the string: ${V}_{1}+{V}_{2}+{V}_{1}.{V}_{1}$ would specify a model with the main effects of variables ${V}_{1}$ and ${V}_{2}$ and the interaction term between them. The modelling language is explained in detail in Section 3 in g22yaf.
• g02jff pre-processes the dataset prior to calling the model fitting routine.
• g02jhf fits the model and returns the parameter estimates etc.
In addition to the routines listed above, the following can also be used:
• g02jgf combines information returned by multiple calls to g02jff. This is useful for large problems as it allows the dataset to be split up into smaller subsets of data, pre-processing each one separately before combining them into a single set of information as if g02jff had been called on the full dataset.
• g22ydf can be used to obtain labels for the parameter estimates returned by g02jhf.
• g22zmf can be used to set any optional arguments.
• g22znf can be used to return the value of any optional arguments.
By default, the model fitting routine, g02jhf, fits the linear mixed effects model using restricted maximum likelihood (REML). In order to fit the model using maximum likelihood (ML) you need to call the optional argument setting routine, g22zmf with optstr set to $\mathit{Likelihood}=\mathit{ML}$, between the call to g02jff and the call to g02jhf.

## G05 – Random Number Generators

### g05caf

Withdrawn at Mark 22.
Replaced by g05saf.

#### Old Code

```DO 20 I = 1, N
X(I) = g05caf(X(I))
20 CONTINUE```

#### New Code

`CALL g05saf(N,STATE,X,IFAIL)`
The integer array state in the call to g05saf contains information on the base generator being used. This array must have been initialized prior to calling g05saf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization. Due to changes in the underlying code the sequence of values produced by g05saf is likely to be different from those produced by g05caf.

### g05cbf

Withdrawn at Mark 22.
Replaced by g05kff.

#### Old Code

`CALL g05cbf(I)`

#### New Code

```LSEED = 1
SEED(1) = I
GENID = 1
SUBID = 1
CALL g05kff(GENID,SUBID,SEED,LSEED,STATE,LSTATE,IFAIL)```
The integer array state in the call to g05kff contains information on the base generator being used. The base generator is chosen via the integer arguments genid and subid. The required length of the array state depends on the base generator chosen. Due to changes in the underlying code a sequence of values produced by using a random number generator initialized via a call to g05kff is likely to be different from a sequence produced by a generator initialized by g05cbf, even if the same value for I is used.

### g05ccf

Withdrawn at Mark 22.
Replaced by g05kgf.

#### Old Code

`CALL g05ccf`

#### New Code

```GENID = 1
SUBID = 1
CALL g05kgf(GENID,SUBID,STATE,LSTATE,IFAIL)```
The integer array state in the call to g05kgf contains information on the base generator being used. The base generator is chosen via the integer arguments genid and subid. The required length of the array state depends on the base generator chosen.

### g05cff

Withdrawn at Mark 22.
Replaced by f06dff.

#### Old Code

`CALL g05cff(IA,NI,XA,NX,IFAIL)`

#### New Code

```LSTATE = STATE(1)
CALL f06dff(LSTATE,STATE,1,CSTATE,1)```
The state of the base generator for the group of routines g05kff, g05kgf, g05khf, g05kjf, g05ncf, g05ndf, g05pdfg05pzf, g05rcfg05rzf, G05S and G05T can be saved by simply creating a local copy of the array state. The first element of the state array contains the number of elements that are used by the random number generating routines, therefore either this number of elements can be copied, or the whole array (as defined in the calling program).

### g05cgf

Withdrawn at Mark 22.
Replaced by f06dff.

#### Old Code

`CALL g05cgf(IA,NI,XA,NX,IFAIL)`

#### New Code

```LSTATE = CSTATE(1)
CALL f06dff(LSTATE,CSTATE,1,STATE,1)```
The state of the base generator for the group of routines g05kff, g05kgf, g05khf, g05kjf, g05ncf, g05ndf, g05pdfg05pzf, g05rcfg05rzf, G05S and G05T can be restored by simply copying back the previously saved copy of the state array. The first element of the state array contains the number of elements that are used by the random number generating routines, therefore either this number of elements can be copied, or the whole array (as defined in the calling program).

### g05daf

Withdrawn at Mark 22.
Replaced by g05sqf.

#### Old Code

```DO 10 I = 1, N
X(I) = g05daf(AA,BB)
10 CONTINUE```

#### New Code

```A = MIN(AA,BB)
B = MAX(AA,BB)
IFAIL = 0
CALL g05sqf(N,A,B,STATE,X,IFAIL)```
The old routine g05daf returns a single variate at a time, whereas the new routine g05sqf returns a vector of n values in one go. In g05sqf the minimum value must be held in the argument a and the maximum in argument b, therefore ${\mathbf{a}}<{\mathbf{b}}$. This was not the case for the equivalent arguments in g05daf.
The integer array state in the call to g05sqf contains information on the base generator being used. This array must have been initialized prior to calling g05sqf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization. Due to changes in the underlying code the sequence of values produced by g05sqf is likely to be different from those produced by g05daf.

### g05dbf

Withdrawn at Mark 22.
Replaced by g05sff.

#### Old Code

```DO 10 I = 1, N
X(I) = g05dbf(AA)
10 CONTINUE```

#### New Code

```A = ABS(AA)
IFAIL = 0
CALL g05sff(N,A,STATE,X,IFAIL)```
The old routine g05dbf returns a single variate at a time, whereas the new routine g05sff returns a vector of n values in one go. In g05sff argument a must be non-negative, this was not the case for the equivalent argument in g05dbf.
The integer array state in the call to g05sff contains information on the base generator being used. This array must have been initialized prior to calling g05sff with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization. Due to changes in the underlying code the sequence of values produced by g05sff is likely to be different from those produced by g05dbf.

### g05dcf

Withdrawn at Mark 22.
Replaced by g05slf.

#### Old Code

```DO 10 I = 1, N
X(I) = g05dcf(A,BB)
10 CONTINUE```

#### New Code

```B = ABS(BB)
IFAIL = 0
CALL g05slf(N,A,B,STATE,X,IFAIL)```
The old routine g05dcf returns a single variate at a time, whereas the new routine g05slf returns a vector of n values in one go. In g05slf the spread (argument a) must be positive, this was not the case for the equivalent arguments in g05dcf.
The integer array state in the call to g05slf contains information on the base generator being used. This array must have been initialized prior to calling g05slf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization. Due to changes in the underlying code the sequence of values produced by g05slf is likely to be different from those produced by g05dcf.

### g05ddf

Withdrawn at Mark 22.
Replaced by g05skf.

#### Old Code

```DO 10 I = 1, N
X(I) = g05ddf(XMU,SD)
10 CONTINUE```

#### New Code

```VAR = SD**2
IFAIL = 0
CALL g05skf(N,XMU,VAR,STATE,X,IFAIL)```
The old routine g05ddf returns a single variate at a time, whereas the new routine g05skf returns a vector of n values in one go. g05skf expects the variance of the Normal distribution (argument var), compared to g05ddf which expected the standard deviation.
The integer array state in the call to g05skf contains information on the base generator being used. This array must have been initialized prior to calling g05skf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization. Due to changes in the underlying code the sequence of values produced by g05skf is likely to be different from those produced by g05ddf.

### g05def

Withdrawn at Mark 22.
Replaced by g05smf.

#### Old Code

```DO 10 I = 1, N
X(I) = g05def(XMU,SD)
10 CONTINUE```

#### New Code

```VAR = SD**2
IFAIL = 0
CALL g05smf(N,XMU,VAR,STATE,X,IFAIL)```
The old routine g05def returns a single variate at a time, whereas the new routine g05smf returns a vector of n values in one go. g05smf expects the variance of the corresponding Normal distribution (argument var), compared to g05def which expected the standard deviation.
The integer array state in the call to g05smf contains information on the base generator being used. This array must have been initialized prior to calling g05smf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization. Due to changes in the underlying code the sequence of values produced by g05smf is likely to be different from those produced by g05def.

### g05dff

Withdrawn at Mark 22.
Replaced by g05scf.

#### Old Code

```DO 10 I = 1, N
X(I) = g05dff(XMED,B)
10 CONTINUE```

#### New Code

```SEMIQR = ABS(B)
IFAIL = 0
CALL g05scf(N,XMED,SEMIQR,STATE,X,IFAIL)```
The old routine g05dff returns a single variate at a time, whereas the new routine g05scf returns a vector of n values in one go. g05scf expects the semi-interquartile range (argument semiqr) to be non-negative, this was not the case for the equivalent argument in g05dff.
The integer array state in the call to g05scf contains information on the base generator being used. This array must have been initialized prior to calling g05scf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization. Due to changes in the underlying code the sequence of values produced by g05scf is likely to be different from those produced by g05dff.

### g05dhf

Withdrawn at Mark 22.
Replaced by g05sdf.

#### Old Code

```DO 10 I = 1, N
X(I) = g05dhf(DF,IFAIL)
10 CONTINUE```

#### New Code

`CALL g05sdf(N,DF,STATE,X,IFAIL)`
The old routine g05dhf returns a single variate at a time, whereas the new routine g05sdf returns a vector of n values in one go.
The integer array state in the call to g05sdf contains information on the base generator being used. This array must have been initialized prior to calling g05sdf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization. Due to changes in the underlying code the sequence of values produced by g05sdf is likely to be different from those produced by g05dhf.

### g05djf

Withdrawn at Mark 22.
Replaced by g05snf.

#### Old Code

```DO 10 I = 1, N
X(I) = g05djf(DF,IFAIL)
10 CONTINUE```

#### New Code

`CALL g05snf(N,DF,STATE,X,IFAIL)`
The old routine g05djf returns a single variate at a time, whereas the new routine g05snf returns a vector of n values in one go.
The integer array state in the call to g05snf contains information on the base generator being used. This array must have been initialized prior to calling g05snf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization. Due to changes in the underlying code the sequence of values produced by g05snf is likely to be different from those produced by g05djf.

### g05dkf

Withdrawn at Mark 22.
Replaced by g05shf.

#### Old Code

```DO 10 I = 1, N
X(I) = g05dkf(DF1,DF2,IFAIL)
10 CONTINUE```

#### New Code

`CALL g05shf(N,DF1,DF2,STATE,X,IFAIL)`
The old routine g05dkf returns a single variate at a time, whereas the new routine g05shf returns a vector of n values in one go.
The integer array state in the call to g05shf contains information on the base generator being used. This array must have been initialized prior to calling g05shf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization. Due to changes in the underlying code the sequence of values produced by g05shf is likely to be different from those produced by g05dkf.

### g05dpf

Withdrawn at Mark 22.
Replaced by g05ssf.

#### Old Code

```DO 10 I = 1, N
X(I) = g05dpf(A,B,IFAIL)
10 CONTINUE```

#### New Code

`CALL g05ssf(N,A,B,STATE,X,IFAIL)`
The old routine g05dpf returns a single variate at a time, whereas the new routine g05ssf returns a vector of n values in one go.
The integer array state in the call to g05ssf contains information on the base generator being used. This array must have been initialized prior to calling g05ssf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization. Due to changes in the underlying code the sequence of values produced by g05ssf is likely to be different from those produced by g05dpf.

### g05drf

Withdrawn at Mark 22.
Replaced by g05tkf.

#### Old Code

```DO 10 I = 1, N
X(I) = g05drf(LAMDA,IFAIL)
10 CONTINUE```

#### New Code

```MODE = 3
CALL g05tjf(MODE,N,LAMBDA,R,LR,STATE,X,IFAIL)```
The old routine g05drf returns a single variate at a time, whereas the new routine g05tjf returns a vector of n values in one go. For efficiency, the new routine can make use of a reference vector, r. If, as in this case, the integer argument mode is set to 3, the real reference vector r is not referenced, and its length, lr, need only be at least one.
The integer array state in the call to g05tjf contains information on the base generator being used. This array must have been initialized prior to calling g05tjf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization. Due to changes in the underlying code the sequence of values produced by g05tjf is likely to be different from those produced by g05drf.

### g05dyf

Withdrawn at Mark 22.
Replaced by g05tlf.

#### Old Code

```DO 10 I = 1, N
X(I) = g05dyf(AA,BB)
10 CONTINUE```

#### New Code

```IFAIL = 0
A = MIN(AA,BB)
B = MAX(AA,BB)
CALL g05tlf(N,A,B,STATE,X,IFAIL)```
The old routine g05dyf returns a single variate at a time, whereas the new routine g05tlf returns a vector of n values in one go. In g05tlf the minimum value must be held in the argument a and the maximum in argument b, therefore ${\mathbf{a}}\le {\mathbf{b}}$. This was not the case for the equivalent arguments in g05dyf.
The integer array state in the call to g05tlf contains information on the base generator being used. This array must have been initialized prior to calling g05tlf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization. Due to changes in the underlying code the sequence of values produced by g05tlf is likely to be different from those produced by g05dyf.

### g05dzf

Withdrawn at Mark 22.
Replaced by g05tbf.

#### Old Code

```DO 20 I = 1, N
X(I) = g05dzf(PP)
20 CONTINUE```

#### New Code

```P = MAX(0.0D0,MIN(PP,1.0D0))
IFAIL = 0
CALL g05tbf(N,P,STATE,X,IFAIL)```
The old routine g05dzf returns a single variate at a time, whereas the new routine g05tbf returns a vector of n values in one go. The real argument p in g05tbf must not be less than zero or greater than one, this was not the case for the equivalent argument in g05dzf.
The integer array state in the call to g05tbf contains information on the base generator being used. This array must have been initialized prior to calling g05tbf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization. Due to changes in the underlying code the sequence of values produced by g05tbf is likely to be different from those produced by g05dzf.

### g05eaf

Withdrawn at Mark 22.
Replaced by g05rzf.

#### Old Code

`CALL g05eaf(XMU,M,C,LDC,EPS,R1,LR1,IFAIL)`

#### New Code

```MODE = 0
CALL g05rzf(MODE,N,M,XMU,C,LDC,R,LR,STATE,X,LDX,IFAIL)```
The old routine g05eaf sets up a reference vector for use by g05ezf. The functionality of both these routines has been combined into the single new routine g05rzf. Setting ${\mathbf{mode}}=0$ in the call to g05rzf only sets up the real reference vector r and hence mimics the functionality of g05eaf.
The length of the real reference vector, r, in g05rzf must be at least ${\mathbf{m}}×\left({\mathbf{m}}+1\right)+1$. In contrast to the equivalent argument in g05eaf, this array must be allocated in the calling program.

### g05ebf

Withdrawn at Mark 22.
Replaced by g05tlf.
There is no direct replacement for routine g05ebf. g05ebf sets up a reference vector for use by g05eyf, this reference vector is no longer required. The replacement routine for g05eyf is g05tlf.

### g05ecf

Withdrawn at Mark 22.
Replaced by g05tjf.

#### Old Code

```CALL g05ecf(LAMBDA,R1,LR1,IFAIL)
DO 10 I = 1, N
X(I) = g05eyf(R1,LR1)
10 CONTINUE```

#### New Code

```MODE = 2
CALL g05tjf(MODE,N,LAMBDA,R,LR,STATE,X,IFAIL)```
The old routine g05ecf sets up a reference vector for use by g05eyf. The replacement routine g05tjf is now used to both set up a reference vector and generate the required variates. Setting ${\mathbf{mode}}=0$ in the call to g05tjf sets up the real reference vector r and hence mimics the functionality of g05ecf. Setting ${\mathbf{mode}}=1$ generates a series of variates from a reference vector mimicking the functionality of g05eyf for this particular distribution. Setting ${\mathbf{mode}}=2$ initializes the reference vector and generates the variates in one go.
The routine g05eyf returns a single variate at a time, whereas the new routine g05tjf returns a vector of n values in one go.
The length of the real reference vector, r, in g05tjf, must be allocated in the calling program in contrast to the equivalent argument in g05ecf, see the documentation for more details.
The integer array state in the call to g05tjf contains information on the base generator being used. This array must have been initialized prior to calling g05tjf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization. Due to changes in the underlying code the sequence of values produced by g05tjf is likely to be different from those produced by a combination of g05ecf and g05eyf.

### g05edf

Withdrawn at Mark 22.
Replaced by g05taf.

#### Old Code

```CALL g05edf(M,P,R1,LR1,IFAIL)
DO 10 I = 1, N
X(I) = g05eyf(R1,LR1)
10 CONTINUE```

#### New Code

```MODE = 2
CALL g05taf(MODE,N,M,P,R,LR,STATE,X,IFAIL)```
The old routine g05edf sets up a reference vector for use by g05eyf. The replacement routine g05taf is now used to both set up a reference vector and generate the required variates. Setting ${\mathbf{mode}}=0$ in the call to g05taf sets up the real reference vector r and hence mimics the functionality of g05edf. Setting ${\mathbf{mode}}=1$ generates a series of variates from a reference vector mimicking the functionality of g05eyf for this particular distribution. Setting ${\mathbf{mode}}=2$ initializes the reference vector and generates the variates in one go.
The routine g05eyf returns a single variate at a time, whereas the new routine g05taf returns a vector of n values in one go.
The length of the real reference vector, r, in g05taf, needs to be a different length from the equivalent argument in g05edf, see the documentation for more details.
The integer array state in the call to g05taf contains information on the base generator being used. This array must have been initialized prior to calling g05taf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization. Due to changes in the underlying code the sequence of values produced by g05taf is likely to be different from those produced by a combination of g05edf and g05eyf.

### g05eef

Withdrawn at Mark 22.
Replaced by g05thf.

#### Old Code

```CALL g05eef(M,P,R1,LR1,IFAIL)
DO 10 I = 1, N
X(I) = g05eyf(R1,LR1)
10 CONTINUE```

#### New Code

```MODE = 2
CALL g05thf(MODE,N,M,P,R,LR,STATE,X,IFAIL)```
The old routine g05eef sets up a reference vector for use by g05eyf. The replacement routine g05thf is now used to both set up a reference vector and generate the required variates. Setting ${\mathbf{mode}}=0$ in the call to g05thf sets up the real reference vector r and hence mimics the functionality of g05eef. Setting ${\mathbf{mode}}=1$ generates a series of variates from a reference vector mimicking the functionality of g05eyf for this particular distribution. Setting ${\mathbf{mode}}=2$ initializes the reference vector and generates the variates in one go.
The routine g05eyf returns a single variate at a time, whereas the new routine g05thf returns a vector of n values in one go.
The length of the real reference vector, r, in g05thf, needs to be a different length from the equivalent argument in g05eef, see the documentation for g05thf for more details.
The integer array state in the call to g05thf contains information on the base generator being used. This array must have been initialized prior to calling g05thf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization. Due to changes in the underlying code the sequence of values produced by g05thf is likely to be different from those produced by a combination of g05eef and g05eyf.

### g05eff

Withdrawn at Mark 22.
Replaced by g05tef.

#### Old Code

```CALL g05eff(NS,M,NP,R1,LR1,IFAIL)
DO 10 I = 1, N
X(I) = g05eyf(R1,LR1)
10 CONTINUE```

#### New Code

```MODE = 2
CALL g05tef(MODE,N,NS,NP,M,R,LR,STATE,X,IFAIL)```
The old routine g05eff sets up a reference vector for use by g05eyf. The replacement routine g05tef is now used to both set up a reference vector and generate the required variates. Setting ${\mathbf{mode}}=0$ in the call to g05tef sets up the real reference vector r and hence mimics the functionality of g05eff. Setting ${\mathbf{mode}}=1$ generates a series of variates from a reference vector mimicking the functionality of g05eyf for this particular distribution. Setting ${\mathbf{mode}}=2$ initializes the reference vector and generates the variates in one go.
The routine g05eyf returns a single variate at a time, whereas the new routine g05tef returns a vector of n values in one go.
The length of the real reference vector, r, in g05tef, needs to be a different length from the equivalent argument in g05eff, see the documentation for more details.
The integer array state in the call to g05tef contains information on the base generator being used. This array must have been initialized prior to calling g05tef with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization. Due to changes in the underlying code the sequence of values produced by g05tef is likely to be different from those produced by a combination of g05eff and g05eyf.

### g05egf

Withdrawn at Mark 22.
Replaced by g05phf.

#### Old Code

`CALL g05egf(E,A,NA,B,NB,R,NR,VAR,IFAIL)`

#### New Code

```AVAR = B(1)**2
IQ = NB - 1
IF (AVAR.GT.0.0D0) THEN
DO 10 I = 1, IQ
THETA(I) = -B(I+1)/B(1)
10    CONTINUE
ELSE
DO 20 I = 1, IQ
THETA(I) = 0.0D0
20    CONTINUE
END IF
MODE = 0
CALL g05phf(MODE,N,E,NA,A,IQ,THETA,AVAR,R,LR,STATE,VAR,X,IFAIL)```
The real vector theta must be of length at least ${\mathbf{iq}}=\mathbf{nb}-1$.
The old routine g05egf sets up a reference vector for use by g05ewf. The replacement routine g05phf is now used to both set up a reference vector and generate the required variates. Setting ${\mathbf{mode}}=0$ in the call to g05phf sets up the real reference vector r and hence mimics the functionality of g05egf. When ${\mathbf{mode}}=0$, the integer array state in the call to g05phf need not be set.

### g05ehf

Withdrawn at Mark 22.
Replaced by g05ncf.

#### Old Code

`CALL g05ehf(INDEX,N,IFAIL)`

#### New Code

`CALL g05ncf(INDEX,N,STATE,IFAIL)`
The integer array state in the call to g05ncf contains information on the base generator being used. This array must have been initialized prior to calling g05ncf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization. Due to changes in the underlying code the sequence of values produced by g05ncf is likely to be different from those produced by g05ehf.

### g05ejf

Withdrawn at Mark 22.
Replaced by g05ndf.

#### Old Code

`CALL g05ejf(IA,N,IZ,M,IFAIL)`

#### New Code

`CALL g05ndf(IA,N,IZ,M,STATE,IFAIL)`
The integer array state in the call to g05ndf contains information on the base generator being used. This array must have been initialized prior to calling g05ndf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization. Due to changes in the underlying code the sequence of values produced by g05ndf is likely to be different from those produced by g05ejf.

### g05ewf

Withdrawn at Mark 22.
Replaced by g05phf.

#### Old Code

```CALL g05egf(E,A,NA,B,NB,R,NR,VAR,IFAIL)
DO 10 I = 1, N
X(I) = g05ewf(R,NR,IFAIL)
10 CONTINUE```

#### New Code

```AVAR = B(1)**2
IQ = NB - 1
IF (AVAR.GT.0.0D0) THEN
DO 10 I = 1, IQ
THETA(I) = -B(I+1)/B(1)
10    CONTINUE
ELSE
DO 20 I = 1, IQ
THETA(I) = 0.0D0
20    CONTINUE
END IF
MODE = 2
CALL g05phf(MODE,N,E,NA,A,NB-1,THETA,AVAR,VAR,R,LR,STATE,X,IFAIL)```
The real vector theta must be of length at least ${\mathbf{iq}}=\mathbf{nb}-1$.
The old routine g05egf sets up a reference vector for use by g05ewf. The replacement routine g05phf is now used to both set up a reference vector and generate the required variates. Setting the integer argument mode to 0 in the call to g05phf sets up the real reference vector r and hence mimics the functionality of g05egf. Setting mode to 1 generates a series of variates from a reference vector mimicking the functionality of g05ewf. Setting mode to 2 initializes the reference vector and generates the variates in one go.
The integer array state in the call to g05phf contains information on the base generator being used. This array must have been initialized prior to calling g05phf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization. Due to changes in the underlying code the sequence of values produced by g05phf is likely to be different from those produced by g05egf.

### g05exf

Withdrawn at Mark 22.
Replaced by g05tdf.

#### Old Code

```CALL g05exf(P,NP,IP1,ITYPE,R1,LR1,IFAIL)
DO 10 I = 1, N
X(I) = g05eyf(R1,LR1)
10 CONTINUE```

#### New Code

```MODE = 2
CALL g05tdf(MODE,N,P,NP,IP1,ITYPE,R,LR,STATE,X,IFAIL)```
The old routine g05exf sets up a reference vector for use by g05eyf. The replacement routine g05tdf is now used to both set up a reference vector and generate the required variates. Setting ${\mathbf{mode}}=0$ in the call to g05tdf sets up the real reference vector r and hence mimics the functionality of g05exf. Setting ${\mathbf{mode}}=1$ generates a series of variates from a reference vector mimicking the functionality of g05eyf for this particular distribution. Setting ${\mathbf{mode}}=2$ initializes the reference vector and generates the variates in one go.
The routine g05eyf returns a single variate at a time, whereas the new routine g05tdf returns a vector of n values in one go.
The length of the real reference vector, r, in g05tdf must be allocated in the calling program in contrast to the equivalent argument in g05exf, see the documentation for more details.
The integer array state in the call to g05tdf contains information on the base generator being used. This array must have been initialized prior to calling g05tdf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization. Due to changes in the underlying code the sequence of values produced by g05tdf is likely to be different from those produced by a combination of g05exf and g05eyf.

### g05eyf

Withdrawn at Mark 22.
Replaced by g05tdf.
There is no direct replacement routine for g05eyf.
g05eyf is designed to generate random draws from a distribution defined by a reference vector. These reference vectors are created by other routines in Chapter G05, for example g05ebf, which have themselves been superseded. In order to replace a call to g05eyf you must identify which NAG routine generated the reference vector being used and look up its replacement. For example, to replace a call to g05eyf preceded by a call to g05ebf, as in:
```    CALL g05ebf(M,IB,R,NR,IFAIL)
X = g05eyf(R,NR)```
you would need to look at the replacement routine for g05ebf.

### g05ezf

Withdrawn at Mark 22.
Replaced by g05rzf.

#### Old Code

```CALL g05eaf(XMU,N,C,LDC,EPS,R1,LR1,IFAIL)
DO 20 I = 1, N
CALL g05ezf(CX,M,R,NR,IFAIL)
DO 30 J = 1, M
X(I,J) = CX(J)
30    CONTINUE
20 CONTINUE```

#### New Code

```MODE = 2
CALL g05rzf(MODE,N,M,XMU,C,LDC,R,LR,STATE,X,LDX,IFAIL)```
The old routine g05eaf sets up a reference vector for use by g05ezf. The functionality of both these routines has been combined into the single new routine g05rzf. Setting ${\mathbf{mode}}=2$ in the call to g05rzf sets up the real reference vector r and generates the draws from the multivariate Normal distribution in one go.
The old routine g05ezf returns a single (m-dimensional vector) draw from the multivariate Normal distribution at a time, whereas the new routine g05rzf returns an n by m matrix of n draws in one go.
The integer array state in the call to g05rzf contains information on the base generator being used. This array must have been initialized prior to calling g05rzf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization. Due to changes in the underlying code the sequence of values produced by g05rzf is likely to be different from those produced by g05ezf.

### g05faf

Withdrawn at Mark 22.
Replaced by g05sqf.

#### Old Code

`CALL g05faf(AA,BB,N,X)`

#### New Code

```A = MIN(AA,BB)
B = MAX(AA,BB)
IFAIL = 0
CALL g05sqf(N,A,B,STATE,X,IFAIL)```
In g05sqf the minimum value must be held in the argument A and the maximum in argument b, therefore ${\mathbf{a}}\le {\mathbf{b}}$. This was not the case for the equivalent arguments in g05faf.
The integer array state in the call to g05sqf contains information on the base generator being used. This array must have been initialized prior to calling g05sqf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization. Due to changes in the underlying code the sequence of values produced by g05sqf is likely to be different from those produced by g05faf.

### g05fbf

Withdrawn at Mark 22.
Replaced by g05sff.

#### Old Code

`CALL g05fbf(AA,N,X)`

#### New Code

```A = ABS(AA)
IFAIL = 0
CALL g05sff(N,A,STATE,X,IFAIL)```
In g05sff argument a must be non-negative, this was not the case for the equivalent argument in g05fbf.
The integer array state in the call to g05sff contains information on the base generator being used. This array must have been initialized prior to calling g05sff with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization. Due to changes in the underlying code the sequence of values produced by g05sff is likely to be different from those produced by g05fbf.

### g05fdf

Withdrawn at Mark 22.
Replaced by g05skf.

#### Old Code

`CALL g05fdf(XMU,SD,N,X)`

#### New Code

```VAR = SD**2
IFAIL = 0
CALL g05skf(N,XMU,VAR,STATE,X,IFAIL)```
g05skf expects the variance of the Normal distribution (argument var), compared to g05fdf which expected the standard deviation.
The integer array state in the call to g05skf contains information on the base generator being used. This array must have been initialized prior to calling g05skf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization. Due to changes in the underlying code the sequence of values produced by g05skf is likely to be different from those produced by g05fdf.

### g05fef

Withdrawn at Mark 22.
Replaced by g05sbf.

#### Old Code

`CALL g05fef(A,B,N,X,IFAIL)`

#### New Code

`CALL g05sbf(N,A,B,STATE,X,IFAIL)`
The integer array state in the call to g05sbf contains information on the base generator being used. This array must have been initialized prior to calling g05sbf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization. Due to changes in the underlying code the sequence of values produced by g05sbf is likely to be different from those produced by g05fef.

### g05fff

Withdrawn at Mark 22.
Replaced by g05sjf.

#### Old Code

`CALL g05fff(A,B,N,X,IFAIL)`

#### New Code

`CALL g05sjf(N,A,B,STATE,X,IFAIL)`
The integer array state in the call to g05sjf contains information on the base generator being used. This array must have been initialized prior to calling g05sjf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization. Due to changes in the underlying code the sequence of values produced by g05sjf is likely to be different from those produced by g05fff.

### g05fsf

Withdrawn at Mark 22.
Replaced by g05srf.

#### Old Code

`CALL g05fsf(VK,N,X,IFAIL)`

#### New Code

`CALL g05srf(N,VK,STATE,X,IFAIL)`
The integer array state in the call to g05srf contains information on the base generator being used. This array must have been initialized prior to calling g05srf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization. Due to changes in the underlying code the sequence of values produced by g05srf is likely to be different from those produced by g05fsf.

### g05gaf

Withdrawn at Mark 22.
Replaced by g05pxf.

#### Old Code

`CALL g05gaf(SIDE,INIT,M,N,A,LDA,WK,IFAIL)`

#### New Code

`CALL g05pxf(SIDE,INIT,M,N,STATE,A,LDA,IFAIL)`
The integer array state in the call to g05pxf contains information on the base generator being used. This array must have been initialized prior to calling g05pxf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization. Due to changes in the underlying code the sequence of values produced by g05pxf is likely to be different from those produced by g05gaf.

### g05gbf

Withdrawn at Mark 22.
Replaced by g05pyf.

#### Old Code

`CALL g05gbf(N,D,C,LDC,EPS,WK,IFAIL)`

#### New Code

`CALL g05pyf(N,D,EPS,STATE,C,LDC,IFAIL)`
The integer array state in the call to g05pyf contains information on the base generator being used. This array must have been initialized prior to calling g05pyf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization. Due to changes in the underlying code the sequence of values produced by g05pyf is likely to be different from those produced by g05gbf.

### g05hdf

Withdrawn at Mark 22.
Replaced by g05pjf.

#### Old Code

```CALL g05hdf(MODE,K,IP,IQ,MEAN,PAR,LPAR,QQ,LDQQ,N,W,REF,LREF,  &
IWORK,LIWORK,IFAIL)```

#### New Code

```IF (MODE.EQ.'S') THEN
IMODE = 0
ELSE IF (MODE.EQ.'C') THEN
IMODE = 1
ELSE IF (MODE.EQ.'R') THEN
IMODE = 3
END IF
LL = 0
DO 30 L = 1, IP
DO 20 I = 1, K
DO 10 J = 1, K
LL = LL + 1
PHI(I,J,L) = PAR(LL)
10       CONTINUE
20    CONTINUE
30 CONTINUE
DO 60 L = 1, IQ-1
DO 50 I = 1, K
DO 40 J = 1, K
LL = LL + 1
THETA(I,J,L) = PAR(LL)
40       CONTINUE
50    CONTINUE
60 CONTINUE
IF (MEAN.EQ.'M') THEN
DO 70 I = 1, K
LL = LL + 1
XMEAN(I) = PAR(LL)
70    CONTINUE
ELSE
DO 80 I = 1, K
XMEAN(I) = 0.0D0
80    CONTINUE
END IF
LDW = N
CALL g05pjf(IMODE,N,K,XMEAN,IP,PHI,IQ,THETA,QQ,LDQQ,REF,LREF, &
STATE,W,LDW,IWORK,LIWORK,IFAIL)```
The integer argument IMODE should be set to 0, 1 or 3 in place of the argument mode having settings of 'S', 'C' or 'R' respectively. The real array phi should have length at least $\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left(1,{\mathbf{ip}}×\left({\mathbf{k}}×{\mathbf{k}}\right)\right)$; if dimensioned as ${\mathbf{phi}}\left({\mathbf{k}},{\mathbf{k}},{\mathbf{ip}}\right)$ (as in the above example) then ${\mathbf{phi}}\left(i,j,l\right)$ will contain the element $\mathbf{par}\left(\left(l-1\right)×k×k+\left(i-1\right)×k+j\right)$. The real array theta should have length at least $\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left(1,{\mathbf{iq}}×\left({\mathbf{k}}×{\mathbf{k}}\right)\right)$; if dimensioned as ${\mathbf{theta}}\left({\mathbf{k}},{\mathbf{k}},{\mathbf{iq}}\right)$ (as in the above example) then ${\mathbf{theta}}\left(i,j,l\right)$ will contain the element $\mathbf{par}\left(\mathbf{ip}×k×k+\left(l-1\right)×k×k+\left(i-1\right)×k+j\right)$. The real array xmean should have length at least k; if $\mathbf{mean}=\text{'M'}$ then ${\mathbf{xmean}}\left(i\right)$ will contain the element $\mathbf{par}\left(\mathbf{ip}+\mathbf{iq}×k×k+i\right)$, otherwise xmean should contain an array of zero values.
The integer array state in the call to g05pjf contains information on the base generator being used. This array must have been initialized prior to calling g05pjf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization. Due to changes in the underlying code the sequence of values produced by g05pjf is likely to be different from those produced by g05hdf.

### g05hkf

Withdrawn at Mark 24.
Replaced by g05pdf.

#### Old Code

```CALL g05hkf(DIST,NUM,IP,IQ,THETA,GAMMA,DF,HT,ET,FCALL,RVEC,IGEN,  &
ISEED,RWSAV,IFAIL)```

#### New Code

```CALL g05pdf(DIST,NUM,IP,IQ,THETA,GAMMA,DF,HT,ET,FCALL,R,LR,STATE, &
IFAIL)```
The integer array state in the call to g05pdf contains information on the base generator being used. This array must have been initialized prior to calling g05pdf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization. Due to changes in the underlying code the sequence of values produced by g05pdf is likely to be different from those produced by g05hkf.

### g05hlf

Withdrawn at Mark 24.
Replaced by g05pef.

#### Old Code

```CALL g05hlf(DIST,NUM,IP,IQ,THETA,GAMMA,DF,HT,ET,FCALL,RVEC,IGEN,  &
ISEED,RWSAV,IFAIL)```

#### New Code

```CALL g05pef(DIST,NUM,IP,IQ,THETA,GAMMA,DF,HT,ET,FCALL,R,LR,STATE, &
IFAIL)```
The integer array state in the call to g05pef contains information on the base generator being used. This array must have been initialized prior to calling g05pef with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization. Due to changes in the underlying code the sequence of values produced by g05pef is likely to be different from those produced by g05hlf.

### g05hmf

Withdrawn at Mark 24.
Replaced by g05pff.

#### Old Code

```CALL g05hmf(DIST,NUM,IP,IQ,THETA,GAMMA,DF,HT,ET,FCALL,RVEC,IGEN,  &
ISEED,RWSAV,IFAIL)```

#### New Code

```CALL g05pff(DIST,NUM,IP,IQ,THETA,GAMMA,DF,HT,ET,FCALL,R,LR,STATE, &
IFAIL)```
The integer array state in the call to g05pff contains information on the base generator being used. This array must have been initialized prior to calling g05pff with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization. Due to changes in the underlying code the sequence of values produced by g05pff is likely to be different from those produced by g05hmf.

### g05hnf

Withdrawn at Mark 24.
Replaced by g05pgf.

#### Old Code

```CALL g05hnf(DIST,NUM,IP,IQ,THETA,DF,HT,ET,FCALL,RVEC,IGEN,ISEED, &
RWSAV,IFAIL)```

#### New Code

```CALL g05pgf(DIST,NUM,IP,IQ,THETA,DF,HT,ET,FCALL,RVEC,STATE,      &
IFAIL)```
The integer array state in the call to g05pgf contains information on the base generator being used. This array must have been initialized prior to calling g05pgf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization.

### g05kaf

Withdrawn at Mark 24.
Replaced by g05saf.

#### Old Code

```DO 20 I = 1, N
X(I) = g05kaf(IGEN,ISEED)
20 CONTINUE```

#### New Code

`CALL g05saf(N,STATE,X,IFAIL)`
The old routine g05kaf returns a single variate at a time, whereas the new routine g05saf returns a vector of n values in one go.
The integer array state in the call to g05saf contains information on the base generator being used. This array must have been initialized prior to calling g05saf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization.

### g05kbf

Withdrawn at Mark 24.
Replaced by g05kff.

#### Old Code

`g05kbf(IGEN,ISEED)`

#### New Code

```IF (IGEN.EQ.0) THEN
CALL g05kff(1,1,ISEED,LSEED,STATE,LSTATE,IFAIL)
ELSE
CALL g05kff(2,IGEN,ISEED,LSEED,STATE,LSTATE,IFAIL)
END IF```

### g05kcf

Withdrawn at Mark 24.
Replaced by g05kgf.

#### Old Code

`CALL g05kcf(IGEN,ISEED)`

#### New Code

```IF (IGEN.EQ.0) THEN
CALL g05kgf(1,1,STATE,LSTATE,IFAIL)
ELSE
CALL g05kgf(2,IGEN,STATE,LSTATE,IFAIL)
END IF```

### g05kef

Withdrawn at Mark 24.
Replaced by g05tbf.

#### Old Code

```DO 20 I = 1, N
X(I) = g05kef(P,IGEN,ISEED,IFAIL)
20 CONTINUE```

#### New Code

`CALL g05tbf(N,P,STATE,X,IFAIL)`
The old routine g05kef returns a single variate at a time, whereas the new routine g05tbf returns a vector of n values in one go.
The integer array state in the call to g05tbf contains information on the base generator being used. This array must have been initialized prior to calling g05tbf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization.

### g05laf

Withdrawn at Mark 24.
Replaced by g05skf.

#### Old Code

`CALL g05laf(XMU,VAR,N,X,IGEN,ISEED,IFAIL)`

#### New Code

`CALL g05skf(N,XMU,VAR,STATE,X,IFAIL)`
The integer array state in the call to g05skf contains information on the base generator being used. This array must have been initialized prior to calling g05skf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization.

### g05lbf

Withdrawn at Mark 24.
Replaced by g05snf.

#### Old Code

`CALL g05lbf(DF,N,X,IGEN,ISEED,IFAIL)`

#### New Code

`CALL g05snf(N,DF,STATE,X,IFAIL)`
The integer array state in the call to g05snf contains information on the base generator being used. This array must have been initialized prior to calling g05snf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization.

### g05lcf

Withdrawn at Mark 24.
Replaced by g05sdf.

#### Old Code

`CALL g05lcf(DF,N,X,IGEN,ISEED,IFAIL)`

#### New Code

`CALL g05sdf(N,DF,STATE,X,IFAIL)`
The integer array state in the call to g05sdf contains information on the base generator being used. This array must have been initialized prior to calling g05sdf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization.

### g05ldf

Withdrawn at Mark 24.
Replaced by g05shf.

#### Old Code

`CALL g05ldf(DF1,DF2,N,X,IGEN,ISEED,IFAIL)`

#### New Code

`CALL g05shf(N,DF1,DF2,STATE,X,IFAIL)`
The integer array state in the call to g05shf contains information on the base generator being used. This array must have been initialized prior to calling g05shf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization.

### g05lef

Withdrawn at Mark 24.
Replaced by g05sbf.

#### Old Code

`CALL g05lef(A,B,N,X,IGEN,ISEED,IFAIL)`

#### New Code

`CALL g05sbf(N,A,B,STATE,X,IFAIL)`
The integer array state in the call to g05sbf contains information on the base generator being used. This array must have been initialized prior to calling g05sbf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization.

### g05lff

Withdrawn at Mark 24.
Replaced by g05sjf.

#### Old Code

`CALL g05lff(A,B,N,X,IGEN,ISEED,IFAIL)`

#### New Code

`CALL g05sjf(N,A,B,STATE,X,IFAIL)`
The integer array state in the call to g05sjf contains information on the base generator being used. This array must have been initialized prior to calling g05sjf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization.

### g05lgf

Withdrawn at Mark 24.
Replaced by g05sqf.

#### Old Code

`CALL g05lgf(A,B,N,X,IGEN,ISEED,IFAIL)`

#### New Code

`CALL g05sqf(N,A,B,STATE,X,IFAIL)`
The integer array state in the call to g05sqf contains information on the base generator being used. This array must have been initialized prior to calling g05sqf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization.

### g05lhf

Withdrawn at Mark 24.
Replaced by g05spf.

#### Old Code

`CALL g05lhf(XMIN,XMAX,XMED,N,X,IGEN,ISEED,IFAIL)`

#### New Code

`CALL g05spf(N,XMIN,XMED,XMAX,STATE,X,IFAIL)`
The integer array state in the call to g05spf contains information on the base generator being used. This array must have been initialized prior to calling g05spf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization.

### g05ljf

Withdrawn at Mark 24.
Replaced by g05sff.

#### Old Code

`CALL g05ljf(A,N,X,IGEN,ISEED,IFAIL)`

#### New Code

`CALL g05sff(N,A,STATE,X,IFAIL)`
The integer array state in the call to g05sff contains information on the base generator being used. This array must have been initialized prior to calling g05sff with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization.

### g05lkf

Withdrawn at Mark 24.
Replaced by g05smf.

#### Old Code

`CALL g05lkf(XMU,VAR,N,X,IGEN,ISEED,IFAIL)`

#### New Code

`CALL g05smf(N,XMU,VAR,STATE,X,IFAIL)`
The integer array state in the call to g05smf contains information on the base generator being used. This array must have been initialized prior to calling g05smf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization.

### g05llf

Withdrawn at Mark 24.
Replaced by g05sjf.

#### Old Code

`CALL g05llf(XMED,SEMIQR,N,X,IGEN,ISEED,IFAIL)`

#### New Code

`CALL g05scf(N,XMED,SEMIQR,STATE,X,IFAIL)`
The integer array state in the call to g05scf contains information on the base generator being used. This array must have been initialized prior to calling g05scf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization.

### g05lmf

Withdrawn at Mark 24.
Replaced by g05ssf.

#### Old Code

`CALL g05lmf(A,B,N,X,IGEN,ISEED,IFAIL)`

#### New Code

`CALL g05ssf(N,A,B,STATE,X,IFAIL)`
The integer array state in the call to g05ssf contains information on the base generator being used. This array must have been initialized prior to calling g05ssf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization.

### g05lnf

Withdrawn at Mark 24.
Replaced by g05slf.

#### Old Code

`CALL g05lnf(A,B,N,X,IGEN,ISEED,IFAIL)`

#### New Code

`CALL g05slf(N,A,B,STATE,X,IFAIL)`
The integer array state in the call to g05slf contains information on the base generator being used. This array must have been initialized prior to calling g05slf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization.

### g05lpf

Withdrawn at Mark 24.
Replaced by g05srf.

#### Old Code

`CALL g05lpf(VK,N,X,IGEN,ISEED,IFAIL)`

#### New Code

`CALL g05srf(N,VK,STATE,X,IFAIL)`
The integer array state in the call to g05srf contains information on the base generator being used. This array must have been initialized prior to calling g05srf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization.

### g05lqf

Withdrawn at Mark 24.
Replaced by g05sgf.

#### Old Code

`CALL g05lqf(NMIX,A,WGT,N,X,IGEN,ISEED,IFAIL)`

#### New Code

`CALL g05sgf(N,NMIX,A,WGT,STATE,X,IFAIL)`
The integer array state in the call to g05sgf contains information on the base generator being used. This array must have been initialized prior to calling g05sgf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization.

### g05lxf

Withdrawn at Mark 24.
Replaced by g05ryf.

#### Old Code

`CALL g05lxf(MODE,DF,M,XMU,C,LDC,N,X,LDX,IGEN,ISEED,R,LR,IFAIL)`

#### New Code

```IF (MODE == 0) THEN
NMODE = 1
ELSE IF (MODE == 1) THEN
NMODE = 0
ELSE
NMODE = MODE
END IF
CALL g05ryf(NMODE,N,DF,M,XMU,C,LDC,R,LR,STATE,X,LDX,IFAIL)```
The integer array state in the call to g05ryf contains information on the base generator being used. This array must have been initialized prior to calling g05ryf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization.

### g05lyf

Withdrawn at Mark 24.
Replaced by g05rzf.

#### Old Code

`g05lyf(MODE,M,XMU,C,LDC,N,X,LDX,IGEN,ISEED,R,LR,IFAIL)`

#### New Code

```IF (MODE == 0) THEN
NMODE = 1
ELSE IF (MODE == 1) THEN
NMODE = 0
ELSE
NMODE = MODE
END IF
CALL g05rzf(NMODE,N,M,XMU,C,LDC,R,LR,STATE,X,LDX,IFAIL)```
The integer array state in the call to g05rzf contains information on the base generator being used. This array must have been initialized prior to calling g05rzf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization.

### g05lzf

Withdrawn at Mark 24.
Replaced by g05rzf.

#### Old Code

`CALL g05lzf(MODE,M,XMU,C,LDC,X,IGEN,ISEED,R,LR,IFAIL)`

#### New Code

```IF (MODE == 0) THEN
NMODE = 1
ELSE IF (MODE == 1) THEN
NMODE = 0
ELSE
NMODE = MODE
END IF
CALL g05rzf(NMODE,N,M,XMU,C,LDC,R,LR,STATE,X,LDX,IFAIL)```
The integer array state in the call to g05rzf contains information on the base generator being used. This array must have been initialized prior to calling g05rzf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization.

### g05maf

Withdrawn at Mark 24.
Replaced by g05tlf.

#### Old Code

`CALL g05maf(A,B,N,X,IGEN,ISEED,IFAIL)`

#### New Code

`CALL g05tlf(N,A,B,STATE,X,IFAIL)`
The integer array state in the call to g05tlf contains information on the base generator being used. This array must have been initialized prior to calling g05tlf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization.

### g05mbf

Withdrawn at Mark 24.
Replaced by g05tcf.

#### Old Code

`CALL g05mbf(MODE,P,N,X,IGEN,ISEED,R,NR,IFAIL)`

#### New Code

```CALL g05tcf(MODE,N,P,R,LR,STATE,X,IFAIL)
DO 20 I = 1, N
X(I) = X(I) + 1
20 CONTINUE```
g05mbf returned the number of trials required to get the first success, whereas g05tcf returns the number of failures before the first success, therefore the value returned by g05tcf is one less than the equivalent value returned from g05mbf.
The integer array state in the call to g05tcf contains information on the base generator being used. This array must have been initialized prior to calling g05tcf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization.

### g05mcf

Withdrawn at Mark 24.
Replaced by g05thf.

#### Old Code

`CALL g05mcf(MODE,M,P,N,X,IGEN,ISEED,R,NR,IFAIL)`

#### New Code

`CALL g05thf(MODE,N,M,P,R,LR,STATE,X,IFAIL)`
The integer array state in the call to g05thf contains information on the base generator being used. This array must have been initialized prior to calling g05thf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization.

### g05mdf

Withdrawn at Mark 24.
Replaced by g05tff.

#### Old Code

`CALL g05mdf(MODE,A,N,X,IGEN,ISEED,R,NR,IFAIL)`

#### New Code

`CALL g05tff(MODE,N,A,R,LR,STATE,X,IFAIL)`
The integer array state in the call to g05tff contains information on the base generator being used. This array must have been initialized prior to calling g05tff with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization.

### g05mef

Withdrawn at Mark 24.
Replaced by g05tkf.

#### Old Code

`CALL g05mef(M,VLAMDA,X,IGEN,ISEED,IFAIL)`

#### New Code

`CALL g05tkf(M,VLAMDA,STATE,X,IFAIL)`
The integer array state in the call to g05tkf contains information on the base generator being used. This array must have been initialized prior to calling g05tkf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization.

### g05mjf

Withdrawn at Mark 24.
Replaced by g05taf.

#### Old Code

`CALL g05mjf(MODE,M,P,N,X,IGEN,ISEED,R,NR,IFAIL)`

#### New Code

`CALL g05taf(MODE,N,M,P,R,LR,STATE,X,IFAIL)`
The integer array state in the call to g05taf contains information on the base generator being used. This array must have been initialized prior to calling g05taf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization.

### g05mkf

Withdrawn at Mark 24.
Replaced by g05tjf.

#### Old Code

`CALL g05mkf(MODE,LAMBDA,N,X,IGEN,ISEED,R,NR,IFAIL)`

#### New Code

`CALL g05tjf(MODE,N,LAMBDA,R,LR,STATE,X,IFAIL)`
The integer array state in the call to g05tjf contains information on the base generator being used. This array must have been initialized prior to calling g05tjf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization.

### g05mlf

Withdrawn at Mark 24.
Replaced by g05tef.

#### Old Code

`CALL g05mlf(MODE,NS,NP,M,N,X,IGEN,ISEED,R,NR,IFAIL)`

#### New Code

`CALL g05tef(MODE,N,NS,NP,M,R,LR,STATE,X,IFAIL)`
The integer array state in the call to g05tef contains information on the base generator being used. This array must have been initialized prior to calling g05tef with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization.

### g05mrf

Withdrawn at Mark 24.
Replaced by g05tgf.

#### Old Code

`CALL g05mrf(MODE,M,K,P,N,X,LDX,IGEN,ISEED,R,NR,IFAIL)`

#### New Code

`CALL g05tgf(MODE,N,M,K,P,R,LR,STATE,X,LDX,IFAIL)`
The integer array state in the call to g05tgf contains information on the base generator being used. This array must have been initialized prior to calling g05tgf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization.

### g05mzf

Withdrawn at Mark 24.
Replaced by g05tdf.

#### Old Code

`CALL g05mzf(MODE,P,NP,IP1,ITYPE,N,X,IGEN,ISEED,R,NR,IFAIL)`

#### New Code

`CALL g05tdf(MODE,N,P,NP,IP1,ITYPE,R,LR,STATE,X,IFAIL)`
The integer array state in the call to g05tdf contains information on the base generator being used. This array must have been initialized prior to calling g05tdf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization.

### g05naf

Withdrawn at Mark 24.
Replaced by g05ncf.

#### Old Code

`CALL g05naf(INDEX,N,IGEN,ISEED,IFAIL)`

#### New Code

`CALL g05ncf(INDEX,N,STATE,IFAIL)`
The integer array state in the call to g05ncf contains information on the base generator being used. This array must have been initialized prior to calling g05ncf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization.

### g05nbf

Withdrawn at Mark 24.
Replaced by g05ndf.

#### Old Code

`CALL g05nbf(IPOP,N,ISAMPL,M,IGEN,ISEED,IFAIL)`

#### New Code

`CALL g05ndf(IPOP,N,ISAMPL,M,STATE,IFAIL)`
The integer array state in the call to g05ndf contains information on the base generator being used. This array must have been initialized prior to calling g05ndf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization.

### g05paf

Withdrawn at Mark 24.
Replaced by g05phf.

#### Old Code

```CALL g05paf(MODE,XMEAN,IP,PHI,IQ,THETA,AVAR,VAR,N,X,IGEN,ISEED,R, &
NR,IFAIL)```

#### New Code

```CALL g05phf(MODE,N,XMEAN,IP,PHI,IQ,THETA,AVAR,R,LR,STATE,VAR,X,   &
IFAIL)```
The integer array state in the call to g05phf contains information on the base generator being used. This array must have been initialized prior to calling g05phf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization.

### g05pcf

Withdrawn at Mark 24.
Replaced by g05pjf.

#### Old Code

```CALL g05pcf(MODE,K,XMEAN,IP,PHI,IQ,THETA,VAR,LDV,N,X,IGEN,ISEED,R,   &
NR,IWORK,LIWORK,IFAIL)```

#### New Code

```CALL g05pjf(MODE,N,K,XMEAN,IP,PHI,IQ,THETA,VAR,LDV,R,LR,STATE,X,LDX, &
IFAIL)```
The integer array state in the call to g05pjf contains information on the base generator being used. This array must have been initialized prior to calling g05pjf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization.

### g05qaf

Withdrawn at Mark 24.
Replaced by g05pxf.

#### Old Code

`CALL g05qaf(SIDE,INIT,M,N,A,LDA,IGEN,ISEED,WK,IFAIL)`

#### New Code

`CALL g05pxf(SIDE,INIT,M,N,STATE,A,LDA,IFAIL)`
The integer array state in the call to g05pxf contains information on the base generator being used. This array must have been initialized prior to calling g05pxf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization.

### g05qbf

Withdrawn at Mark 24.
Replaced by g05pyf.

#### Old Code

`CALL g05qbf(N,D,C,LDC,EPS,IGEN,ISEED,WK,IFAIL)`

#### New Code

`CALL g05pyf(N,D,EPS,STATE,C,LDC,IFAIL)`
The integer array state in the call to g05pyf contains information on the base generator being used. This array must have been initialized prior to calling g05pyf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization.

### g05qdf

Withdrawn at Mark 24.
Replaced by g05pzf.

#### Old Code

```CALL g05qdf(MODE,NROW,NCOL,TOTR,TOTC,X,LDX,IGEN,ISEED,R,NR,IW,LIW, &
IFAIL)```

#### New Code

`CALL g05pzf(MODE,NROW,NCOL,TOTR,TOTC,R,LR,STATE,X,LDX,IFAIL)`
The integer array state in the call to g05pzf contains information on the base generator being used. This array must have been initialized prior to calling g05pzf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization.

### g05raf

Withdrawn at Mark 24.
Replaced by g05rdf.

#### Old Code

`CALL g05raf(MODE,M,C,LDC,N,X,LDX,IGEN,ISEED,R,LR,IFAIL)`

#### New Code

```IF (MODE == 0) THEN
NMODE = 1
ELSE IF (MODE == 1) THEN
NMODE = 0
ELSE
NMODE = MODE
END IF
CALL CALL g05rdf(NMODE,N,M,C,LDC,R,LR,STATE,X,LDX,IFAIL)```
The integer array state in the call to g05rdf contains information on the base generator being used. This array must have been initialized prior to calling g05rdf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization.

### g05rbf

Withdrawn at Mark 24.
Replaced by g05rcf.

#### Old Code

`CALL g05rbf(MODE,DF,M,C,LDC,N,X,LDX,IGEN,ISEED,R,LR,IFAIL)`

#### New Code

```IF (MODE == 0) THEN
NMODE = 1
ELSE IF (MODE == 1) THEN
NMODE = 0
ELSE
NMODE = MODE
END IF
CALL CALL g05rcf(NMODE,N,DF,M,C,LDC,R,LR,STATE,X,LDX,IFAIL)```
The integer array state in the call to g05rcf contains information on the base generator being used. This array must have been initialized prior to calling g05rcf with a call to either g05kff or g05kgf. The required length of the array state will depend on the base generator chosen during initialization.

### g05yaf

Withdrawn at Mark 23.
Replaced by g05ylf and g05ymf.
Faure quasi-random numbers

#### Old Code

`CALL g05yaf(.TRUE.,'F',ISKIP,IDIM,QUAS,IREF,IFAIL)`

#### New Code

`CALL g05ylf(4,IDIM,IREF,LIREF,ISKIP,IFAIL)`

#### Old Code

`CALL g05yaf(.FALSE.,'F',ISKIP,IDIM,QUAS,IREF,IFAIL)`

#### New Code

`CALL g05ymf(1,2,QUAS,LDQUAS,IREF,IFAIL)`
Sobol quasi-random numbers

#### Old Code

`CALL g05yaf(.TRUE.,'S',ISKIP,IDIM,QUAS,IREF,IFAIL)`

#### New Code

`CALL g05ylf(2,IDIM,IREF,LIREF,ISKIP,IFAIL)`

#### Old Code

`CALL g05yaf(.FALSE.,'S',ISKIP,IDIM,QUAS,IREF,IFAIL)`

#### New Code

`CALL g05ymf(1,2,QUAS,LDQUAS,IREF,IFAIL)`
Neiderreiter quasi-random numbers

#### Old Code

`CALL g05yaf(.TRUE.,'N',ISKIP,IDIM,QUAS,IREF,IFAIL)`

#### New Code

`CALL g05ylf(3,IDIM,IREF,LIREF,ISKIP,IFAIL)`

#### Old Code

`CALL g05yaf(.FALSE.,'N',ISKIP,IDIM,QUAS,IREF,IFAIL)`

#### New Code

`CALL g05ymf(1,2,QUAS,LDQUAS,IREF,IFAIL)`

### g05ybf

Withdrawn at Mark 23.
Replaced by g05ylf and either g05yjf or g05ykf.
This routine has been replaced by a suite of routines consisting of the relevant initialization routine followed by one of two possible generator routines.
Faure quasi-random numbers with Gaussian probability:

#### Old Code

`CALL g05ybf(.TRUE.,'F',.FALSE.,MEAN,STD,ISKIP,IDIM,QUASI,IREF,IFAIL)`

#### New Code

`CALL g05ylf(4,IDIM,IREF,LIREF,ISKIP,IFAIL)`

#### Old Code

`CALL g05ybf(.FALSE.,'F',.FALSE.,MEAN,STD,ISKIP,IDIM,QUASI,IREF,IFAIL)`

#### New Code

`CALL g05yjf(MEAN,STD,N,QUASI,IREF,IFAIL)`
Sobol quasi-random numbers with Gaussian probability:

#### Old Code

`CALL g05ybf(.TRUE.,'S',.FALSE.,MEAN,STD,ISKIP,IDIM,QUASI,IREF,IFAIL)`

#### New Code

`CALL g05ylf(2,IDIM,IREF,LIREF,ISKIP,IFAIL)`

#### Old Code

`CALL g05ybf(.FALSE.,'S',.FALSE.,MEAN,STD,ISKIP,IDIM,QUASI,IREF,IFAIL)`

#### New Code

`CALL g05yjf(MEAN,STD,N,QUASI,IREF,IFAIL)`
Neiderreiter quasi-random numbers with Gaussian probability:

#### Old Code

`CALL g05ybf(.TRUE.,'N',.FALSE.,MEAN,STD,ISKIP,IDIM,QUASI,IREF,IFAIL)`

#### New Code

`CALL g05ylf(3,IDIM,IREF,LIREF,ISKIP,IFAIL)`

#### Old Code

`CALL g05ybf(.FALSE.,'N',.FALSE.,MEAN,STD,ISKIP,IDIM,QUASI,IREF,IFAIL)`

#### New Code

`CALL g05yjf(MEAN,STD,N,QUASI,IREF,IFAIL)`
Faure quasi-random numbers with log Normal probability:

#### Old Code

`CALL g05ybf(.TRUE.,'F',.TRUE.,MEAN,STD,ISKIP,IDIM,QUASI,IREF,IFAIL)`

#### New Code

`CALL g05ylf(4,IDIM,IREF,LIREF,ISKIP,IFAIL)`

#### Old Code

`CALL g05ybf(.FALSE.,'F',.TRUE.,MEAN,STD,ISKIP,IDIM,QUASI,IREF,IFAIL)`

#### New Code

`CALL g05ykf(MEAN,STD,N,QUASI,IREF,IFAIL)`
Sobol quasi-random numbers with log Normal probability:

#### Old Code

`CALL g05ybf(.TRUE.,'S',.TRUE.,MEAN,STD,ISKIP,IDIM,QUASI,IREF,IFAIL)`

#### New Code

`CALL g05ylf(2,IDIM,IREF,LIREF,ISKIP,IFAIL)`

#### Old Code

`CALL g05ybf(.FALSE.,'S',.TRUE.,MEAN,STD,ISKIP,IDIM,QUASI,IREF,IFAIL)`

#### New Code

`CALL g05ykf(MEAN,STD,N,QUASI,IREF,IFAIL)`
Neiderreiter quasi-random numbers with log Normal probability:

#### Old Code

`CALL g05ybf(.TRUE.,'N',.TRUE.,MEAN,STD,ISKIP,IDIM,QUASI,IREF,IFAIL)`

#### New Code

`CALL g05ylf(3,IDIM,IREF,LIREF,ISKIP,IFAIL)`

#### Old Code

`CALL g05ybf(.FALSE.,'N',.TRUE.,MEAN,STD,ISKIP,IDIM,QUASI,IREF,IFAIL)`

#### New Code

`CALL g05ykf(MEAN,STD,N,QUASI,IREF,IFAIL)`

### g05ycf

Withdrawn at Mark 24.
Replaced by g05ylf.

#### Old Code

`CALL g05ycf(IDIM,IREF,IFAIL)`

#### New Code

```GENID = 4
CALL g05ylf(GENID,IDIM,IREF,LIREF,ISKIP,IFAIL)```

### g05ydf

Withdrawn at Mark 24.
Replaced by g05ymf.

#### Old Code

`CALL g05ydf(N,QUASI,IREF,IFAIL)`

#### New Code

`CALL g05ymf(N,QUAS,LDQUAS,IREF,IFAIL)`

### g05yef

Withdrawn at Mark 24.
Replaced by g05ylf.

#### Old Code

`CALL g05yef(IDIM, IREF, ISKIP, IFAIL)`

#### New Code

```GENID = 2
CALL g05ylf(GENID,IDIM,IREF,LIREF,ISKIP,IFAIL)```

### g05yff

Withdrawn at Mark 24.
Replaced by g05ymf.

#### Old Code

`CALL g05yff(N,QUASI,IREF,IFAIL)`

#### New Code

`CALL g05ymf(N,QUAS,LDQUAS,IREF,IFAIL)`

### g05ygf

Withdrawn at Mark 24.
Replaced by g05ylf.

#### Old Code

`CALL g05ygf(IDIM,IREF,ISKIP,IFAIL)`

#### New Code

```GENID = 3
CALL g05ylf(GENID,IDIM,IREF,LIREF,ISKIP,IFAIL)```

### g05yhf

Withdrawn at Mark 24.
Replaced by g05ymf.

#### Old Code

`CALL g05yhf(N,QUASI,IREF,IFAIL)`

#### New Code

`CALL g05ymf(N,RCORD,QUAS,LDQUAS,IREF,IFAIL)`

### g05zaf

Withdrawn at Mark 22.
There is no replacement for this routine.
g05zaf was used to select the underlying generator for the old style random number generation routines. These routines are no longer available and hence no direct replacement routine for g05zaf is required. See g05kff and g05kgf for details on how to select the underlying generator for the newer style routines.

## G10 – Smoothing in Statistics

### g10baf

Withdrawn at Mark 27.
Replaced by g10bbf.
Withdrawn primarily due to threadsafety. The replacement routine also introduces new functionality with respect to the automatic selection of a suitable window width.

#### Old Code

`CALL g10baf(N,X,WINDOW,SLO,SHI,NS,SMOOTH,T,USEFFT,FFT,IFAIL)`

#### New Code

```ALLOCATE(RCOMM,NS+20)
CALL g10bbf(N,X,1,WINDOW,SLO,SHI,NS,SMOOTH,T,USEFFT,RCOMM,IFAIL)
! the next step is only required if the information in FFT
! was being used outside another call to G10BAF
FFT(1:NS) = RCOMM(21:NS+20)```

## G13 – Time Series Analysis

### g13dcf

Withdrawn at Mark 24.
Replaced by g13ddf.

#### Old Code

```CALL g13dcf(K,N,IP,IQ,MEAN,PAR,NPAR,QQ,KMAX,W,PARHLD,EXACT,IPRINT,  &
CGETOL,MAXCAL,ISHOW,NITER,RLOGL,V,G,CM,LDCM,WORK,LWORK, &
IW,LIW,IFAIL)```

#### New Code

```CALL g13ddf(K,N,IP,IQ,MEAN,PAR,NPAR,QQ,KMAX,W,PARHLD,EXACT,IPRINT,  &
CGETOL,MAXCAL,ISHOW,NITER,RLOGL,V,G,CM,LDCM,IFAIL)```
The workspace arguments work, lwork, iw and liw are no longer required in the call to g13ddf.

## P01 – Error Trapping

### p01abf

Withdrawn at Mark 24.
There is no replacement for this routine.

## X02 – Machine Constants

### x02daf

Withdrawn at Mark 24.
There is no replacement for this routine.

### x02djf

Withdrawn at Mark 24.
There is no replacement for this routine.