Why constrained minimization fails in python

I'm trying to numerically solve some mechanism design question using the following code. However, the output is not optimal(just stay at my initial guess). How should I make it work?

I tried different initial guesses ( I can calculate myself so the initial guess is valid, each x0 may give different solutions).

I also switched methods, but some results (like Nelder-Mead) generate really invalid solutions (violating my constraints).

import numpy as np
from scipy.optimize import minimize

p0=0.3
r1=0.4
r2=0.8


p1s=r1/(r1+1)
p2s=r2/(r2+1)

q1=(r1+r1*r2)/(r2+r1*r2)

q2=r1/r2

q=(q1+q2)/2

def v_1(k,p0):
    if p0>p1s and k>p0:
        return r1-p0*r1/k
    elif p0<p1s and k>p1s:
        return p0+p0*r1-p0*r1/k
    return 0

def v_2(k,p0):
    if p0>p2s and k>p0:
        return r2-p0*r2/k
    elif p0<p2s and k>p2s:
        return p0+p0*r2-p0*r2/k
    return 0

def V(p):
    if p<p1s:
        return p
    if p>p2s:
        return q*r2*(1-p)
    A=r1-r1/(r2-r1)*(p+p*r2-r1)
    return q*p+(1-q)*A

def objective(x):
    k1,k2,p1,p2 = x
    return -q*(p2+V(k2))+(1-q)*(p1+V(k1))  # Convert to minimization

def IR1(x):
    k1,k2,p1,p2 = x
    return v_1(k1,p0)-p1

def IR2(x):
    k1,k2,p1,p2 = x
    return v_2(k2,p0)-p2

def IC1(x):
    k1,k2,p1,p2 = x
    return v_1(k1,p0)-p1-v_1(k2,p0)-p2

def IC2(x):
    k1,k2,p1,p2 = x
    return v_2(k2,p0)-p2-v_2(k1,p0)-p1

constraints = [
    {'type': 'ineq', 'fun': lambda x: IR1(x)},  
    {'type': 'ineq', 'fun': lambda x: IR2(x)},  
    {'type': 'ineq', 'fun': lambda x: IC1(x)},  
    {'type': 'ineq', 'fun': lambda x: IC2(x)}  
]

# Bounds for x1, x2
bounds = [(p0, 1), (p0, 1), (0,p0), (0,p0)]  #k1,k2>=p0,p0>=p1,p2>=0

# Initial guess
x0 = [p0, 1, 0,p0 ]

# Solve the problem
result = minimize(objective, x0, bounds=bounds, constraints=constraints,method='SLSQP')

# Print results
if result.success:
    print("Optimal solution found:")
    print(f"k1 = {result.x[0]:.4f}, k2 = {result.x[1]:.4f}")
    print(f"p1 = {result.x[2]:.4f}, p2 = {result.x[3]:.4f}")
    print(f"Maximum function value: {-result.fun:.4f}")
else:
    print("Optimization failed:", result.message)**strong text**

Upvotes: 0

Views: 22

Answers (0)

Related Questions