NAG Library Manual, Mark 29.3
```!   G13NBF Example Program Text
!   Mark 29.3 Release. NAG Copyright 2023.

Module g13nbfe_mod

!     G13NBF Example Program Module:
!            Parameters and User-defined Routines

!     .. Use Statements ..
Use nag_library, Only: nag_wp
!     .. Implicit None Statement ..
Implicit None
!     .. Accessibility Statements ..
Private
Public                           :: costfn, get_data
Contains
Subroutine costfn(ts,nr,r,c,y,iuser,ruser,info)
!       Cost function, C. This cost function is based on the likelihood of
!       the gamma distribution

!       .. Scalar Arguments ..
Integer, Intent (Inout)        :: info
Integer, Intent (In)           :: nr, ts
!       .. Array Arguments ..
Real (Kind=nag_wp), Intent (Out) :: c(nr)
Real (Kind=nag_wp), Intent (Inout) :: ruser(*), y(0:*)
Integer, Intent (Inout)        :: iuser(*)
Integer, Intent (In)           :: r(nr)
!       .. Local Scalars ..
Real (Kind=nag_wp)             :: dn, shape, si
Integer                        :: i
!       .. Intrinsic Procedures ..
Intrinsic                      :: log, real
!       .. Executable Statements ..

!       RUSER(1) holds the shape parameter (a) for the gamma distribution
shape = ruser(1)

!       Test which way around TS and R are (only needs to be done once)
If (ts<r(1)) Then
Do i = 1, nr
si = y(r(i)) - y(ts)
dn = real(r(i)-ts,kind=nag_wp)
c(i) = 2.0_nag_wp*dn*shape*(log(si)-log(dn*shape))
End Do

Else
Do i = 1, nr
si = y(ts) - y(r(i))
dn = real(ts-r(i),kind=nag_wp)
c(i) = 2.0_nag_wp*dn*shape*(log(si)-log(dn*shape))
End Do
End If

!       Set info nonzero to terminate execution for any reason
info = 0
End Subroutine costfn

Subroutine get_data(nin,n,k,y,iuser,ruser)
!       Read in data that is specific to the cost function

!       .. Scalar Arguments ..
Real (Kind=nag_wp), Intent (Out) :: k
Integer, Intent (In)           :: n, nin
!       .. Array Arguments ..
Real (Kind=nag_wp), Allocatable, Intent (Out) :: ruser(:), y(:)
Integer, Allocatable, Intent (Out) :: iuser(:)
!       .. Local Scalars ..
Real (Kind=nag_wp)             :: shape
Integer                        :: i
!       .. Executable Statements ..

!       Read in the series of interest
!       NB: we are starting Y allocation at 0 as we manipulate
!       the data in Y in a moment
Allocate (y(0:n))

!       Read in the shape parameter for the Gamma distribution

!       Store the shape parameter in RUSER. IUSER is not used
Allocate (ruser(1),iuser(0))
ruser(1) = shape

!       The cost function is a function of the sum of Y, so for
!       efficiency we will calculate the cumulative sum
!       It should be noted that this may introduce some rounding issues
!       with very extreme data
y(0) = 0.0_nag_wp
Do i = 1, n
y(i) = y(i-1) + y(i)
End Do

!       The value of K is defined by the cost function being used
!       in this example a value of 0.0 is the required value
k = 0.0_nag_wp

Return
End Subroutine get_data
End Module g13nbfe_mod

Program g13nbfe

!     .. Use Statements ..
Use g13nbfe_mod, Only: costfn, get_data
Use nag_library, Only: g13nbf, nag_wp
!     .. Implicit None Statement ..
Implicit None
!     .. Parameters ..
Integer, Parameter               :: nin = 5, nout = 6
!     .. Local Scalars ..
Real (Kind=nag_wp)               :: beta, k
Integer                          :: i, ifail, minss, n, ntau
!     .. Local Arrays ..
Real (Kind=nag_wp), Allocatable  :: ruser(:), y(:)
Integer, Allocatable             :: iuser(:), tau(:)
!     .. Intrinsic Procedures ..
Intrinsic                        :: repeat
!     .. Executable Statements ..
Write (nout,*) 'G13NBF Example Program Results'
Write (nout,*)

!     Skip heading in data file

!     Read in the problem size, penalty and minimum segment size

!     Read in the rest of the data, that (may be) dependent on the cost
!     function
Call get_data(nin,n,k,y,iuser,ruser)

!     Allocate output arrays
Allocate (tau(n))

!     Call routine to detect change points
ifail = 0
Call g13nbf(n,beta,minss,k,costfn,ntau,tau,y,iuser,ruser,ifail)

!     Display the results
Write (nout,99999) ' -- Change Points --'
Write (nout,99999) ' Number     Position'
Write (nout,99999) repeat('=',21)
Do i = 1, ntau
Write (nout,99998) i, tau(i)
End Do

99999 Format (1X,A)
99998 Format (1X,I4,7X,I6)
End Program g13nbfe
```