python數(shù)學(xué)建模 (1)規(guī)劃問(wèn)題

線性規(guī)劃


x自變量,c是題目中給出的求解的參數(shù)(有多少個(gè)自變量就有多少參數(shù)),A不等式約束條件,c^T就是c的轉(zhuǎn)置。
求解前應(yīng)轉(zhuǎn)化為標(biāo)準(zhǔn)形式

scipy庫(kù)求解

from scipy import optimize
import numpy as np
#求解函數(shù)
res = optimize.linprog(c, A, b, Aeq, beq, LB, UB, X0, OPTIONS)
#目標(biāo)函數(shù)最小值
print(res.fun)
#最優(yōu)解
print(res.x)

樣例

樣例1
#導(dǎo)入包
from scipy import optimize
import numpy as np
#確定c, A, b, Aeq, beq
c = np.array([2, 3, -5])
A = np.array([[-2, 5, -1], [1,3,1]]) #方向是反的,取負(fù)
B = np.array([-10, 12])  
Aeq = np.array([[1,1,1]])  # A和Aeq都是二維數(shù)組
Beq = np.array([7])
#求解
res = optimize.linprog(-c, A, B, Aeq, Beq) #求最大值,加個(gè)負(fù)號(hào)
print(res)

注意A和Aeq里面有兩個(gè)括號(hào)(矩陣運(yùn)算的規(guī)則,兩個(gè)括號(hào)表示是在一行的)

結(jié)果


只需關(guān)注funx,fun是函數(shù)值,x是最優(yōu)解

pulp庫(kù)求解

樣例2
import pulp
#目標(biāo)函數(shù)的系數(shù)
z = [2, 3, 1]
#約束
a = [[1, 4, 2], [3, 2, 0]]
b = [8, 6]
aeq = [[1,2,4]]
beq = [101]
#確定最大化最小化問(wèn)題,最大化只要把Min改成Max即可
m = pulp.LpProblem(sense = pulp.LpMaximize)
#定義三個(gè)變量放到列表中
x = [pulp.LpVariable(f'x{i}', lowBound = 0) for i in [1, 2, 3]]
#定義目標(biāo)函數(shù), lpDot可以將兩個(gè)列表的對(duì)應(yīng)位相乘再加和
#相當(dāng)于z[0] * x[0] + z[1] * x[1] + z[2] * x[2]
m += pulp.lpDot(z, x)   
#設(shè)置約束條件
for i in range(len(a)):
  m += (pulp.lpDot(a[i], x) >= b[i] )
for i in range(len(aeq)):
  m += (pulp.lpDot(a[i], x) == beq[i])
#求解
m.solve()
#輸出結(jié)果
print(f'優(yōu)化結(jié)果:{pulp.value(m.objective)} ')
print(f'參數(shù)取值: {[pulp.value(var) for var in x]}')
結(jié)果

運(yùn)輸問(wèn)題

樣例

import pulp
import numpy as np
from pprint import pprint
def transportation_problem(costs, x_max, y_max):
  row = len(costs)
  col = len(costs[0])
  prob = pulp.LpProblem('Transportation Problem', sense = pulp.LpMaximize)
  var = [[pulp.LpVariable(f'x{i}{j}', lowBound = 0, cat = pulp.LpInteger) for j in range(col)] for i in range(row)]
  flatten = lambda x : [y for l in x for y in flatten(l)] if type(x) is list else [x]  
  prob += pulp.lpDot(flatten(var), costs.flatten())
  for i in range(row):
    prob +=  (pulp.lpSum(var[i]) <= x_max[i])
  for j in range(col):
    prob += (pulp.lpSum([var[i][j] for i in range(row)]) <= y_max[j])
  prob.solve()
  return {'objective' : pulp.value(prob.objective), 'var': [[pulp.value(var[i][j]) for j in range(col)] for i in range(row)]}

if __name__ == '__main__':
  costs = np.array([[500, 550, 630, 1000, 800, 700], [800,700,600,950,900,930], [1000,960,840,650,600,700], [1200,1040,980,860,880,780]])
  max_plant = [76,88,96,40]
  max_cultivation = [42,56,44,39,60,59]
  res = transportation_problem(costs, max_plant, max_cultivation)
  print(f'最大值為{res["objective"]}')
  print('各變量的取值為:')
  pprint(res['var'])

整數(shù)規(guī)劃

分支定界代碼

import math
from scipy.optimize import linprog
import sys
def integerPro(c, A, b, Aeq, beq, t=1.0E-12):
  res = linprog(c, A_ub = A, b_ub = b, A_eq = Aeq, b_eq = beq)
  if (type(res.x) is float):
    bestX = [sys.maxsize]*len(c)
  else:
    bestX = res.x
  bestVal = sum([x*y for x,y in zip(c, bestX)])
  if all(((x-math.floor(x)) < t or (math.ceil(x)-x) < t) for x in bestX):
    return (bestVal, bestX) 
  else:
    ind = [i for i, x in enumerate(bestX) if (x-math.floor(x)) > t and (math.ceil(x) - x) > t][0]
    newCon1 = [0]*len(A[0])
    newCon2 = [0]*len(A[0])
    newCon1[ind] = -1
    newCon2[ind] = 1
    newA1 = A.copy()
    newA2 = A.copy()
    newA1.append(newCon1)
    newA2.append(newCon2)
    newB1 = b.copy()
    newB2 = b.copy()
    newB1.append(-math.ceil(bestX[ind]))    
    newB2.append(math.floor(bestX[ind]))
    r1 = integerPro(c, newA1, newB1, Aeq, beq)
    r2 = integerPro(c, newA2, newB2, Aeq, beq)
    if r1[0] < r2[0]:
      return r1
    else:
      return r2    

非線性規(guī)劃

樣例

1.計(jì)算\frac{1}{x} + x的最小值

from scipy.optimize import minimize
import numpy as np
#計(jì)算1/x+x的最小值
def fun(args):
  a = args
  v = lambda x : a / x[0] + x[0]
  return v

if __name__ == "__main__":
  args = (1)    #a
  x0 = np.asarray((2))  #初始猜測(cè)值
  res = minimize(fun(args), x0, method = 'SLSQP')
  print(res.fun)
  print(res.success)
  print(res.x)   

2.計(jì)算\frac {2+x_1}{1+x_2} - 3x_1 + 4x_3的最小值,其中x_1、x_2、x_3范圍在0.10.9之間

from scipy.optimize import minimize
import numpy as np
def fun(args):
  a,b,c,d = args
  v = lambda x : (a + x[0]) / (b + x[1]) - c*x[0] + d*x[2]
  return v

def con(args):
  #約束條件分為eq和ineq
  #eq表示函數(shù)結(jié)果等于0;ineq表示表達(dá)式大于等于0
  x1min, x1max, x2min, x2max, x3min, x3max = args
  cons = ({'type': 'ineq', 'fun': lambda x : x[0] - x1min}, 
          {'type': 'ineq', 'fun': lambda x : -x[0] + x1max}, 
          {'type': 'ineq', 'fun': lambda x : x[1] - x2min}, 
          {'type': 'ineq', 'fun': lambda x : -x[1] + x2max},  
          {'type': 'ineq', 'fun': lambda x : x[2] - x3min}, 
          {'type': 'ineq', 'fun': lambda x : -x[2] + x3max})
  return cons

if __name__ == "__main__":
  #定義常量值
  args = (2, 1, 3, 4)  #abcd
  #設(shè)置參數(shù)范圍/約束條件
  argsl = (0.1, 0.9, 0.1, 0.9, 0.1, 0.9) #x1min, x1max, x2min, x2max……
  cons = con(argsl)
  #設(shè)置初始猜測(cè)值
  x0 = np.asarray((0.5, 0.5, 0.5))
  res = minimize(fun(args), x0, method = 'SLSQP', constraints = cons)
  print(res.fun)
  print(res.success)
  print(res.x)
?著作權(quán)歸作者所有,轉(zhuǎn)載或內(nèi)容合作請(qǐng)聯(lián)系作者
【社區(qū)內(nèi)容提示】社區(qū)部分內(nèi)容疑似由AI輔助生成,瀏覽時(shí)請(qǐng)結(jié)合常識(shí)與多方信息審慎甄別。
平臺(tái)聲明:文章內(nèi)容(如有圖片或視頻亦包括在內(nèi))由作者上傳并發(fā)布,文章內(nèi)容僅代表作者本人觀點(diǎn),簡(jiǎn)書系信息發(fā)布平臺(tái),僅提供信息存儲(chǔ)服務(wù)。

相關(guān)閱讀更多精彩內(nèi)容

友情鏈接更多精彩內(nèi)容