VB.Net Program to Find a Function Minimum

Using the Powell Search Method

by Namir Shammas

The following program calculates the minimum point of a multi-variable function using the Powell search method.

Click here to download a ZIP file containing the project files for this program.

The program prompts you to either use the predefined default input values or to enter the following for each variable:

1. The function tolerance value

2. The step tolerance value.

3. The maximum number of search cycles

4. The values for the initial set of variables

In case you choose the default input values, the program displays these values and proceeds to find the optimum point. In the case you select being prompted, the program displays the name of each input variable along with its default value. You can then either enter a new value or simply press Enter to use the default value. This approach allows you to quickly and efficiently change only a few input values if you so desire.

The program displays the following final results:

1. The coordinates of the minimum value.

2. The minimum function value.

3. The number of iterations

The current code finds the minimum for the following function:

f(x1,x2) = x1 - x2 + 2 * x1 ^ 2 + 2 * x1 * x2 + x2 ^ 2

Using an initial value of 0 for each variable, a function tolerance of 1e-7, and a maximum number of 100 search cycles. Here is the sample console screen:

Here is the listing for the main module. The module contains several test functions:

Module Module1

  Sub Main()
    Dim nNumVars As Integer = 2
    Dim fX() As Double = {0, 0}
    Dim fParam() As Double = {0, 0}
    Dim nIter As Integer = 0
    Dim nMaxCycles As Integer = 100
    Dim fEpsFx As Double = 0.0000001
    Dim fEpsStep As Double = 0.0000001
    Dim I As Integer
    Dim fBestF As Double
    Dim sAnswer As String, sErrorMsg As String = ""
    Dim oOpt As CPowellSearch1
    Dim MyFx As MyFxDelegate = AddressOf Fx3
    Dim SayFx As SayFxDelegate = AddressOf SayFx3

    oOpt = New CPowellSearch1

    Console.WriteLine("Powell Search Optimization (version 1)")
    Console.WriteLine("Finding the minimum of function:")
    Console.WriteLine(SayFx())
    Console.Write("Use default input values? (Y/N) ")
    sAnswer = Console.ReadLine()
    If sAnswer.ToUpper() = "Y" Then
      For I = 0 To nNumVars - 1
        Console.WriteLine("X({0}) = {1}", I + 1, fX(I))
      Next
      Console.WriteLine("Function tolerance = {0}", fEpsFx)
      Console.WriteLine("Step tolerance = {0}", fEpsStep)
      Console.WriteLine("Maxumum cycles = {0}", nMaxCycles)
    Else
      For I = 0 To nNumVars - 1
        fX(I) = GetIndexedDblInput("X", I + 1, fX(I))
      Next
      fEpsFx = GetDblInput("Function tolerance", fEpsFx)
      fEpsStep = GetDblInput("Function tolerance", fEpsStep)
      nMaxCycles = GetDblInput("Maxumum cycles", nMaxCycles)
    End If

    Console.WriteLine("******** FINAL RESULTS *************")
    fBestF = oOpt.CalcOptim(nNumVars, fX, fParam, fEpsFx, fEpsStep, nMaxCycles, nIter, sErrorMsg, MyFx)
    If sErrorMsg.Length > 0 Then
      Console.WriteLine("**ERROR: {0} ***", sErrorMsg)
    End If
    Console.WriteLine("Optimum at")
    For I = 0 To nNumVars - 1
      Console.WriteLine("X({0}) = {1}", I + 1, fX(I))
    Next
    Console.WriteLine("Function value = {0}", fBestF)
    Console.WriteLine("Number of iterations = {0}", nIter)
    Console.WriteLine()
    Console.Write("Press Enter to end the program ...")
    Console.ReadLine()
  End Sub

  Function GetDblInput(ByVal sPrompt As String, ByVal fDefInput As Double) As Double
    Dim sInput As String

    Console.Write("{0}? ({1}): ", sPrompt, fDefInput)
    sInput = Console.ReadLine()
    If sInput.Trim().Length > 0 Then
      Return Double.Parse(sInput)
    Else
      Return fDefInput
    End If
  End Function

  Function GetIntInput(ByVal sPrompt As String, ByVal nDefInput As Integer) As Integer
    Dim sInput As String

    Console.Write("{0}? ({1}): ", sPrompt, nDefInput)
    sInput = Console.ReadLine()
    If sInput.Trim().Length > 0 Then
      Return Double.Parse(sInput)
    Else
      Return nDefInput
    End If
  End Function

  Function GetIndexedDblInput(ByVal sPrompt As String, ByVal nIndex As Integer, ByVal fDefInput As Double) As Double
    Dim sInput As String

    Console.Write("{0}({1})? ({2}): ", sPrompt, nIndex, fDefInput)
    sInput = Console.ReadLine()
    If sInput.Trim().Length > 0 Then
      Return Double.Parse(sInput)
    Else
      Return fDefInput
    End If
  End Function

  Function GetIndexedIntInput(ByVal sPrompt As String, ByVal nIndex As Integer, ByVal nDefInput As Integer) As Integer
    Dim sInput As String

    Console.Write("{0}({1})? ({2}): ", sPrompt, nIndex, nDefInput)
    sInput = Console.ReadLine()
    If sInput.Trim().Length > 0 Then
      Return Double.Parse(sInput)
    Else
      Return nDefInput
    End If
  End Function

  Function SayFx1() As String
    Return "F(X) = 10 + (X(1) - 2) ^ 2 + (X(2) + 5) ^ 2"
  End Function

  Function Fx1(ByVal N As Integer, ByRef X() As Double, ByRef fParam() As Double) As Double
    Return 10 + (X(0) - 2) ^ 2 + (X(1) + 5) ^ 2
  End Function

  Function SayFx2() As String
    Return "F(X) = 100 * (X(1) - X(2) ^ 2) ^ 2 + (X(2) - 1) ^ 2"
  End Function

  Function Fx2(ByVal N As Integer, ByRef X() As Double, ByRef fParam() As Double) As Double
    Return 100 * (X(0) - X(1) ^ 2) ^ 2 + (X(1) - 1) ^ 2
  End Function

  Function SayFx3() As String
    Return "F(X) = X(1) - X(2) + 2 * X(1) ^ 2 + 2 * X(1) * X(2) + X(2) ^ 2"
  End Function

  Function Fx3(ByVal N As Integer, ByRef X() As Double, ByRef fParam() As Double) As Double
    Return X(0) - X(1) + 2 * X(0) ^ 2 + 2 * X(0) * X(1) + X(1) ^ 2
  End Function

  ' X(0) - X(1) + 2 * X(0) ^ 2 + 2 * X(0) * X(1) + X(1) ^ 2
End Module

Notice that the user-defined functions have accompanying helper functions to display the mathematical expression of the function being optimized. For example, function Fx1 has the helper function SayFx1 to list the function optimized in Fx1. Please observe the following rules::

When you add new user-defined functions to optimize, you must also add an accompanying helper function that states the function mathematical make up.
When you edit an existing user-defined functions to optimize, you must also edit its accompanying helper function to reflect the changes you make in the function mathematical make up
When you select a new user-defined function you must alter both variables MyFx, and SayFx to specify the new function and its accompanying helper function.

The program uses the following class to optimize the objective function:

Public Delegate Function MyFxDelegate(ByVal nNumVars As Integer, ByRef fX() As Double, ByRef fParam() As Double) As Double
Public Delegate Function SayFxDelegate() As String

Public Class CPowellSearch1

  Dim m_MyFx As MyFxDelegate

  Protected Function MyFxEx(ByVal nNumVars As Integer, ByRef fX() As Double, ByRef fParam() As Double, _
                  ByRef fSmat(,) As Double, ByVal nII As Integer, _
                  ByVal fAlpha As Double) As Double
    Dim I As Integer
    Dim fXX(nNumVars) As Double

    For I = 0 To nNumVars - 1
      fXX(I) = fX(I) + fAlpha * fSmat(nII, I)
    Next I

    MyFxEx = m_MyFx(nNumVars, fXX, fParam)
  End Function

  Protected Function LinSearch_DirectSearch( _
            ByVal nNumVars As Integer, ByRef fX() As Double, ByRef fParam() As Double, _
            ByRef fAlpha As Double, ByRef fSmat(,) As Double, ByVal nII As Integer, _
            ByVal fInitStep As Double, ByVal fMinStep As Double) As Boolean
    Dim F1, F2 As Double

    F1 = MyFxEx(nNumVars, fX, fParam, fSmat, nII, fAlpha)

    Do
      F2 = MyFxEx(nNumVars, fX, fParam, fSmat, nII, fAlpha + fInitStep)
      If F2 < F1 Then
        F1 = F2
        fAlpha += fInitStep
      Else
        F2 = MyFxEx(nNumVars, fX, fParam, fSmat, nII, fAlpha - fInitStep)
        If F2 < F1 Then
          F1 = F2
          fAlpha -= fInitStep
        Else
          ' reduce search step size
          fInitStep /= 10
        End If
      End If
    Loop Until fInitStep < fMinStep

    Return True

  End Function


  Public Function CalcOptim(ByVal nNumVars As Integer, ByRef fX() As Double, ByRef fParam() As Double, _
                    ByVal fEpsFx As Double, ByVal fEpsStep As Double, ByVal nMaxCycles As Integer, _
                    ByRef nIter As Integer, ByRef sErrorMsg As String, _
                    ByVal MyFx As MyFxDelegate) As Double

    Dim I, J, K, L As Integer
    Dim fAlpha, fSum, F1, F2 As Double
    Dim fX1(nNumVars) As Double, fX2(nNumVars) As Double
    Dim fSmat(nNumVars + 1, nNumVars) As Double

    m_MyFx = MyFx
    F1 = MyFx(nNumVars, fX, fParam)
    ' Build fSmat matrix
    For K = 0 To nNumVars - 1
      fX1(K) = fX(K)
      For L = 0 To nNumVars - 1
        fSmat(K, L) = IIf(K = L, 1, 0)
      Next L
    Next K

    sErrorMsg = ""
    nIter = 1
    J = 1

    Do

      nIter += 1

      ' reset row fSmat(nNumVars+1,:)
      For K = 0 To nNumVars - 1
        fSmat(nNumVars, K) = 0
      Next K

      For I = 0 To nNumVars - 1
        fAlpha = 0
        If Not LinSearch_DirectSearch(nNumVars, fX, fParam, fAlpha, fSmat, I, 0.1, 0.0001) Then
          sErrorMsg = "Linear Search failed"
          Exit Do
        End If
        For K = 0 To nNumVars - 1
          fX(K) += fAlpha * fSmat(I, K)
          fSmat(nNumVars, K) = fSmat(nNumVars, K) + fAlpha * fSmat(I, K)
        Next K
      Next I

      fAlpha = 0
      If Not LinSearch_DirectSearch(nNumVars, fX, fParam, fAlpha, fSmat, nNumVars, 0.1, 0.0001) Then
        sErrorMsg = "Linear Search failed"
        Exit Do
      End If

      For K = 0 To nNumVars - 1
        fX(K) += fAlpha * fSmat(nNumVars, K)
        fX2(K) = fX(K)
      Next K
      F2 = MyFx(nNumVars, fX2, fParam)

      ' test end of iterations criteria
      If Math.Abs(F2 - F1) < fEpsFx Then Exit Do
      fSum = 0
      For K = 0 To nNumVars - 1
        fSum += (fX2(K) - fX1(K)) ^ 2
      Next K
      fSum = Math.Sqrt(fSum)
      If fSum < fEpsStep Then Exit Do
      If J >= nMaxCycles Then Exit Do

      J += 1

      ' rotate data
      For K = 0 To nNumVars - 1
        fX1(K) = fX2(K)
      Next K
      F1 = F2

      ' copy fSmat matrix
      For K = 0 To nNumVars - 1
        For L = 0 To nNumVars - 1
          fSmat(K, L) = fSmat(K + 1, L)
        Next L
      Next K

    Loop

    Return F1

  End Function

End Class

BACK