Model Infeasibility

I am running a known feasible model, but I'm not sure what I have done incorrectly to make it state that the model is infeasible. Removing the binary constraints makes it feasible, but these constraints should not be the issue.

import numpy as np

from numpy import sin, cos, tan, pi

import matplotlib.pyplot as plt

import gurobipy as gp

from gurobipy import GRB

import scipy

# General Vehicle Parameters

g = 9.8 # Gravity Constant

m = 3600/2.2 # Vehicle Mass

L = 2.1336 # Wheelbase

l_fl_r = 0.85 # A/B Ratio

alfaf_max = 7*pi/180 # Design Limit

alat_max = 0.7*g # Design Limit

K_u = 1.1 # Understeer

b = L/(1+l_fl_r) # CG Distance to Rear Wheel Center

a = L-b # CG Distance to Front Wheel Center

Iz = .5*m*a*b # Vehicle Yaw Inertia

mu = 0.7

# Tire parameters

Cr = 5.9649e+04 # Rear tire stiffness

Cf = 6.1138e+04 # Front tire stiffness

# Constant longitudinal velocity

U = 15 # Longitudinal Velocity [m/s]

# Linear Bycicle model Dynamics

# [V(m/s), r(rad/s), Psi(deg)]

A = np.array([

[-(Cf+Cr)/(m*U),(b*Cr-a*Cf)/(m*U)-U, 0],

[(b*Cr-a*Cf)/(Iz*U),-(b**2*Cr+a**2*Cf)/(Iz*U),0],

[0, 1, 0] ])

#B matrix with steering inputs

B = np.array([

[Cf/m],

[a*Cf/Iz],

[0]])

# obstacle parameters

vertexL = np.array([[-1,11,25],[-4,0,-4]]) #[xmin;ymin]

vertexU = np.array([[1,13,27],[1.5,8,2]]) #[xmax;ymax]

xminO = vertexL[0]

xmaxO = vertexU[0]

yminO = vertexL[1]

ymaxO = vertexU[1]

# M = [45,35,20,25,35] # big-M

M = [52,52,9,9,15]; # big-M

nobs = vertexL.shape[1] # number of obstacles

N = 30

T = 3.5

# initial condition

xinit = [0,0,0]

XYinit = [-15,0]

Z0= np.array([xinit[0],xinit[1],xinit[2],XYinit[0],XYinit[1]])

# Create variables

Z = np.array([])

bv1 = np.array([])

# bv2 = np.array([])

# bv3 = np.array([])

u = np.array([])

# Bounds

umin = -35*pi/180

umax = 35*pi/180

# Create a new model

mod = gp.Model("MPC")

Z = mod.addMVar(shape=(N+1,5), lb=-GRB.INFINITY, ub=GRB.INFINITY, name='Z')

bv1 = mod.addMVar(shape=(N,4), lb=0, ub =1, vtype=GRB.INTEGER, name='bv1')

# bv2 = mod.addMVar(shape=(N,4), lb=0, ub =1, type=GRB.INTEGER, name='bv2')

# bv3 = mod.addMVar(shape=(N,4), lb=0, ub =1, type=GRB.INTEGER, name='bv3')

u = mod.addMVar(N, lb=umin, ub=umax, name='u')

L = mod.addMVar(N+1, lb =0, ub = GRB.INFINITY, name='L')

cosx = mod.addMVar(N, name='cosx')

sinx = mod.addMVar(N, name='sinx')

# MPC Formulation

mod.addConstr(Z[0] == Z0)

for k in range(N):

mod.addConstr( Z[k+1, :3] == Z[k,:3] + T/N*(A @ Z[k,:3] + B *u[k]) )

mod.addConstr( Z[k+1,3] == Z[k,3] + T/N * ( U*cosx[k] - Z[k,0]*sinx[k] ) )

mod.addConstr( Z[k+1,4] == Z[k,4] + T/N * ( U*sinx[k] + Z[k,0]*cosx[k] ) )

mod.addGenConstrCos(Z[k,2], cosx[k])

mod.addGenConstrSin(Z[k,2], sinx[k])

mod.addConstr( Z[k,3] <= xminO[0] + M[0]*bv1[k,0] )

mod.addConstr( -Z[k,3] <= -xmaxO[0] + M[1]*bv1[k,1] )

mod.addConstr( Z[k,4] <= yminO[0] + M[2]*bv1[k,2] )

mod.addConstr( -Z[k,4] <= -ymaxO[0] + M[3]*bv1[k,3] )

mod.addConstr( sum(bv1[k,:]) <= 3 )

# if k!=0:

# m.addConstr( Z[k-1,3] <= xminO[0] + M[0]*bv1[k,0] )

# m.addConstr( -Z[k-1,3] <= -xmaxO[0] + M[1]*bv1[k,1] )

# m.addConstr( Z[k-1,4] <= yminO[0] + M[2]*bv1[k,2] )

# m.addConstr( -Z[k-1,4] <= -ymaxO[0] + M[3]*bv1[k,3])

mod.addConstr( Z[k,4] <= L[k] )

mod.addConstr( -Z[k,4] <= L[k] )

mod.addConstr( Z[N,4] <= L[N] )

mod.addConstr( -Z[N,4] <= L[N] )

obj = 3*sum(L)

mod.setObjective(obj, GRB.MINIMIZE)

mod.params.FuncNonlinear=1

mod.write("mod_minlp.lp")

# mod_presolve = mod.presolve()

# mod_presolve.write("mod_minlp_presolve.lp")

mod.optimize()

answer from BONMIN solver with casadi interface for Z[3:,:] :

Model Infeasibility

コメント

お探しのものを見つけられませんでしたか？