pen01.f90

Go to the documentation of this file.

!> @file pen01.f90
!! @ingroup FORT1THREAD_EXAMPLES
!!
!! Test model pen01.gms -- a model with penalty pairs.
!!
!! \f$\max x1\f$
!! <table>
!! <tr><th>     <th>x1  <th>x2  <th>x3   <th>x4  <th>x5  <th>x6  <th>x7 <th>
!! <tr><td>1:   <td>1   <td>-1  <td>-1   <td>    <td>    <td>    <td>   <td>=   0
!! <tr><td>2:   <td>    <td>1   <td> -10 <td>-5  <td>    <td>    <td>   <td>=   0
!! <tr><td>3:   <td>    <td>    <td>     <td>-1  <td>    <td>1   <td>1  <td>\f$\leq\f$  10
!! <tr><td>4:   <td>    <td>    <td>     <td>    <td>-1  <td>1   <td>-1 <td>\f$\leq\f$ -10
!! <tr><td>5:   <td>    <td>    <td>1    <td>    <td>    <td>    <td>   <td>\f$\leq\f$   8
!! <tr><td>     <td>    <td>x2  <td>x3   <td>x4  <td>x5  <td>x6  <td>x7 <td>\f$\geq\f$ 0
!! </table>
!!
!! The base-case has a penalty-pair with two structural variables and
!! one with a structural and a slack.
!!
!! In case 2 variable x6 and x7 have an upper bound and can therefore
!! not be used a a penalty variable. However, the effect is that x5
!! becomes implied free and can be used as post-triangular.
!!
!! In case 3 variable x6 and x7 have different upper bounds and can again
!! not be used a a penalty variable, and this time the slack is used
!! so we have two structural/slack-pairs.
!!
!!
!! 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
 
!> Main program. A simple setup and call of CONOPT
!!
Program pen01
 
   Use proginfo
   Use conopt
   Use casedata_num
   implicit None
!
!  Declare the user callback routines as Integer, External:
!
   Integer, External :: pen01_readmatrix ! Mandatory Matrix definition routine defined below
   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
#ifdef dec_directives_win32
!DEC$ ATTRIBUTES STDCALL, REFERENCE, NOMIXED_STR_LEN_ARG :: Pen01_ReadMatrix
!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
#endif
!
!  Control vector
!
   INTEGER, Dimension(:), Pointer :: cntvect
   INTEGER :: coi_error
 
!
!  Create and initialize a Control Vector
!
   call startup
 
   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, 7 ) )  ! # variables
   coi_error = max( coi_error, coidef_numcon( cntvect, 5 ) )  ! # constraints
   coi_error = max( coi_error, coidef_numnz( cntvect,13 ) )  ! # nonzeros in the Jacobian
   coi_error = max( coi_error, coidef_numnlnz( cntvect, 0 ) )  ! # 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 1
   coi_error = max( coi_error, coidef_optfile( cntvect, 'pen01.opt' ) )
!
!  Tell CONOPT about the callback routines:
!
   coi_error = max( coi_error, coidef_readmatrix( cntvect, pen01_readmatrix ) )
   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     ) )
 
#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 CONOPT4 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.
!
!  Start CONOPT:
!
   casenum = 1
   coi_error = coi_solve( cntvect )
   If ( coi_error /= 0 ) then
      call flog( "Case 1: Errors encountered during solution", 1 )
   elseif ( stacalls == 0 .or. solcalls == 0 ) then
      call flog( "Case 1: Status or Solution routine was not called", 1 )
   elseif ( sstat /= 1 ) then
      call flog( "Case 1: Solver Status was 1 not as expected.", 1 )
   endif
!
!  Solve the second case:
!
   stacalls = 0; solcalls = 0;
   casenum = 2
   coi_error = coi_solve( cntvect )
   If ( coi_error /= 0 ) then
      call flog( "Case 2: Errors encountered during solution", 1 )
   elseif ( stacalls == 0 .or. solcalls == 0 ) then
      call flog( "Case 2: Status or Solution routine was not called", 1 )
   elseif ( sstat /= 1 ) then
      call flog( "Case 2: Solver Status was 1 not as expected.", 1 )
   endif
!
!  And the third:
!
   stacalls = 0; solcalls = 0;
   casenum = 3
   coi_error = coi_solve( cntvect )
   If ( coi_error /= 0 ) then
      call flog( "Case 3: Errors encountered during solution", 1 )
   elseif ( stacalls == 0 .or. solcalls == 0 ) then
      call flog( "Case 3: Status or Solution routine was not called", 1 )
   elseif ( sstat /= 1 ) then
      call flog( "Case 3: Solver Status was 1 not as expected.", 1 )
   endif
 
   write(*,*)
   write(*,*) 'End of pen01 example. Return code=',coi_error
 
   if ( coi_free(cntvect) /= 0 ) call flog( "Error while freeing control vector",1)
 
   call flog( "Successful Solve", 0 )
!
!  Free solution memory
!
   call finalize
 
End Program pen01
!
! ============================================================================
!  Define information about the model:
!
 
!>  Define information about the model
!!
!! @include{doc} readMatrix_params.dox
Integer Function pen01_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 :: Pen01_ReadMatrix
#endif
   use casedata_num
   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.
!
   lower(2) = 0.0d0
   lower(3) = 0.0d0
   lower(4) = 0.0d0
   lower(5) = 0.0d0
   lower(6) = 0.0d0
   lower(7) = 0.0d0
   if ( casenum == 2 ) Then
      upper(6) = 12.0d0
      upper(7) =  7.0d0
   Endif
   if ( casenum == 3 ) Then
      upper(6) =  7.0d0
      upper(7) = 12.0d0
   Endif
!
!  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.
!
   type(1) = 0
   type(2) = 0
   type(3) = 2
   type(4) = 2
   type(5) = 2
   rhs(3) =  10.d0
   rhs(4) = -10.d0
   rhs(5) =   8.d0
!
!  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
!       x1  x2  x3  x4  x5  x6  x7
!  1:    1   2   4
!  2:        3       6   8
!  3:                7      10  12
!  4:                    9  11  13
!  5:            5
!
   colsta(1) = 1
   colsta(2) = 2
   colsta(3) = 4
   colsta(4) = 6
   colsta(5) = 8
   colsta(6) =10
   colsta(7) =12
   colsta(8) =14
   rowno(1)  = 1
   rowno(2)  = 1
   rowno(3)  = 2
   rowno(4)  = 1
   rowno(5)  = 5
   rowno(6)  = 2
   rowno(7)  = 3
   rowno(8)  = 2
   rowno(9)  = 4
   rowno(10) = 3
   rowno(11) = 4
   rowno(12) = 3
   rowno(13) = 4
!
!  Nonlinearity Structure: Not needed for a linear model
!
!  Value (Linear only)
!       x1  x2  x3  x4  x5  x6  x7
!  1:    1  -1  -1
!  2:        1     -10  -5
!  3:               -1       1   1
!  4:                   -1   1  -1
!  5:            1
!
   value(1)  =  1.d0
   value(2)  = -1.d0
   value(3)  =  1.d0
   value(4)  = -1.d0
   value(5)  =  1.d0
   value(6)  =-10.d0
   value(7)  = -1.d0
   value(8)  = -5.d0
   value(9)  = -1.d0
   value(10) =  1.d0
   value(11) =  1.d0
   value(12) =  1.d0
   value(13) = -1.d0
 
  pen01_readmatrix = 0 ! Return value means OK
 
end Function pen01_readmatrix