NewtonMin2D Function

public function NewtonMin2D(tol, a, b, x, fobj) result(fmin)

Newton's iterative algorithm to find the minimimum of function f(x,y) in the interval x(i)[a(i),b(i)], using grad(f(x)=0 and the derivative

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(in) :: tol

abort tolerance
real(kind=wp), intent(in) :: a(2)

search interval (2D)
real(kind=wp), intent(in) :: b(2)

search interval (2D)
real(kind=wp), intent(inout) :: x(2)

initial guess on input, result on output
class(c_newton_Min2D) :: fobj

functional f(x,y) to minimize with f, f', f''

Return Value real(kind=wp)

on output =f(x,y)


Calls

proc~~newtonmin2d~~CallsGraph proc~newtonmin2d NewtonMin2D FR FR proc~newtonmin2d->FR dFR dFR proc~newtonmin2d->dFR ddFR ddFR proc~newtonmin2d->ddFR

Source Code

FUNCTION NewtonMin2D(tol,a,b,x,fobj) RESULT (fmin)
  ! MODULES
  IMPLICIT NONE
  ! INPUT VARIABLES -------------------------!
  REAL(wp),INTENT(IN)    :: tol        !! abort tolerance
  REAL(wp),INTENT(IN)    :: a(2),b(2)  !! search interval (2D)
  ! OUTPUT VARIABLES -------------------------!
  REAL(wp),INTENT(INOUT) :: x(2) !! initial guess on input, result on output
  CLASS(c_newton_Min2D)  :: fobj !! functional f(x,y) to minimize with f, f', f''
  REAL(wp)               :: fmin !! on output =f(x,y)
  ! LOCAL VARIABLES --------------------------!
  INTEGER             :: iter,maxiter
  REAL(wp)            :: dx(2)
  REAL(wp)            :: det_Hess
  REAL(wp)            :: gradF(2),Hess(2,2),HessInv(2,2)
  LOGICAL             :: converged
  ! CODE --------------------------------------------------------------------------------------------------------------------------!
  converged=.FALSE.
  maxiter=50
  DO iter=1,maxiter
    Hess=fobj%ddFR(x)
    det_Hess = Hess(1,1)*Hess(2,2)-Hess(1,2)*Hess(2,1)
    IF(det_Hess.LT.1.0E-12) CALL abort(__STAMP__,&
                                      'det Hessian=0 in NewtonMin')
    HessInv(1,1)= Hess(2,2)
    HessInv(1,2)=-Hess(1,2)
    HessInv(2,1)=-Hess(2,1)
    HessInv(2,2)= Hess(1,1)
    HessInv=HessInv/det_Hess
    gradF=fobj%dFR(x)
    dx=-MATMUL(HessInv,gradF)
    dx = MAX(-(x-a),MIN(b-x,dx)) !respect bounds
    x = x+dx
    converged=(SQRT(SUM(dx*dx)).LT.tol).AND.ALL(x.GT.a).AND.ALL(x.LT.b)
    IF(converged) EXIT
  END DO !iter
  IF(.NOT.converged) CALL abort(__STAMP__,&
                                'NewtonMin2D not converged')
  fmin=fobj%FR(x)

END FUNCTION NewtonMin2D