代码之家 › 专栏 › 技术社区 › Jeff Tilton

Scipy在df中最小化约束

scipy optimization python

Jeff Tilton · 技术社区 · 6 年前

我正在尝试运行中发现的最小化过程 this publication 。方程式见pdf文件16第6页。

我有一个如下所示的数据框

 df = pd.DataFrame({'h_t':[7.06398,6.29948,5.04570,6.20774,4.80106],           
'p_atm':[101057.772801,101324.416001,101857.702401,101724.380801,101991.024001],
    'q_p':[5.768132,3.825600,2.772215,5.830429,2.619304],
    'q_s':[2.684433,3.403679,2.384275,1.008078,2.387106],
    'tdg_f':[117.678100,110.131579,108.376963,103.669725,113.594771],
    'tdg_tw':[121.052635,119.710907,114.921463,112.156868,115.444900],
    'temp_water':[11.92,19.43,16.87,7.45,11.83]})

我有一个约束条件,下面的函数必须是正的,其中b1和b3是我正在优化的系数。

def q_ge(q_p,q_s,b1,b3):   
    return min(q_p,(b1*q_s+b3))

我在下面写下了我的限制,但我不确定它是否正确。

def constraint_q_ge(x):
    b1,b2,b3=x
    power_flow = df.apply(lambda x:q_ge(x['q_p'],x['q_s'],b1,b3), axis = 1)
    const = power_flow<0
    return -const.sum()

对吗?我在所有行上运行该函数,检查是否有任何行小于0,然后求和。该和的负数应大于或等于0。如果甚至有一个值小于0,则不满足此约束。

编辑:

下面是完整的问题。

from scipy.constants import g as gravity
from sklearn.metrics import mean_squared_error
from math import sqrt
from scipy.optimize import minimize
import warnings

try:
    from numpy import any as _any
except ImportError:
    def _any(arg):
        if arg is True:
            return True
        if arg is False:
            return False
        return any(arg)

def water_density(T=None, T0=None, units=None, a=None,
                  just_return_a=False, warn=True):
    if units is None:
        K = 1
        m = 1
        kg = 1
    else:
        K = units.Kelvin
        m = units.meter
        kg = units.kilogram
    if T is None:
        T = 298.15*K
    m3 = m**3
    if a is None:
        a = (-3.983035*K,  # C
             301.797*K,  # C
             522528.9*K*K,  # C**2
             69.34881*K,  # C
             999.974950*kg/m3)
    if just_return_a:
        return a
    if T0 is None:
        T0 = 273.15*K
    t = T - T0
    if warn and (_any(t < 0*K) or _any(t > 40*K)):
        warnings.warn("Temperature is outside range (0-40 degC)")
    return a[4]*(1-((t + a[0])**2*(t + a[1]))/(a[2]*(t + a[3])))


def celsius_to_kelvin(t_celsius):
    return t_celsius+273.15


def tailwater(h_t, temp_water, p_atm):
    t_water_kelvin = celsius_to_kelvin(temp_water)
    rho = water_density(t_water_kelvin)
    g = gravity
    return (1+(rho*g*h_t)/(2*p_atm))

def tailwater_tdg(q_s,q_p,x, h_t,temp_water,p_atm,tdg_f):
    b1,b2,b3=x
    A = ((q_s+b1*q_s+b3)/(q_s+q_p))
    B = tailwater(h_t, temp_water, p_atm)
    C = ((q_p-b1*q_s-b3)/(q_s+q_p))

    return 100*A*B*b2+tdg_f*C

def q_ge(q_p,q_s,b1,b3):
    return min(q_p,(b1*q_s+b3))


def rmse(y, y_hat):
    return sqrt(mean_squared_error(y,y_hat))

def objective(x):
    y_hat = df.apply(lambda r:tailwater_tdg(q_s=r['q_s'],q_p=r['q_p'],x=x, h_t=r['h_t'],temp_water=r['temp_water'],p_atm=r['p_atm'],tdg_f=r['tdg_f']), axis = 1)
    y = df['tdg_tw']
    return rmse(y, y_hat)

#constraints and bounds for optimization model.  See reference for more information
def constraint_q_ge(x):
    b1,b2,b3=x
    power_flow = df.apply(lambda x:q_ge(x['q_p'],x['q_s'],b1,b3), axis = 1)
    const = power_flow<0
    return -const.sum()

constraints = [{'type':'ineq', 'fun':constraint_q_ge}]          
bounds = [(-500,10000),(.00001,10000),(-500,10000)]
x0=[1,1,1]
sol = minimize(objective, x0, method = 'SLSQP',constraints = constraints, bounds = bounds,options={'disp':True, 'maxiter':100})

0 回复 | 直到 6 年前