Solver, Minimization, 2 Variables

Question

Solver, Minimization, 2 Variables

ikkebins (TechnicalUser)

(OP)

28 May 14 07:58

Hey there,

following problem:

A set of non linear equations, that can be solved by determining two variables.

So far i have some experience using the NEQNF Solver, that minimizes an objective function (or several functions) by varying one variable.

Now i have 7 objective functions (non linear) and i need to vary and determine two variables to solve minimize the functions. Is there an solver that you can suggest me to use for this problem?

I`m using fortrann 77.

Thanks for your help,
ikkebins

salgerman (Programmer) · Answer 1 · 2014-05-28 19:48:41

I don't know; but, have you looked for Optimization in netlib?

ikkebins (TechnicalUser) · Answer 2 · 2014-06-05 03:59:25

RE: Solver, Minimization, 2 Variables

ikkebins (TechnicalUser)

(OP)

5 Jun 14 03:59

Yes, but i couldnt find anything

mikrom (Programmer) · Answer 3 · 2014-06-05 16:31:07

If you cannot find nothing appropriate, you can try to code it yourself.
For solving the nonlinear equation system of this form

f(x) = 0

we can use Newton method:

x₁ = x₀ - J^-1(f, x₀) * f(x₀)

That means:

1) For the given starting point x₀ we need to do:

1.a) Compute the Jacobian matrix J(f, x₀) of the vector function f(x) in the point x₀

1.b) Solve the linear equation system
J(f, x₀) * d = - f(x₀)
For this step we can use a linear equation solver from Lapack.

2) finaly compute the next approximation point
x₁ = x₀ + d

All steps should be repeated until the desired accuracy is reached, it means for given eps -> 0 or delta -> 0 until
||x_n - x_n-1|| < delta or ||f(x_n)|| < eps

IMO it doesn't seem very complicated...

mikrom (Programmer) · Answer 4 · 2014-07-02 12:42:55

I looked at the description of NEQNF method here
http://faculty.ksu.edu.sa/Almutaz/Documents/Summer...
I took as example the 3x3 system of nonlinear equations from the document above and tried to solve it with simple Newton method.
For solving of linear equation systems I first used self written GEM method, but then - to make the code shorter for this posting - I took a method from LAPACk:

nlsolve.f95

CODE

module linear_algebra
  implicit none
  ! This module contains the tools
  ! for Linear Algebra
contains
  real function norm(X)
    real, dimension(:) :: X
    ! short notation of Fortran 90
    norm = sqrt(sum(X*X))
  end function norm

  subroutine Lsolve(M, V, X, n)
    ! LAPACK routine for solving linear equation system
    ! M * X = V
    integer, intent(in) :: n
    real, dimension(n,n), intent(inout) :: M
    real, dimension(n), intent(inout) :: V
    real, dimension(n), intent(out) :: X 
    !
    real :: A(n,n), b(n)
    integer :: i, pivot(n), ok
    A = M
    b = V
    call SGESV(n, 1, A, n, pivot, b, 3, ok)    
    X = b
  end subroutine Lsolve
end module linear_algebra

module solving_methods
  use linear_algebra
  implicit none
  ! This module contains the solving method
contains
  subroutine jacobian(F, P, Jcb, n)
    integer :: n
    ! interface for functional argument
    interface
      function F(X, n)
        integer :: n
        real, dimension(n), intent(in) :: x
        real, dimension(n) :: f
      end function f
    end interface    
    real, dimension(n), intent(in) :: P
    real, dimension(n,n), intent(out) :: Jcb
    !
    integer, parameter :: info = 0 ! view informations
    real, parameter :: h=0.1   ! step    
    integer :: i, j
    real, dimension(n) :: diffj, T1, T2
    do i=1, n   ! nr of functions
      do j=1, n ! nr of variables
        T1 = P
        T2 = P
        T1(j) = P(j) + h
        T2(j) = P(j) - h

        diffj = (F(T1, n) - F(T2, n)) / (2*h)
        if (info .ne. 0) then
           write(*,*) 'T1=', T1
           write(*,*) 'T2=', T2
           write(*,*) 'diffj = ', diffj
        end if
        ! difference of i-th function on j-th variable
        Jcb(i,j) = diffj(i)
      end do
    end do
  end subroutine jacobian

  subroutine newton(F, X, n)
    ! Newton-Raphson Method for solving nonlinear systems
    integer :: n
    ! interface for functional argument
    interface
      function F(X, n)
        integer :: n
        real, dimension(n), intent(in) :: x
        real, dimension(n) :: f
      end function f
    end interface    
    real, dimension(n), intent(inout) :: X
    !
    integer, parameter :: info = 1 ! view informations
    integer, parameter :: max_steps = 100
    real, parameter :: delta=1e-5, eps=1e-8
    real, dimension(n,n) :: JF
    real, dimension(n) :: T, D, X_old
    real :: dnorm, fnorm
    integer :: k
    
    do k=1, max_steps
      X_old = X
      T = F(X, n)
      call Jacobian(F,X,JF,n)
      call Lsolve(JF,T,D,n)
      X = X - D
      dnorm = norm(D)
      fnorm = norm(T)
      if ((fnorm <= eps) .or. (dnorm <= delta)) then
        ! exit loop 
        exit
      end if      
    end do
    if (info .ne. 0) then
      write(*,*) 'fnorm =', fnorm
      write(*,*) 'dnorm =', dnorm
      write(*,*) 'nr steps =', k
    end if    
  end subroutine newton
end module solving_methods

module user_functions
  implicit none
  ! Define here your own right-hand function
  ! for the non linear system
  !    F(X) = 0
contains
  function nls1(X, n)
    ! starting point: 4.0, 4.0, 4.0
    ! solution      : 1.0 2.0 3.0
    integer :: n
    real, dimension(n), intent(in) :: X
    real, dimension(n) :: nls1
    !
    NLS1(1) = X(1) + EXP(X(1)-1.0) + (X(2)+X(3))*(X(2)+X(3)) - 27.0
    NLS1(2) = EXP(X(2)-2.0)/X(1) + X(3)*X(3) - 10.0
    NLS1(3) = X(3) + SIN(X(2)-2.0) + X(2)*X(2) - 7.0
  end function nls1
end module user_functions

program nlsys
  use solving_methods
  use user_functions
  implicit none
  
  real, dimension(3) :: X

  write(*,*) 'Solving Nonlinear System:'
  x =(/4,4,4/)
  write(*,*) 'Starting Point =', X
  call newton(nls1, X, 3)
  write(*,*) 'Solution =', X
  write(*,*)
end program nlsys

Here are the results:

CODE

$ gfortran -o nlsolve nlsolve.f95 -L. liblapack.dll

$ nlsolve
 Solving Nonlinear System:
 Starting Point =   4.00000000       4.00000000       4.00000000
 fnorm =   1.33280039E-07
 dnorm =   5.93307874E-08
 nr steps =           7
 Solution =   1.00000000       2.00000000       3.00000000

Come Join Us!

Posting Guidelines

Jobs

Solver, Minimization, 2 Variables