NAG Library Manual, Mark 30.1
Interfaces:  FL   CL   CPP   AD 

NAG FL Interface Introduction
Example description
    Program g02mafe

!     G02MAF Example Program Text

!     Mark 30.1 Release. NAG Copyright 2024.

!     .. Use Statements ..
      Use nag_library, Only: g02maf, nag_wp
!     .. Implicit None Statement ..
      Implicit None
!     .. Parameters ..
      Integer, Parameter               :: nin = 5, nout = 6
!     .. Local Scalars ..
      Integer                          :: i, ifail, ip, k, ldb, ldd, lisx,     &
                                          lropt, m, mnstep, mtype, n, nstep,   &
                                          pred, prey
!     .. Local Arrays ..
      Real (Kind=nag_wp), Allocatable  :: b(:,:), d(:,:), fitsum(:,:),         &
                                          ropt(:), y(:)
      Integer, Allocatable             :: isx(:)
!     .. Intrinsic Procedures ..
      Intrinsic                        :: count, floor, max, repeat
!     .. Executable Statements ..
      Write (nout,*) 'G02MAF Example Program Results'
      Write (nout,*)

!     Skip heading in data file
      Read (nin,*)

!     Read in the problem size
      Read (nin,*) n, m

!     Read in the model specification
      Read (nin,*) mtype, pred, prey, mnstep, lisx

!     Read in variable inclusion flags (if specified) and calculate IP
      Allocate (isx(lisx))
      If (lisx==m) Then
        Read (nin,*) isx(1:lisx)
        ip = count(isx(1:m)==1)
      Else
        ip = m
      End If

!     Optional arguments (using defaults)
      lropt = 0
      Allocate (ropt(lropt))

!     Read in the data
      ldd = n
      Allocate (y(n),d(ldd,m))
      Read (nin,*)(d(i,1:m),y(i),i=1,n)

!     Allocate output arrays
      ldb = ip
      Allocate (b(ldb,mnstep+2),fitsum(6,mnstep+1))

!     Call the model fitting routine
      ifail = -1
      Call g02maf(mtype,pred,prey,n,m,d,ldd,isx,lisx,y,mnstep,ip,nstep,b,ldb,  &
        fitsum,ropt,lropt,ifail)
      If (ifail/=0) Then
        If (ifail/=112 .And. ifail/=161 .And. ifail/=162 .And. ifail/=163)     &
          Then
!         IFAIL = 112, 161, 162 and 163 are warnings, so no need to terminate
!         if they occur
          Go To 100
        End If
      End If

!     Display the parameter estimates
      Write (nout,*) ' Step ', repeat(' ',max((ip-2),0)*5),                    &
        ' Parameter Estimate'
      Write (nout,*) repeat('-',5+ip*10)
      Do k = 1, nstep
        Write (nout,99998) k, b(1:ip,k)
      End Do
      Write (nout,*)
      Write (nout,99999) 'alpha: ', fitsum(1,nstep+1)
      Write (nout,*)
      Write (nout,*)                                                           &
        ' Step     Sum      RSS       df       Cp       Ck     Step Size'
      Write (nout,*) repeat('-',64)
      Do k = 1, nstep
        Write (nout,99997) k, fitsum(1:2,k), floor(fitsum(3,k)+0.5_nag_wp),    &
          fitsum(4:6,k)
      End Do
      Write (nout,*)
      Write (nout,99999) 'sigma^2: ', fitsum(5,nstep+1)

100   Continue
99999 Format (1X,A,F9.3)
99998 Format (2X,I3,10(1X,F9.3))
99997 Format (2X,I3,2(1X,F9.3),1X,I6,1X,3(1X,F9.3))
    End Program g02mafe