Program e04rdfe

!     E04RDF Example Program Text
!     Mark 28.3 Release. NAG Copyright 2022.

!     Load a linear semidefinite programming problem from a sparse SDPA
!     file, formulate the problem via a handle, pass it to the solver
!     and print both primal and dual variables.

!     .. Use Statements ..
      Use, Intrinsic                   :: iso_c_binding, Only: c_ptr
      Use nag_library, Only: e04raf, e04rdf, e04ref, e04rnf, e04ryf, e04rzf,   &
                             e04svf, e04zmf, nag_wp, x04acf, x04adf, x04ccf
!     .. Implicit None Statement ..
      Implicit None
!     .. Parameters ..
      Integer, Parameter               :: filelst = 0, infile = 42, nout = 6
      Character (*), Parameter         :: fname_default = 'e04rdfe.opt'
!     .. Local Scalars ..
      Type (c_ptr)                     :: handle
      Integer                          :: idblk, idx, ifail, ifail_e04rd,      &
                                          inform, k, maxnblk, maxnnz, maxnvar, &
                                          nblk, nnz, nnzu, nnzua, nnzuc,       &
                                          ntests, nvar
      Character (256)                  :: fname
      Character (60)                   :: title
!     .. Local Arrays ..
      Real (Kind=nag_wp), Allocatable  :: a(:), cvec(:), u(:), ua(:), uc(:),   &
                                          x(:)
      Real (Kind=nag_wp)               :: rinfo(32), stats(32)
      Integer, Allocatable             :: blksizea(:), icola(:), irowa(:),     &
                                          nnza(:)
!     .. Intrinsic Procedures ..
      Intrinsic                        :: command_argument_count,              &
                                          get_command_argument, sum, trim
!     .. Executable Statements ..

      Write (nout,*) 'E04RDF Example Program Results'
      Write (nout,*)

!     Use the first command line argument as the filename or
!     choose default hard-coded filename in 'fname_default'.
      ntests = command_argument_count()
      If (ntests==0) Then
!       Assume the default filename.
        fname = fname_default
      Else
        Call get_command_argument(1,fname)
      End If

      Write (nout,*) 'Reading SDPA file: ', trim(fname)
      Flush (nout)

!     Open the input file.
      ifail = 0
      Call x04acf(infile,fname,0,ifail)

!     Go through the file and find the dimension of the problem.
!     Unless the file format is wrong, the routine should finish
!     with IFAIL = 1 (not enough space).
      maxnvar = 0
      maxnblk = 0
      maxnnz = 0
      Allocate (cvec(maxnvar),nnza(maxnvar+1),irowa(maxnnz),icola(maxnnz),     &
        a(maxnnz),blksizea(maxnblk))
      ifail_e04rd = -1
      Call e04rdf(infile,maxnvar,maxnblk,maxnnz,filelst,nvar,nblk,nnz,cvec,    &
        nnza,irowa,icola,a,blksizea,ifail_e04rd)
      Deallocate (cvec,nnza,irowa,icola,a,blksizea)

!     Close the file, it will need to be reopened later.
      ifail = 0
      Call x04adf(infile,ifail)

      If (ifail_e04rd/=1) Then
!       Possible problem with formatting, etc.
        Write (nout,99999) 'Reading the SDPA file failed with IFAIL = ', ifail
99999   Format (1X,A,I3)
        Write (nout,*) 'Terminating the example program.'
        Go To 100
      End If

!     Allocate the right size of arrays for the data.
      Write (nout,*) 'Allocating space for the problem.'
      Write (nout,Fmt=99998) 'NVAR = ', nvar
      Write (nout,Fmt=99998) 'NBLK = ', nblk
      Write (nout,Fmt=99998) 'NNZ  = ', nnz
99998 Format (6X,A,I7)
      Flush (nout)

      maxnvar = nvar
      maxnblk = nblk
      maxnnz = nnz

      Allocate (cvec(maxnvar),nnza(maxnvar+1),irowa(maxnnz),icola(maxnnz),     &
        a(maxnnz),blksizea(maxnblk))

!     Reopen the file.
      ifail = 0
      Call x04acf(infile,fname,0,ifail)

!     Read the problem data, there should be enough space this time.
      ifail = 0
      Call e04rdf(infile,maxnvar,maxnblk,maxnnz,filelst,nvar,nblk,nnz,cvec,    &
        nnza,irowa,icola,a,blksizea,ifail)

!     Close the file.
      ifail = 0
      Call x04adf(infile,ifail)

!     Problem was successfully decoded.
      Write (nout,*)                                                           &
        'Linear SDP problem was read, start formulating the problem'
      Flush (nout)

!     Initialize the handle of the problem.
      ifail = 0
      Call e04raf(handle,nvar,ifail)

!     Add the linear objective function to the formulation.
      ifail = 0
      Call e04ref(handle,nvar,cvec,ifail)

!     Add all linear matrix constraints to the formulation.
      idblk = 0
      ifail = 0
      Call e04rnf(handle,nvar,sum(blksizea(1:nblk)),nnza,nnz,irowa,icola,a,    &
        nblk,blksizea,idblk,ifail)

      Write (nout,*) 'The problem formulation in a handle is completed.'
      Write (nout,*)
      Flush (nout)

!     Print overview of the handle.
      ifail = 0
      Call e04ryf(handle,nout,'Overview',ifail)

!     Set optional arguments.
      ifail = 0
      Call e04zmf(handle,'DIMACS Measures = Check',ifail)
      ifail = 0
      Call e04zmf(handle,'Initial X = Automatic',ifail)

!     Compute memory needed for primal & dual variables.

!     There are no box constraints or linear constraints set
!     by E04RHF or E04RJF, neither second order cone constraints.
      nnzu = 0
      nnzuc = 0

!     Count size of the matrix multipliers, stored as packed
!     triangle respecting the block structure.
      nnzua = 0
      Do k = 1, nblk
        nnzua = nnzua + blksizea(k)*(blksizea(k)+1)/2
      End Do

      Allocate (x(nvar),ua(nnzua),u(nnzu),uc(nnzuc))

!     Call the solver.
      ifail = 0
      Call e04svf(handle,nvar,x,nnzu,u,nnzuc,uc,nnzua,ua,rinfo,stats,inform,   &
        ifail)

!     Print results.

      Write (nout,*)
      Write (nout,*) 'Optimal solution:'
      Write (nout,99997) x(1:nvar)
99997 Format (1X,'X = ',2F9.2)
      Flush (nout)

!     Print packed lower triangles of the Lagrangian multipliers.
      idx = 1
      Do k = 1, nblk
        Write (title,99996) 'Lagrangian multiplier for A_', k
99996   Format (A,I0)
        nnz = blksizea(k)*(blksizea(k)+1)/2
        ifail = 0
        Call x04ccf('Lower','N',blksizea(k),ua(idx:idx+nnz-1),title,ifail)
        idx = idx + nnz
      End Do

!     Deallocate memory within the handle.
      ifail = 0
      Call e04rzf(handle,ifail)

100   Continue
    End Program e04rdfe