# Path: pybobyqa/hessian.py
# def to_upper_triangular_vector(A):
#     n = A.shape[0]
#     hq = np.zeros((n*(n+1)//2,))
#     ih = -1
#     for j in range(n):  # j = 0, ..., n-1
#         for i in range(j + 1):  # i = 0, ..., j
#             ih += 1
#             hq[ih] = A[i,j]
#     return hq
# 
# Path: pybobyqa/trust_region.py
# def trsbox_geometry(xbase, c, g, H, lower, upper, Delta, use_fortran=USE_FORTRAN):
#     # Given a Lagrange polynomial defined by: L(x) = c + g' * (x - xbase) + 0.5*(x-xbase)*H*(x-xbase)
#     # Maximise |L(x)| in a box + trust region - that is, solve:
#     #   max_x  abs(c + g' * (x - xbase) + 0.5*(x-xbase)*H*(x-xbase))
#     #    s.t.  lower <= x <= upper
#     #          ||x-xbase|| <= Delta
#     # Setting s = x-xbase (or x = xbase + s), this is equivalent to:
#     #   max_s  abs(c + g' * s + 0.5*s*H*s)
#     #   s.t.   lower <= xbase + s <= upper
#     #          ||s|| <= Delta
#     smin, gmin, crvmin = trsbox(xbase, g, H, lower, upper, Delta, use_fortran=use_fortran)  # minimise L(x)
#     smax, gmax, crvmax = trsbox(xbase, -g, -H, lower, upper, Delta, use_fortran=use_fortran)  # maximise L(x)
#     if abs(c + model_value(g, H, smin)) >= abs(c + model_value(g, H, smax)):  # take largest abs value
#         return xbase + smin
#     else:
#         return xbase + smax
# 
# Path: pybobyqa/util.py
# def sumsq(x):
#     # There are several ways to calculate sum of squares of a vector:
#     #   np.dot(x,x)
#     #   np.sum(x**2)
#     #   np.sum(np.square(x))
#     #   etc.
#     # Using the timeit routine, it seems like dot(x,x) is ~3-4x faster than other methods
#     return np.dot(x, x)
# 
# def model_value(g, H, s):
#     # Calculate model value (s^T * g + 0.5* s^T * H * s) = s^T * (gopt + 0.5 * H*s)
#     assert g.shape == s.shape, "g and s have incompatible sizes"
#     Hs = H.dot(s)
#     return np.dot(s, g + 0.5*Hs)


import logging
import numpy as np
import scipy.linalg as LA
from math import sqrt
from .hessian import to_upper_triangular_vector
from .trust_region import trsbox_geometry
from .util import sumsq, model_value
            # Apply bounds and convert back to absolute coordinates
            return self.xbase + np.minimum(np.maximum(self.sl, self.points[k, :]), self.su)

    def fval(self, k):
        assert 0 <= k < self.npt(), "Invalid index %g" % k
        return self.f_values[k]

    def as_absolute_coordinates(self, x):
        # If x were an interpolation point, get the absolute coordinates of x
        return self.xbase + np.minimum(np.maximum(self.sl, x), self.su)

    def xpt_directions(self, include_kopt=True):
        if include_kopt:
            ndirs = self.npt()
        else:
            ndirs = self.npt() - 1

        dirns = np.zeros((ndirs, self.n()))  # vector of directions xpt - xopt, excluding for xopt
        xopt = self.xopt()
        for k in range(self.npt()):
            if not include_kopt and k == self.kopt:
                continue  # skipt
            idx = k if include_kopt or k < self.kopt else k - 1
            dirns[idx, :] = self.xpt(k) - xopt
        return dirns

    def distances_to_xopt(self):
        sq_distances = np.zeros((self.npt(),))
        xopt = self.xopt()
        for k in range(self.npt()):

# Path: pybobyqa/hessian.py
# def to_upper_triangular_vector(A):
#     n = A.shape[0]
#     hq = np.zeros((n*(n+1)//2,))
#     ih = -1
#     for j in range(n):  # j = 0, ..., n-1
#         for i in range(j + 1):  # i = 0, ..., j
#             ih += 1
#             hq[ih] = A[i,j]
#     return hq
# 
# Path: pybobyqa/trust_region.py
# def trsbox_geometry(xbase, c, g, H, lower, upper, Delta, use_fortran=USE_FORTRAN):
#     # Given a Lagrange polynomial defined by: L(x) = c + g' * (x - xbase) + 0.5*(x-xbase)*H*(x-xbase)
#     # Maximise |L(x)| in a box + trust region - that is, solve:
#     #   max_x  abs(c + g' * (x - xbase) + 0.5*(x-xbase)*H*(x-xbase))
#     #    s.t.  lower <= x <= upper
#     #          ||x-xbase|| <= Delta
#     # Setting s = x-xbase (or x = xbase + s), this is equivalent to:
#     #   max_s  abs(c + g' * s + 0.5*s*H*s)
#     #   s.t.   lower <= xbase + s <= upper
#     #          ||s|| <= Delta
#     smin, gmin, crvmin = trsbox(xbase, g, H, lower, upper, Delta, use_fortran=use_fortran)  # minimise L(x)
#     smax, gmax, crvmax = trsbox(xbase, -g, -H, lower, upper, Delta, use_fortran=use_fortran)  # maximise L(x)
#     if abs(c + model_value(g, H, smin)) >= abs(c + model_value(g, H, smax)):  # take largest abs value
#         return xbase + smin
#     else:
#         return xbase + smax
# 
# Path: pybobyqa/util.py
# def sumsq(x):
#     # There are several ways to calculate sum of squares of a vector:
#     #   np.dot(x,x)
#     #   np.sum(x**2)
#     #   np.sum(np.square(x))
#     #   etc.
#     # Using the timeit routine, it seems like dot(x,x) is ~3-4x faster than other methods
#     return np.dot(x, x)
# 
# def model_value(g, H, s):
#     # Calculate model value (s^T * g + 0.5* s^T * H * s) = s^T * (gopt + 0.5 * H*s)
#     assert g.shape == s.shape, "g and s have incompatible sizes"
#     Hs = H.dot(s)
#     return np.dot(s, g + 0.5*Hs)


# Path: pybobyqa/model.py
import logging
import numpy as np
import scipy.linalg as LA
from math import sqrt
from .hessian import to_upper_triangular_vector
from .trust_region import trsbox_geometry
from .util import sumsq, model_value

            # Apply bounds and convert back to absolute coordinates
            return self.xbase + np.minimum(np.maximum(self.sl, self.points[k, :]), self.su)

    def fval(self, k):
        assert 0 <= k < self.npt(), "Invalid index %g" % k
        return self.f_values[k]

    def as_absolute_coordinates(self, x):
        # If x were an interpolation point, get the absolute coordinates of x
        return self.xbase + np.minimum(np.maximum(self.sl, x), self.su)

    def xpt_directions(self, include_kopt=True):
        if include_kopt:
            ndirs = self.npt()
        else:
            ndirs = self.npt() - 1

        dirns = np.zeros((ndirs, self.n()))  # vector of directions xpt - xopt, excluding for xopt
        xopt = self.xopt()
        for k in range(self.npt()):
            if not include_kopt and k == self.kopt:
                continue  # skipt
            idx = k if include_kopt or k < self.kopt else k - 1
            dirns[idx, :] = self.xpt(k) - xopt
        return dirns

    def distances_to_xopt(self):
        sq_distances = np.zeros((self.npt(),))
        xopt = self.xopt()
        for k in range(self.npt()):

            sq_distances[k] = sumsq(self.points[k, :] - xopt)