```    Program f08kgfe

!     F08KGF Example Program Text

!     Mark 26.2 Release. NAG Copyright 2017.

!     .. Use Statements ..
Use nag_library, Only: dgebrd, dgelqf, dgeqrf, dorglq, dorgqr, dormbr,   &
f06qff, f06qhf, nag_wp, x04caf
!     .. Implicit None Statement ..
Implicit None
!     .. Parameters ..
Real (Kind=nag_wp), Parameter    :: zero = 0.0E0_nag_wp
Integer, Parameter               :: nin = 5, nout = 6
!     .. Local Scalars ..
Integer                          :: i, ic, ifail, info, lda, ldpt, ldu,  &
lwork, m, n
!     .. Local Arrays ..
Real (Kind=nag_wp), Allocatable  :: a(:,:), d(:), e(:), pt(:,:), tau(:), &
taup(:), tauq(:), u(:,:), work(:)
!     .. Executable Statements ..
Write (nout,*) 'F08KGF Example Program Results'
!     Skip heading in data file
Do ic = 1, 2
lda = m
ldpt = n
ldu = m
lwork = 64*(m+n)
Allocate (a(lda,n),d(n),e(n-1),pt(ldpt,n),tau(n),taup(n),tauq(n),      &
u(ldu,n),work(lwork))

!       Read A from data file

If (m>=n) Then

!         Compute the QR factorization of A
!         The NAG name equivalent of dgeqrf is f08aef
Call dgeqrf(m,n,a,lda,tau,work,lwork,info)

!         Copy A to U
Call f06qff('Lower',m,n,a,lda,u,ldu)

!         Form Q explicitly, storing the result in U
!         The NAG name equivalent of dorgqr is f08aff
Call dorgqr(m,n,n,u,ldu,tau,work,lwork,info)

!         Copy R to PT (used as workspace)
Call f06qff('Upper',n,n,a,lda,pt,ldpt)

!         Set the strictly lower triangular part of R to zero
Call f06qhf('Lower',n-1,n-1,zero,zero,pt(2,1),ldpt)

!         Bidiagonalize R
!         The NAG name equivalent of dgebrd is f08kef
Call dgebrd(n,n,pt,ldpt,d,e,tauq,taup,work,lwork,info)

!         Update Q, storing the result in U
!         The NAG name equivalent of dormbr is f08kgf
Call dormbr('Q','Right','No transpose',m,n,n,pt,ldpt,tauq,u,ldu,     &
work,lwork,info)

!         Print bidiagonal form and matrix Q

Write (nout,*)
Write (nout,*) 'Example 1: bidiagonal matrix B'
Write (nout,*) 'Diagonal'
Write (nout,99999) d(1:n)
Write (nout,*) 'Superdiagonal'
Write (nout,99999) e(1:n-1)
Write (nout,*)
Flush (nout)

!         ifail: behaviour on error exit
!                =0 for hard exit, =1 for quiet-soft, =-1 for noisy-soft
ifail = 0
Call x04caf('General',' ',m,n,u,ldu,'Example 1: matrix Q',ifail)

Else

!         Compute the LQ factorization of A
!         The NAG name equivalent of dgelqf is f08ahf
Call dgelqf(m,n,a,lda,tau,work,lwork,info)

!         Copy A to PT
Call f06qff('Upper',m,n,a,lda,pt,ldpt)

!         Form Q explicitly, storing the result in PT
!         The NAG name equivalent of dorglq is f08ajf
Call dorglq(n,n,m,pt,ldpt,tau,work,lwork,info)

!         Copy L to U (used as workspace)
Call f06qff('Lower',m,m,a,lda,u,ldu)

!         Set the strictly upper triangular part of L to zero
Call f06qhf('Upper',m-1,m-1,zero,zero,u(1,2),ldu)

!         Bidiagonalize L
!         The NAG name equivalent of dgebrd is f08kef
Call dgebrd(m,m,u,ldu,d,e,tauq,taup,work,lwork,info)

!         Update P**T, storing the result in PT
!         The NAG name equivalent of dormbr is f08kgf
Call dormbr('P','Left','Transpose',m,n,m,u,ldu,taup,pt,ldpt,work,    &
lwork,info)

!         Print bidiagonal form and matrix P**T

Write (nout,*)
Write (nout,*) 'Example 2: bidiagonal matrix B'
Write (nout,*) 'Diagonal'
Write (nout,99999) d(1:m)
Write (nout,*) 'Superdiagonal'
Write (nout,99999) e(1:m-1)
Write (nout,*)
Flush (nout)

ifail = 0
Call x04caf('General',' ',m,n,pt,ldpt,'Example 2: matrix P**T',      &
ifail)

End If
Deallocate (a,d,e,pt,tau,taup,tauq,u,work)
End Do

99999 Format (3X,(8F8.4))
End Program f08kgfe
```