! *****************************************************************
! *****************************************************************

program analyse

  implicit none

  ! SCALAR PARAMETERS
  integer, parameter :: nmeth = 2, nprob1 = 361, nprob2 = 157, nprob3 = 740, nprob = nprob1 + nprob2 + nprob3
  real(kind=8) :: epsfeas = 1.0d-08
  
  ! LOCAL SCALARS
  integer :: i,j,jmin,jmax,k,idum,pos
  real(kind=8) :: btol,ftol,tau,ddum
  
  ! LOCAL ARRAYS
  logical :: needsrerun(nprob),comparable(nprob),compandeqf(nprob)
  integer :: istop(nmeth,nprob),n(nmeth,nprob),m(nmeth,nprob),howmany(nmeth),vec(0:20)
  real(kind=8) :: f(nmeth,nprob),csupn(nmeth,nprob),bdsvio(nmeth,nprob),time(nmeth,nprob), &
       metric(nmeth,nprob),bestfun(nprob)
  character(len=23) :: mname(nmeth,nprob),pname(nmeth,nprob)

  character(len=500) :: line
  character(len=100) :: sdum
  
  i = 1

  write(*,*) 'reading table-unc-bound-constr.txt ...'
  
  open(unit=10,file='table-unc-bound-constr.txt')

  do j = 1,nprob1
     read(10,'(A)') line
     do while ( line(1:1) .eq. '#' )
        read(10,'(A)') line
     end do
     
     open(unit=20,file='tempo.txt',status='new')
     write(20,*) line
     
     rewind(20)
     
     read(20,*) mname(i,j),pname(i,j),istop(i,j),idum,time(i,j),n(i,j),f(i,j)
     
     csupn(i,j) = 0.0d0
     bdsvio(i,j) = 0.0d0
     
     close(20,status='delete')
  end do
  
  close(10)

  write(*,*) 'reading table-eqsys.txt ...'
  
  open(unit=10,file='table-eqsys.txt')

  do j = nprob1 + 1,nprob1 + nprob2
     read(10,'(A)') line
     do while ( line(1:1) .eq. '#' )
        read(10,'(A)') line
     end do
     
     open(unit=20,file='tempo.txt',status='new')
     write(20,*) line
     
     rewind(20)
     
     read(20,*) mname(i,j),pname(i,j),istop(i,j),idum,time(i,j),n(i,j),m(i,j),csupn(i,j)

     f(i,j) = 0.0d0
     bdsvio(i,j) = 0.0d0

     close(20,status='delete')
  end do
  
  close(10)

  write(*,*) 'reading table-nlp.txt ...'
  
  open(unit=10,file='table-nlp.txt')

  do j = nprob1 + nprob2 + 1,nprob1 + nprob2 + nprob3
     read(10,'(A)') line
     do while ( line(1:1) .eq. '#' )
        read(10,'(A)') line
     end do
     
     open(unit=20,file='tempo.txt',status='new')
     write(20,*) line
     
     rewind(20)
     
     read(20,*) mname(i,j),pname(i,j),istop(i,j),idum,time(i,j),n(i,j),m(i,j),f(i,j),csupn(i,j),ddum,ddum,bdsvio(i,j)

     close(20,status='delete')
  end do
  
  close(10)

  i = 2

  write(*,*) 'reading table-ipopt.txt ...'
  
  open(10,file='../../ipopt-experiments/table-ipopt.txt')

  do j = 1,nprob
     read(10,'(A)') line
     do while ( line(1:1) .eq. '#' )
        read(10,'(A)') line
     end do
     
     open(unit=20,file='tempo.txt',status='new')
     write(20,*) line
     rewind(20)

     needsrerun(j) = .false.
     
     if ( line(18:40) .eq. 'CPU time limit exceeded' ) then
        read(20,*) mname(i,j),pname(i,j)
        needsrerun(j) = .true.

     else
        read(20,*) mname(i,j),pname(i,j),n(i,j),m(i,j),istop(i,j),sdum,f(i,j),csupn(i,j),bdsvio(i,j),time(i,j)
     end if

     close(20,status='delete')
        
     if ( pname(i,j) .ne. pname(1,j) ) then
        write(*,*) 'Error in the input file.',pname(i,j),pname(1,j)
        stop
     end if
  end do

  do j = 1,nprob
     if ( needsrerun(j) ) then
        read(10,'(A)') line
        do while ( line(1:1) .eq. '#' )
           read(10,'(A)') line
        end do
     
        open(unit=20,file='tempo.txt',status='new')
        write(20,*) line
        rewind(20)

        read(20,*) mname(i,j),pname(i,j),n(i,j),m(i,j),istop(i,j),sdum,f(i,j),csupn(i,j),bdsvio(i,j),time(i,j)
        
        istop(i,j) = 99
        time(i,j) = 600.0d0
        
        close(20,status='delete')
        
        if ( pname(i,j) .ne. pname(1,j) ) then
           write(*,*) 'Error in the input file.',pname(i,j),pname(1,j)
           stop
        end if
     end if
  end do
  
  close(10)

  do k = 1,3
     vec(0:12) = 0
     if ( k .eq. 1 ) then
        jmin = 1
        jmax = nprob1
     else if ( k .eq. 2 ) then
        jmin = nprob1 + 1
        jmax = nprob1 + nprob2
     else
        jmin = nprob1 + nprob2 + 1
        jmax = nprob1 + nprob2 + nprob3
     end if
     do j = jmin,jmax
        if ( istop(1,j) .eq. 99 ) then
           vec(12) = vec(12) + 1
        else if ( istop(1,j) .eq. 20 ) then
           vec(11) = vec(11) + 1
        else
           if ( istop(1,j) .lt. 0 .or. istop(1,j) .gt. 10 ) then
              write(*,*) 'In problem ',j,' Algencan istop = ',istop(1,j)
           else
              vec(istop(1,j)) = vec(istop(1,j)) + 1
           end if
        end if
     end do

     write(*,*)
     write(*,*) 'Algencan stopping criteria summary:'
     write(*,*) 'vec = ',vec(0:12),' sum = ',sum(vec(0:12))
  end do
  
  vec(0:11) = 0
  do j = 1,nprob
     if ( istop(2,j) .eq. -13 ) then
        pos = 0
     else if ( istop(2,j) .eq. -10 ) then
        pos = 1
     else if ( istop(2,j) .eq.  -3 ) then
        pos = 2
     else if ( istop(2,j) .eq.  -2 ) then
        pos = 3
     else if ( istop(2,j) .eq.  -1 ) then
        pos = 4
     else if ( istop(2,j) .eq.   0 ) then
        pos = 5
     else if ( istop(2,j) .eq.   1 ) then
        pos = 6
     else if ( istop(2,j) .eq.   2 ) then
        pos = 7
     else if ( istop(2,j) .eq.   3 ) then
        pos = 8
     else if ( istop(2,j) .eq.   4 ) then
        pos = 9
     else if ( istop(2,j) .eq.   6 ) then
        pos = 10
     else if ( istop(2,j) .eq.  99 ) then
        pos = 11
     else
        write(*,*) 'In problem ',j,' Ipopt istop = ',istop(2,j)
     end if
     vec(pos) = vec(pos) + 1
  end do
  
  write(*,*)
  write(*,*) 'Ipopt stopping criteria summary:'
  write(*,*) 'vec = ',vec(0:11),' sum = ',sum(vec(0:11))

  write(*,*)
  write(*,*) 'Satisfaction of constraints with tolerance 1e-8 and bound constraints with tolerance 1e-1, 1e-2, ..., 1e-16, 0'
  do i = 1,nmeth
     do k = 1,17
        if ( k .eq. 1 ) then
           btol = 0.1d0
        else if ( 2 .le. k .and. k .le. 16 ) then
           btol = 0.1d0 * btol
        else ! if ( btol .eq. 17 ) then
           btol = 0.0d0
        end if
        vec(k) = count( csupn(i,1:nprob) .le. epsfeas .and. bdsvio(i,1:nprob) .le. btol )
     end do
     write(*,*) mname(i,1),vec(1:17)
  end do

  comparable(1:nprob) = csupn(1,1:nprob) .le. epsfeas .and. bdsvio(1,1:nprob) .le. 0.0d0 .and. &
                        csupn(2,1:nprob) .le. epsfeas .and. bdsvio(2,1:nprob) .le. 0.0d0

  write(*,*)
  write(*,*) 'We assume that performances are comparable if the methods found a point with equivalent feasibilities.'
  write(*,*) 'And we say feasibilities are equivalent if constraints (as originally formulated) are satisfied with '
  write(*,*) 'precision 1e-8 and bounds are satisfied with zero tolerance. There are ',count( comparable(1:nprob) ), &
       ' that satisfy these conditions.'
  
  do j = 1,nprob
     if ( comparable(j) ) then
        bestfun(j) = minval( f(1:nmeth,j) )
     end if
  end do

  write(*,*)
  write(*,*) 'Among comparable problems, number of best solutions found as a function o ftol: 1e-1, 1e-2, 1e-3, ..., 1e-8, 0'

  do i = 1,nmeth
     ftol = 1.0d0
     do k = 1,9
        if ( k .eq. 1 ) then
           ftol = 0.1d0
        else if ( 2 .le. k .and. k .le. 8 ) then
           ftol = 0.1d0 * ftol
        else ! if ( k .eq. 9 )
           ftol = 0.0d0
        end if
        vec(k) = count( comparable(1:nprob) .and. &
             f(i,1:nprob) .le. bestfun(1:nprob) + ftol * max( 1.0d0, abs( bestfun(1:nprob) ) ) )
     end do
     write(*,*) mname(i,1),vec(1:9)
  end do
  
  ftol = 1.0d-01

  write(*,*)
  write(*,*) 'Problems with comparable feasibility in which f differs (CONSIDERING FTOL=1e-1).'
  do j = 1,nprob
     if ( comparable(j) ) then
        if ( all( f(1:nmeth,j) .le. bestfun(j) + ftol * max( 1.0d0, abs( bestfun(j) ) ) ) ) then
           compandeqf(j) = .true.
        else
           compandeqf(j) = .false.
!!$           write(*,100) pname(1,j),f(1,j),csupn(1,j),bdsvio(1,j),f(2,j),csupn(2,j),bdsvio(2,j)
        end if
     else
        compandeqf(j) = .false.
     end if
  end do

  write(*,*)
  write(*,*) 'CONSIDERING FTOL=1e-1, all methods found comparable points with equivalent f in ', &
       count( compandeqf(1:nprob) ),' problems.'
  write(*,*) 'Efficiency will be evaluated considering these problems only.'

  write(*,*)
  write(*,*) 'This is the performance profile in terms of number of CPU time'
  
  do j = 1,nprob
     if ( compandeqf(j) ) then
        metric(1:nmeth,j) = ( 1.0d0 * time(1:nmeth,j) ) / minval( time(1:nmeth,j) )
        if ( all( time(1:nmeth,j) .ge. 1.0d0 ) .and. any( metric(1:nmeth,j) .ge. 10.0d0 ) ) then
           write(*,*) pname(1,j),time(1,j),metric(1,j),time(2,j),metric(2,j)
        end if
     else
        metric(1:nmeth,j) = huge( 1.0d0 )
     end if
  end do
  
  open(900,file='pptime.dat')
  
  tau = 1.0d0
  do while ( tau .le. 2.0d+04 )
     do i = 1,nmeth
        howmany(i) = count( metric(i,1:nprob) .le. tau )
     end do

     write(900,200) tau,( 1.0d0 * howmany(1:nmeth) ) / count( compandeqf(1:nprob) )

     if ( abs( tau -     1.0d0 ) .le. 5.0d-02 .or. &
          abs( tau -     1.1d0 ) .le. 5.0d-02 .or. &
          abs( tau -     1.2d0 ) .le. 5.0d-02 .or. &
          abs( tau -     1.3d0 ) .le. 5.0d-02 .or. &
          abs( tau -     1.4d0 ) .le. 5.0d-02 .or. &
          abs( tau -     1.5d0 ) .le. 5.0d-02 .or. &
          abs( tau -     1.6d0 ) .le. 5.0d-02 .or. &
          abs( tau -     1.7d0 ) .le. 5.0d-02 .or. &
          abs( tau -     1.8d0 ) .le. 5.0d-02 .or. &
          abs( tau -     1.9d0 ) .le. 5.0d-02 .or. &
          abs( tau -     2.0d0 ) .le. 5.0d-02 .or. &
          abs( tau -     3.0d0 ) .le. 5.0d-02 .or. &
          abs( tau -     4.0d0 ) .le. 5.0d-02 .or. &
          abs( tau -     5.0d0 ) .le. 5.0d-02 .or. &
          abs( tau -    10.0d0 ) .le. 5.0d-02 .or. &
          abs( tau -    15.0d0 ) .le. 5.0d-02 .or. &
          abs( tau -    20.0d0 ) .le. 5.0d-02 .or. &
          abs( tau -    30.0d0 ) .le. 5.0d-02 .or. &
          abs( tau -    40.0d0 ) .le. 5.0d-02 .or. &
          abs( tau -    50.0d0 ) .le. 5.0d-02 .or. &
          abs( tau -   100.0d0 ) .le. 5.0d-02 .or. &
          abs( tau -   500.0d0 ) .le. 5.0d-02 .or. &
          abs( tau -  1000.0d0 ) .le. 5.0d-02 .or. &
          abs( tau -  1500.0d0 ) .le. 5.0d-02 .or. &
          abs( tau -  5000.0d0 ) .le. 5.0d-02 .or. &
          abs( tau -  7500.0d0 ) .le. 5.0d-02 .or. &
          abs( tau - 10000.0d0 ) .le. 5.0d-02 .or. &
          abs( tau - 15000.0d0 ) .le. 5.0d-02 .or. &
          abs( tau - 20000.0d0 ) .le. 5.0d-02 ) then
        write(*,200) tau,( 1.0d0 * howmany(1:nmeth) ) / count( compandeqf(1:nprob) )
     end if
     
     tau = tau + 1.0d-01
  end do

  close(900)

  stop

!!$100 format(1X,A10,2(1X,1P,D24.16,2(1X,1P,D7.1)))
200 format(F8.2,16(1X,F10.6))
  
end program analyse