ogincl Subroutine

private subroutine ogincl(me, Inc)

ADDS CONSTRAINT TO ACTIVE SET AND REDUCES DERIVATIVES

2008/01/16 | J. SCHOENMAEKERS | NEW

Type Bound

optgra

Arguments

Type IntentOptional Attributes Name
class(optgra), intent(inout) :: me
integer(kind=ip), intent(in) :: Inc

CONSTRAINT TO BE INCLUDED

Calls

proc~~ogincl~~CallsGraph proc~ogincl optgra%ogincl proc~ogwrit optgra%ogwrit proc~ogincl->proc~ogwrit

Called by

proc~~ogincl~~CalledByGraph proc~ogincl optgra%ogincl proc~ogcorr optgra%ogcorr proc~ogcorr->proc~ogincl proc~ogopti optgra%ogopti proc~ogopti->proc~ogincl proc~ogexec optgra%ogexec proc~ogexec->proc~ogcorr proc~ogexec->proc~ogopti

Source Code

   subroutine ogincl(me,Inc)

      !! ADDS CONSTRAINT TO ACTIVE SET AND REDUCES DERIVATIVES
      !!
      !! 2008/01/16 | J. SCHOENMAEKERS | NEW

      class(optgra),intent(inout) :: me
      integer(ip),intent(in) :: Inc !! CONSTRAINT TO BE INCLUDED

      real(wp) :: val , fac , gam , sav , max
      integer(ip) :: row , col , ind , lst
      character(len=str_len) :: str

      ! GENERAL
      me%Numact = me%Numact + 1

      ! PERMUTATION TO GET ZERO DERIVATIVES AT END FOR NEW ACTIVE CONSTRAINT
      lst = me%Numvar
      do col = me%Numvar , me%Numact , -1
         if ( me%Conred(Inc,col)==0.0_wp ) then
            if ( col/=lst ) then
               do row = 1 , me%Numcon + 3
                  if ( me%Conact(row)<=0 ) then
                     sav = me%Conred(row,col)
                     me%Conred(row,col) = me%Conred(row,lst)
                     me%Conred(row,lst) = sav
                  end if
               end do
            end if
            lst = lst - 1
         end if
      end do

      ! PERMUTATION TO GET MAXIMUM PIVOT
      ind = me%Numact
      max = abs(me%Conred(Inc,ind))
      do col = me%Numact + 1 , lst
         val = abs(me%Conred(Inc,col))
         if ( val>max ) then
            ind = col
            max = val
         end if
      end do

      if ( ind/=me%Numact ) then
         do row = 1 , me%Numcon + 3
            if ( me%Conact(row)<=0 ) then
               sav = me%Conred(row,ind)
               me%Conred(row,ind) = me%Conred(row,me%Numact)
               me%Conred(row,me%Numact) = sav
            end if
         end do
      end if

      ! UPDATE LIST OF ACTIVE CONSTRAINTS
      me%Actcon(me%Numact) = Inc
      me%Conact(Inc) = me%Numact

      ! REDUCE FOR NEW ACTIVE CONSTRAINT
      if ( abs(me%Conred(Inc,me%Numact))<1.0e-12_wp ) then
         write (str,*) "OGINCL-WARNING: CONSTRAINT SINGULAR"
         call me%ogwrit(2,str)
         write (str,*) "INC=" , Inc
         call me%ogwrit(2,str)
         write (str,*) "PIV=" , me%Conred(Inc,me%Numact)
         call me%ogwrit(2,str)
         me%Numact = me%Numact - 1
         me%Conact(Inc) = 0
         return
      end if

      val = sqrt(sum(me%Conred(Inc,me%Numact:lst)**2))
      if ( me%Conred(Inc,me%Numact)>0.0_wp ) val = -val

      me%Conred(Inc,me%Numact) = me%Conred(Inc,me%Numact) - val

      sav = me%Conred(Inc,me%Numact)
      fac = 1.0_wp/sav
      me%Conred(Inc,me%Numact:lst) = me%Conred(Inc,me%Numact:lst)*fac

      fac = sav/val
      do row = 1 , me%Numcon + 3
         if ( me%Conact(row)<=0 ) then
            gam = dot_product(me%Conred(row,me%Numact:lst),me%Conred(Inc,me%Numact:lst))
            if ( gam/=0.0_wp ) then
               gam = gam*fac
               me%Conred(row,me%Numact:lst) = me%Conred(row,me%Numact:lst) + &
                                              me%Conred(Inc,me%Numact:lst)*gam
            end if
         end if
      end do

      me%Conred(Inc,me%Numact) = val
      me%Conred(Inc,me%Numact+1:lst) = 0.0_wp

   end subroutine ogincl