2016-01-12 13 views
18

私は、人が使う最適なN単位を見つけることによって、効用関数を最大化しようとしています。制約の1つは、彼らが有限の金額を持っていることです。m。だから私は長さ3の配列N倍の価格、Pも長さ3の配列は、mより大きいことはできません制約を設定しようとしている。2つの引数を使用するScipy Optimizer制約

例と同様に以下のように、この最適化のため

P = np.array([3,4,5]) 
N = np.array([1,2,1]) 
m = 50 
sum(P*N) > m 

を、Pが以前の最適化に基づいて与えられます。

cons_c = [{'type':'ineq', 'fun': lambda N: 10 - sum(np.round(N)*P)},{'type':'ineq', 'fun': lambda N: 24 - sum(N*T)}] 
bnds = [(0.,None) for x in range(len(N))] 
optimized_c = scipy.optimize.minimize(utility_c, N, (P,Q,T), method='SLSQP', bounds=bnds, constraints=cons_c) 

機能:

def utility_c(N,P,Q,T): 

    print "N: {0}".format(N) 
    print "P: {0}".format(P) 
    print "Q: {0}".format(Q) 
    print "T: {0}".format(T) 
    N = np.round(N) 
    m = 10 - sum(N*P) 
    b = sum(N*Q) 
    t = 24 - sum(N*T) 
    print "m in C: {0}".format(m) 
    print "b: {0}".format(b) 
    print "t: {0}".format(t) 
    # if m < 0 or t < 0: 
    #  return 0 
    return 1/ ((b**0.3)*(t**0.7))+(5*(m**0.5)) 

問題は、私はまだ負mを取得している今、ここに私の実際のコードです!だから、私は制約を適切に渡していません。 Pが正しく使用されていないためだと思いますか?

出力:

N: [ 1. 1. 1.] 
P: [ 5. 14. 4.] 
Q: [ 1. 3. 1.] 
T: [ 1. 1. 1.01] 
m in C: -13.0 

私が試したもの:

を私もそうのように、引数でPを渡して試してみた:

cons_c = [{'type':'ineq', 'fun': lambda N,P: 10 - sum(np.round(N)*P), 'args':P},{'type':'ineq', 'fun': lambda N: 24 - sum(N*T)}] 

しかし、それは告げます私はラムダが2つの議論を望んで4を受け取った。

**更新:'args'(F,)を使用して**

は、プログラムがエラーを出さずに実行することはできません、しかし制約がまだホールドアップに失敗しました。

また、mが負の値として定義された後にnanが返されます。これはもちろん、scipyの最適化を完全に無効にします。

**全プロジェクトコード:**

import scipy.optimize 
import numpy as np 
import sys 

def solve_utility(P,Q,T): 
    """ 
    Here we are given the pricing already (P,Q,T), but solve for the quantities each type 
    would purchase in order to maximize their utility (N). 
    """ 

    def utility_a(N,P,Q,T): 
     N = np.round(N) 
     m = 50 - sum(N*P) 
     b = sum(N*Q) 
     t = 8 - sum(N*T) 
     return 1/ ((b**0.5)*(t**0.5))+(5*(m**0.5)) 

    def utility_b(N,P,Q,T): 
     N = np.round(N) 
     m = 50 - sum(N*P) 
     b = sum(N*Q) 
     t = 8 - sum(N*T) 
     return 1/ ((b**0.7)*(t**0.3))+(5*(m**0.5)) 

    def utility_c(N,P,Q,T): 
     N = np.round(N) 
     print "N: {0}".format(N) 
     print "P: {0}".format(P) 
     print "Q: {0}".format(Q) 
     print "T: {0}".format(T) 
     m = 10 - sum(N*P) 
     b = sum(N*Q) 
     t = 24 - sum(N*T) 
     print "m in C: {0}".format(m) 
     print "b: {0}".format(b) 
     print "t: {0}".format(t) 
     return 1/ ((b**0.3)*(t**0.7))+(5*(m**0.5)) 



    # Establishing constraints so no negative money or time: 
    N = np.array([2,2,1]) 
    cons_a = [{'type':'ineq', 'fun': lambda N, P: 50 - sum(np.round(N)*P), 'args':(P,)},{'type':'ineq', 'fun': lambda N: 8 - sum(N*T)}] 
    cons_b = [{'type':'ineq', 'fun': lambda N, P: 50 - sum(np.round(N)*P), 'args':(P,)},{'type':'ineq', 'fun': lambda N: 8 - sum(N*T)}] 
    cons_c = [{'type':'ineq', 'fun': lambda N, P: 10 - sum(np.round(N)*P), 'args':(P,)},{'type':'ineq', 'fun': lambda N: 24 - sum(N*T)}] 

    maxes = P/50 
    bnds = [(0.,None) for x in range(len(N))] 
    b = [()] 

    optimized_a = scipy.optimize.minimize(utility_a, N, (P,Q,T), method='SLSQP', constraints=cons_a) 
    optimized_b = scipy.optimize.minimize(utility_b, N, (P,Q,T), method='SLSQP', constraints=cons_b) 
    optimized_c = scipy.optimize.minimize(utility_c, N, (P,Q,T), method='SLSQP', constraints=cons_c) 

    if not optimized_a.success: 
     print "Solving Utilities A didn't work..." 
     return None 
    if not optimized_b.success: 
     print "Solving Utilities B didn't work..." 
     return None 
    if not optimized_c.success: 
     print "Solving Utilities C didn't work..." 
     return None 
    else: 
     print "returning N: {0}".format(np.array([optimized_a.x,optimized_b.x,optimized_c.x])) 
     return np.array([optimized_a.x,optimized_b.x,optimized_c.x]) 


# solve_utility(P,Q,T,N) 


def solve_profits(): 
    """ 
    Here we build the best pricing strategy to maximize solve_profits 
    """ 

    P = np.array([ 3,  10.67,  2.30]) # Pricing 
    Q = np.array([ 1,  4,  1]) # Quantity of beer for each unit 
    T = np.array([ 1,  1, 4]) # Time cost per unit 
    N = np.array([ 1,  0,  1]) # Quantities of unit taken by customer 

    def profit(X): 
     P,Q,T = X[0:3], X[3:6], X[6:9] 
     Q[1] = round(Q[1]) # needs to be an integer 
     N = solve_utility(P,Q,T) 
     print "N: {0}".format(N) 
     N = np.sum(N,axis=1) 
     # print "P: {0}".format(P) 
     # print "Q: {0}".format(Q) 
     # print "T: {0}".format(T) 
     denom = sum(N*P*Q) - sum(Q*N) 
     return 1/ (sum(N*P*Q) - sum(Q*N)) 

    cons = [{'type':'ineq', 'fun': lambda X: X[8] - X[6] - 0.01 }, # The time expense for a coupon must be 0.01 greater than regular 
      {'type':'ineq', 'fun': lambda X: X[4] - 2 }, # Packs must contain at least 2 beers 
      {'type':'eq', 'fun': lambda X: X[3] - 1}, # Quantity has to be 1 for single beer 
      {'type':'eq', 'fun': lambda X: X[5] - 1}, # same with coupons 
      {'type':'ineq', 'fun': lambda X: X[6] - 1}, # Time cost must be at least 1 
      {'type':'ineq', 'fun': lambda X: X[7] - 1}, 
      {'type':'ineq', 'fun': lambda X: X[8] - 1}, 
      ] 
    X = np.concatenate([P,Q,T]) 
    optimized = scipy.optimize.minimize(profit, X, method='L-BFGS-B', constraints=cons) 

    if not optimized.success: 
     print "Solving Profits didn't work..." 
    else: 
     return optimized.x, N 

X, N = solve_profits() 
print "X: {0} N {1}".format(X,N) 
P,Q,T = X[0:3], X[3:6], X[6:9] 
rev = sum(P * Q * N) 
cost = sum(Q * N) 
profit = (rev-cost)*50 
print "N: {0}".format(N) 
print "P: {0}".format(P) 
print "Q: {0}".format(Q) 
print "T: {0}".format(T) 
print "profit = {0}".format(profit) 
+1

をフォローする素敵な導関数を持っていない場合

scipyのダウンロード最適化は素晴らしいではありません、あなたの制約でラウンドを使用する理由はありますか? – Moritz

+0

そして、私はそれがPの3つの要素でリストを見ているので文句を言うと思う。 'args =(P、)' – Moritz

+0

@Moritzはまだ残念ながら働かなかった – thefoxrocks

答えて

3

あなたはoptimized_aおよび実行のために、あなたはそれがスローエラーが表示隔離した場合はエラー8である - これは正の派生エラーです。

BFGSとSLSQPはどちらも、最初の推測を取り、グラジエントとその派生関数を評価し、最適な方向を見つけ、常に下り坂に降りて停止し、値が設定した許容値を下回っているか、最小値に達しています。

エラーは、少なくともあなたの最初の推測では、問題に強い派生がないことを示唆しています。一般に、SQLSPは、正方形の合計として定式化できる問題で最もよく使用されます。おそらくより現実的な初期推測を試みることが役に立ちます。私は間違いなくほとんどのコードを破棄して、まずoptimize_aを使って最小限の例を実行します。

おそらく非勾配ベースのソルバーが動作するか、問題のサイズとパラメータの現実的な境界に応じて、グローバルな最適化が実行可能になる場合があります。あなたは

関連する問題