Python array optimization with two constraints

Question

I have an optimization problem where I'm trying to find an array that needs to optimize two functions simultaneously.

In the minimal example below I have two known arrays w and x and an unknown array y. I initialize array y to contains only 1s.

I then specify function np.sqrt(np.sum((x-np.array)**2) and want to find the array y where

np.sqrt(np.sum((x-y)**2) approaches 5

np.sqrt(np.sum((w-y)**2) approaches 8

The code below can be used to successfully optimize y with respect to a single array, but I would like to find that the solution that optimizes y with respect to both x and y simultaneously, but am unsure how to specify the two constraints.

y should only consist of values greater than 0.

Any ideas on how to go about this ?

w = np.array([6, 3, 1, 0, 2])
x = np.array([3, 4, 5, 6, 7])
y = np.array([1, 1, 1, 1, 1])

def func(x, y):

    z = np.sqrt(np.sum((x-y)**2)) - 5
    return  np.zeros(x.shape[0],) + z

r = opt.root(func, x0=y, method='hybr')
print(r.x)
# array([1.97522498 3.47287981 5.1943792  2.10120135 4.09593969])

print(np.sqrt(np.sum((x-r.x)**2)))
# 5.0

Answer 1

One option is to use scipy.optimize.minimize instead of root, Here you have multiple solver options and some of them (ie SLSQP) allow you to specify multiple constraints. Note that I changed the variable names so that x is the array you want to optimise and y and z define the constraints.

from scipy.optimize import minimize
import numpy as np

x0 = np.array([1, 1, 1, 1, 1])
y = np.array([6, 3, 1, 0, 2])
z = np.array([3, 4, 5, 6, 7])

constraint_x = dict(type='ineq',
                    fun=lambda x: x)   # fulfilled if > 0
constraint_y = dict(type='eq',
                    fun=lambda x: np.linalg.norm(x-y) - 5)  # fulfilled if == 0
constraint_z = dict(type='eq',
                    fun=lambda x: np.linalg.norm(x-z) - 8)  # fulfilled if == 0

res = minimize(fun=lambda x: np.linalg.norm(x), constraints=[constraint_y, constraint_z], x0=x0,
               method='SLSQP', options=dict(ftol=1e-8))  # default 1e-6

print(res.x)                    # [1.55517124 1.44981672 1.46921122 1.61335466 2.13174483]
print(np.linalg.norm(res.x-y))  # 5.00000000137866
print(np.linalg.norm(res.x-z))  # 8.000000000930026

This is a minimizer so besides the constraints it also wants a function to minimize, I chose just the norm of y, but setting the function to a constant (ie lambda x: 1) would have also worked. Note also that the constraints are not exactly fulfilled, you can increase the accuracy by setting optional argument ftol to a smaller value ie 1e-10. For more information see also the documentation and the corresponding sections for each solver.