mono07.f90

Go to the documentation of this file.

!> @file mono07.f90
!! @ingroup FORT1THREAD_EXAMPLES
!!
!! 
!! Monotone function to bound conversion example 07
!! 
!! This is a CONOPT implementation of the GAMS model:
!! 
!! @verbatim
!! variable x1
!! equation e1;
!! 
!! e1 .. logistic(x) =R= C;
!! where logistic is defined as:
!! x >= 0: 1/(1+exp(-x))
!! x <= 0: exp(x)/(1+exp(x))
!! 
!! model mono / all /;
!! solve mono using nlp maximizing x1;
!! @endverbatim
!! 
!! Where we have the following 9 cases:
!! 1:  =E= -2.0  ! Infeasible
!! 2:  =E=  0.75 ! Feasible, x1 = log(2)
!! 3:  =E=  1.1  ! Infeasible
!! 4:  =L= -2.0  ! Infeasible
!! 5:  =L=  0.75 ! Feasible, x1 = log(2)
!! 6:  =L=  1.1  ! Unbounded
!! 7:  =G= -2.0  ! Unbounded
!! 8:  =G=  0.75 ! Unbounded
!! 9:  =G=  1.1  ! Infeasible
!! 
!!
!! For more information about the individual callbacks, please have a look at the source code.
 
#if defined(_WIN32) && !defined(_WIN64)
#define dec_directives_win32
#endif
 
module mono07data
   Integer, Parameter :: MaxCase = 9
   real*8, Parameter, dimension(MaxCase) :: caserhs = &
     (/ -2.0d0, 0.75d0, 1.1d0, -2.0d0, 0.75d0, 1.1d0, -2.0d0, 0.75d0, 1.1d0 /)
   Integer, Parameter, dimension(MaxCase) :: casetype = &
     (/ 0, 0, 0, 2, 2, 2, 1, 1, 1 /)
   Integer, Parameter, dimension(MaxCase) :: casemstat = &
     (/ 4, 1, 4, 4, 1, 3, 3, 3, 4 /)
   real*8, Parameter, dimension(MaxCase) :: caseobj = &
     (/ 0.d0, log(3.0d0), 0.0d0, 0.0d0, log(3.0d0), 0.0d0, 0.0d0, 0.0d0, 0.0d0 /)
   Integer :: casenum
end module mono07data
 
!> Main program. A simple setup and call of CONOPT
!!
Program mono07
 
   Use proginfo
   Use conopt
   Use mono07data
   implicit None
!
!  Declare the user callback routines as Integer, External:
!
   Integer, External :: mono_readmatrix ! Mandatory Matrix definition routine defined below
   Integer, External :: mono_fdeval     ! Function and Derivative evaluation routine
                                        ! needed a nonlinear model.
   Integer, External :: mono_fdinterval ! Function and Derivative evaluation routine
                                       ! needed a nonlinear model.
   Integer, External :: std_status     ! Standard callback for displaying solution status
   Integer, External :: std_solution   ! Standard callback for displaying solution values
   Integer, External :: std_message    ! Standard callback for managing messages
   Integer, External :: std_errmsg     ! Standard callback for managing error messages
   Integer, External :: std_triord     ! Standard callback for Monongular order
#ifdef dec_directives_win32
!DEC$ ATTRIBUTES STDCALL, REFERENCE, NOMIXED_STR_LEN_ARG :: Mono_ReadMatrix
!DEC$ ATTRIBUTES STDCALL, REFERENCE, NOMIXED_STR_LEN_ARG :: Mono_FDEval
!DEC$ ATTRIBUTES STDCALL, REFERENCE, NOMIXED_STR_LEN_ARG :: Mono_FDInterval
!DEC$ ATTRIBUTES STDCALL, REFERENCE, NOMIXED_STR_LEN_ARG :: Std_Status
!DEC$ ATTRIBUTES STDCALL, REFERENCE, NOMIXED_STR_LEN_ARG :: Std_Solution
!DEC$ ATTRIBUTES STDCALL, REFERENCE, NOMIXED_STR_LEN_ARG :: Std_Message
!DEC$ ATTRIBUTES STDCALL, REFERENCE, NOMIXED_STR_LEN_ARG :: Std_ErrMsg
!DEC$ ATTRIBUTES STDCALL, REFERENCE, NOMIXED_STR_LEN_ARG :: Std_TriOrd
#endif
!
!  Control vector
!
   INTEGER, Dimension(:), Pointer :: cntvect
   INTEGER :: coi_error
 
   call startup
!
!  Create and initialize a Control Vector
!
   coi_error = coi_create( cntvect )
!
!  Tell CONOPT about the size of the model by populating the Control Vector:
!
   coi_error = max( coi_error, coidef_numvar( cntvect, 1 ) )  ! # variables
   coi_error = max( coi_error, coidef_numcon( cntvect, 1 ) )  ! # constraints
   coi_error = max( coi_error, coidef_numnz( cntvect, 1 ) )  ! # nonzeros in the Jacobian
   coi_error = max( coi_error, coidef_numnlnz( cntvect, 1 ) )  ! # of which are nonlinear
   coi_error = max( coi_error, coidef_optdir( cntvect, +1 ) ) ! Maximize
   coi_error = max( coi_error, coidef_objvar( cntvect, 1 ) )  ! Objective is variable 3
   coi_error = max( coi_error, coidef_optfile( cntvect, 'Mono07.opt' ) )
!
!  Tell CONOPT about the callback routines:
!
   coi_error = max( coi_error, coidef_readmatrix( cntvect, mono_readmatrix ) )
   coi_error = max( coi_error, coidef_fdeval( cntvect, mono_fdeval     ) )
   coi_error = max( coi_error, coidef_fdinterval( cntvect, mono_fdinterval ) )
   coi_error = max( coi_error, coidef_status( cntvect, std_status     ) )
   coi_error = max( coi_error, coidef_solution( cntvect, std_solution   ) )
   coi_error = max( coi_error, coidef_message( cntvect, std_message    ) )
   coi_error = max( coi_error, coidef_errmsg( cntvect, std_errmsg     ) )
   coi_error = max( coi_error, coidef_triord( cntvect, std_triord     ) )
 
#if defined(CONOPT_LICENSE_INT_1) && defined(CONOPT_LICENSE_INT_2) && defined(CONOPT_LICENSE_INT_3) && defined(CONOPT_LICENSE_TEXT)
   coi_error = max( coi_error, coidef_license( cntvect, conopt_license_int_1, conopt_license_int_2, conopt_license_int_3, conopt_license_text) )
#endif
 
   If ( coi_error .ne. 0 ) THEN
      write(*,*)
      write(*,*) '**** Fatal Error while loading CONOPT Callback routines.'
      write(*,*)
      call flog( "Skipping Solve due to setup errors", 1 )
   ENDIF
!
!  Save the solution so we can check the duals:
!
   do_allocate = .true.
   DO casenum = 1, maxcase
!
!  Start CONOPT:
!
      coi_error = coi_solve( cntvect )
 
      write(*,*)
      write(*,*) 'End of Mono07 example case',casenum,'. Return code=',coi_error
 
      If ( coi_error /= 0 ) then
         call flog( "Errors encountered during solution", 1 )
      elseif ( stacalls == 0 .or. solcalls == 0 ) then
         call flog( "Status or Solution routine was not called", 1 )
      elseif ( sstat /= 1 .or. mstat /= casemstat(casenum) ) then
         call flog( "Solver and Model Status was not as expected", 1 )
      elseif ( mstat == 1 .and. (abs( obj-caseobj(casenum) ) > 0.000001d0 ) ) then
         call flog( "Incorrect objective returned", 1 )
      Elseif ( mstat == 1 ) Then
         Call checkdual( 'Mono07', maximize )
      Elseif ( mstat == 4 ) Then
         Call checkdual( 'Mono07', infeasible )
      endif
 
   EndDo ! end Casenum loop
 
   if ( coi_free(cntvect) /= 0 ) call flog( "Error while freeing control vector",1)
 
   call flog( "Successful Solve", 0 )
 
End Program mono07
!
! ============================================================================
!  Define information about the model:
!
 
!>  Define information about the model
!!
!! @include{doc} readMatrix_params.dox
Integer Function mono_readmatrix( lower, curr, upper, vsta, type, rhs, esta,      &
                                 colsta, rowno, value, nlflag, n, m, nz,  &
                                 usrmem )
#ifdef dec_directives_win32
!DEC$ ATTRIBUTES STDCALL, REFERENCE, NOMIXED_STR_LEN_ARG :: Mono_ReadMatrix
#endif
   Use mono07data
   implicit none
   integer, intent (in)                      :: n      ! number of variables
   integer, intent (in)                      :: m      ! number of constraints
   integer, intent (in)                      :: nz     ! number of nonzeros
   real*8,  intent (in out), dimension(n)    :: lower  ! vector of lower bounds
   real*8,  intent (in out), dimension(n)    :: curr   ! vector of initial values
   real*8,  intent (in out), dimension(n)    :: upper  ! vector of upper bounds
   integer, intent (in out), dimension(n)    :: vsta   ! vector of initial variable status
                                                       ! (not defined here)
   integer, intent (out),    dimension(m)    :: type   ! vector of equation types
   integer, intent (in out), dimension(m)    :: esta   ! vector of initial equation status
                                                       ! (not defined here)
   real*8,  intent (in out), dimension(m)    :: rhs    ! vector of right hand sides
   integer, intent (in out), dimension(n+1)  :: colsta ! vector with start of column indices
   integer, intent (out),    dimension(nz)   :: rowno  ! vector of row numbers
   integer, intent (in out), dimension(nz)   :: nlflag ! vector of nonlinearity flags
   real*8,  intent (in out), dimension(nz)   :: value  ! vector of matrix values
   real*8   usrmem(*)                                  ! optional user memory
!
!  Information about Variables:
!  Default: Lower = -Inf, Curr = 0, and Upper = +inf.
!  Default: the status information in Vsta is not used.
!
!  The model uses defaults
!
!  Information about Constraints:
!  Default: Rhs = 0
!  Default: the status information in Esta and the function
!           value in FV are not used.
!  Default: Type: There is no default.
!           0 = Equality,
!           1 = Greater than or equal,
!           2 = Less than or equal,
!           3 = Non binding.
!
!  Constraint 1: e1
!  Rhs = 10.0 and type Less than or Equal
!
   rhs(1)  = caserhs(casenum)
   type(1) = casetype(casenum)
!
!  Information about the Jacobian. CONOPT expects a columnwise
!  representation in Rowno, Value, Nlflag and Colsta.
!
!  Colsta = Start of column indices (No Defaults):
!  Rowno  = Row indices
!  Value  = Value of derivative (by default only linear
!                                derivatives are used)
!  Nlflag = 0 for linear and 1 for nonlinear derivative
!           (not needed for completely linear models)
!
!  Indices
!       x(1)
!  1:    1
!
   colsta(1) =  1
   colsta(2) =  2
   rowno(1)  =  1
!
!  Nonlinearity Structure: L = 0 are linear and NL = 1 are nonlinear
!       x(1)
!  1:    NL
!
   nlflag(1) =  1
!
!  Value (Linear only)
!       x(1)
!  1:   NL
!
  mono_readmatrix = 0 ! Return value means OK
 
end Function mono_readmatrix
!
!==========================================================================
!  Compute nonlinear terms and non-constant Jacobian elements
!
 
!>  Compute nonlinear terms and non-constant Jacobian elements
!!
!! @include{doc} fdeval_params.dox
Integer Function mono_fdeval( x, g, jac, rowno, jcnm, mode, ignerr, errcnt, &
                              n, nz, thread, usrmem )
#ifdef dec_directives_win32
!DEC$ ATTRIBUTES STDCALL, REFERENCE, NOMIXED_STR_LEN_ARG :: Mono_FDEval
#endif
   implicit none
   integer, intent (in)                      :: n      ! number of variables
   integer, intent (in)                      :: rowno  ! number of the row to be evaluated
   integer, intent (in)                      :: nz     ! number of nonzeros in this row
   real*8,  intent (in),     dimension(n)    :: x      ! vector of current solution values
   real*8,  intent (in out)                  :: g      ! constraint value
   real*8,  intent (in out), dimension(n)    :: jac    ! vector of derivatives for current constraint
   integer, intent (in),     dimension(nz)   :: jcnm   ! list of variables that appear nonlinearly
                                                       ! in this row. Ffor information only.
   integer, intent (in)                      :: mode   ! evaluation mode: 1 = function value
                                                       ! 2 = derivatives, 3 = both
   integer, intent (in)                      :: ignerr ! if 1 then errors can be ignored as long
                                                       ! as errcnt is incremented
   integer, intent (in out)                  :: errcnt ! error counter to be incremented in case
                                                       ! of function evaluation errors.
   integer, intent (in)                      :: thread
   real*8   usrmem(*)                                  ! optional user memory
   real*8 :: e
!
!  Row 1: e1
!
   if ( rowno .eq. 1 ) then
!
!  Mode = 1 or 3. G = log(x1)
!
      e = exp(-abs(x(1)))
      if ( mode .eq. 1 .or. mode .eq. 3 ) then
         if ( x(1) >= 0.0d0 ) then
            g    = 1.0d0 / (1.0d0+e)
         else
            g    = e / (1.0d0+e)
         endif
      endif
!
!  Mode = 2 or 3: Derivative values:
!
      if ( mode .eq. 2 .or. mode .eq. 3 ) then
         jac(1) = e / (1.0d0+e)**2
      endif
      mono_fdeval = 0
   else
!
!  There are no other rows:
!
      mono_fdeval = 1
   endif
 
end Function mono_fdeval
 
 
!>  Evaluating nonlinear functions and derivatives on an interval. Used in preprocessing
!!
!! @include{doc} fdinterval_params.dox
Integer Function mono_fdinterval( XMIN, XMAX, GMIN, GMAX,  &
                                  JMIN, JMAX, ROWNO, JCNM, &
                                  MODE, PINF, N, NJ, USRMEM )
#ifdef dec_directives_win32
!DEC$ ATTRIBUTES STDCALL, REFERENCE, NOMIXED_STR_LEN_ARG :: Mono_FDInterval
#endif
   Implicit None
   INTEGER, Intent(IN)                   :: rowno, mode, n, nj
   INTEGER, Dimension(NJ), Intent(IN)    :: jcnm
   real*8,  Dimension(N), Intent(IN)     :: xmin, xmax
   real*8,  Intent(IN OUT)               :: gmin, gmax
   real*8,  Dimension(N), Intent(IN OUT) :: jmin, jmax
   real*8,  Intent(IN)                   :: pinf
   real*8,  Intent(IN OUT)               :: usrmem(*)
   real*8 :: emin, emax
!
!  Row 1: e1
!
   write(10,*) 'Enter Mono_FDInterval. Row=',rowno,' Mode=',mode
   write(10,*) 'Xmin=',xmin
   write(10,*) 'Xmax=',xmax
   if ( rowno .eq. 1 ) then
!
!  Mode = 1 or 3. G = log(x1)
!
      emin = exp(-abs(xmin(1)))
      emax = exp(-abs(xmax(1)))
      if ( mode .eq. 1 .or. mode .eq. 3 ) then
         If ( xmin(1) >= 0.0d0 ) then
            gmin = 1.0d0 / (1.0d0+emin)
         else
            gmin = emin / (1.d0+emin)
         endif
         If ( xmax(1) >= 0.0d0 ) then
            gmax = 1.0d0 / (1.0d0+emax)
         else
            gmax = emax / (1.d0+emax)
         endif
         write(10,*) 'Gmin=',gmin,' Gmax=',gmax
      endif
!
!  Mode = 2 or 3: Derivative values:
!
      if ( mode .eq. 2 .or. mode .eq. 3 ) then
         If ( xmin(1) >= 0.0d0 ) then
            jmin(1) = emax/(1.0d0+emax)**2
            jmax(1) = emin/(1.0d0+emin)**2
         else if ( xmax(1) <= 0.0d0 ) then
            jmin(1) = emin/(1.0d0+emin)**2
            jmax(1) = emax/(1.0d0+emax)**2
         else
            jmin(1) = min( emin/(1.0d0+emin)**2, emax/(1.0d0+emax)**2 )
            jmax(1) = 0.25d0
         endif
         write(10,*) 'Jmin=',jmin
         write(10,*) 'Jmax=',jmax
      endif
      mono_fdinterval = 0
   else
!
!  There are no other rows:
!
      mono_fdinterval = 1
   endif
 
end Function mono_fdinterval