! *****************************************************************
! *****************************************************************

subroutine getn(n)

  implicit none

  ! SCALAR ARGUMENTS
  integer, intent(out) :: n

  n = 2

end subroutine getn

! *****************************************************************
! *****************************************************************

subroutine setx0(n,x)

  implicit none

  ! SCALAR ARGUMENTS
  integer, intent(in) :: n

  ! ARRAY ARGUMENTS
  real(kind=8), intent(out) :: x(n)

  x(1) = 0.0d0
  x(2) = 0.0d0

end subroutine setx0

! *****************************************************************
! *****************************************************************

subroutine evalf(n,x,f,flag)

  implicit none

  ! SCALAR ARGUMENTS
  integer, intent(in) :: n
  integer, intent(out) :: flag
  real(kind=8), intent(out) :: f

  ! ARRAY ARGUMENTS
  real(kind=8), intent(in) :: x(n)

  flag = 0

  f = x(1) * x(2) + 0.1d0 * ( x(1) - x(2) ) ** 4 + ( x(1) + x(2) ) ** 4

end subroutine evalf

! *****************************************************************
! *****************************************************************

subroutine evalg(n,x,g,flag)

  implicit none

  ! SCALAR ARGUMENTS
  integer, intent(in) :: n
  integer, intent(out) :: flag

  ! ARRAY ARGUMENTS
  real(kind=8), intent(in) :: x(n)
  real(kind=8), intent(out) :: g(n)

  flag = 0

  g(1) = x(2) + 0.4d0 * ( x(1) - x(2) ) ** 3 + 4.0d0 * ( x(1) + x(2) ) ** 3
  g(2) = x(1) - 0.4d0 * ( x(1) - x(2) ) ** 3 + 4.0d0 * ( x(1) + x(2) ) ** 3

end subroutine evalg

! *****************************************************************
! *****************************************************************

subroutine evalh(n,x,H,flag)

  implicit none

  ! SCALAR ARGUMENTS
  integer, intent(in) :: n
  integer, intent(out) :: flag

  ! ARRAY ARGUMENTS
  real(kind=8), intent(in) :: x(n)
  real(kind=8), intent(out) :: H(n,n)

  ! Upper triangle only

  flag = 0

  H(1,1) =         1.2d0 * ( x(1) - x(2) ) ** 2 + 12.0d0 * ( x(1) + x(2) ) ** 2
  H(1,2) = 1.0d0 - 1.2d0 * ( x(1) - x(2) ) ** 2 + 12.0d0 * ( x(1) + x(2) ) ** 2
  H(2,2) =         1.2d0 * ( x(1) - x(2) ) ** 2 + 12.0d0 * ( x(1) + x(2) ) ** 2

end subroutine evalh
