scipy 如何通过多个参数最小化一个函数？(第一部分)

我想通过4个参数a，b，c和d来最小化函数dchi(a, b, c, d)。
res = sp.optimize.minimize(dchi, first_guess)这一行好像写错了，你能告诉我到底哪里写错了，怎么纠正吗？

import functools
import scipy as sp

def CH(chi: float, T: float, a: float, b: float, c: float, d: float, H: float, TC: float, C:float) -> float:
    return (a*(T - TC)*chi**(1/c)+ b * chi**(1/c) * chi**(1/d) * H**(1/d)- 1)+C/(T-TC)

exp=[[10.0052, 87.9716], [11.0433, 86.3866], [12.359, 84.8015], [13.7127,   83.2164], [14.9939, 82.0276], [16.3367, 80.8388], [17.6282,   80.0463], [18.9464, 78.8575], [18.9934, 78.8575], [18.9969,   78.8575], [18.9962, 78.8575], [19.4221, 78.4612], [20.8034,   77.6687], [22.1306, 76.4799], [23.4352, 75.6873], [24.7503,   75.291], [26.0295, 74.4985], [27.3692, 73.706], [30.0264,   72.9134], [31.3262, 72.5172], [32.6313, 72.1209], [33.9404,   71.7246], [35.2494, 71.3284], [36.551, 70.9321], [37.8445,   70.5358], [40.5813, 69.7433], [41.7696, 69.347], [43.0624,   69.347], [44.3609, 68.9507], [45.6536, 68.5545], [46.9616,   68.5545], [48.254, 68.1582], [49.5428, 67.7619], [52.1298,   67.3657], [53.4147, 67.3657], [54.7091, 66.9694], [56.0001,   66.9694], [57.2883, 66.5731], [58.5865, 66.5731], [59.8911,   66.1769], [62.4856, 66.1769], [63.7665, 65.7806], [65.0689,   65.7806], [66.3588, 65.3843], [67.6581, 65.3843], [68.9546,   65.3843], [70.2572, 64.9881], [72.8482, 64.9881], [74.1637,   64.9881], [75.4546, 64.5918], [76.7495, 64.5918], [78.0347,   64.5918], [79.3376, 64.1955], [80.6437, 64.1955], [83.5861,   64.1955], [84.8957, 63.7993], [86.2123, 63.7993], [87.5317,   63.7993], [88.8468, 63.7993], [90.1674, 63.7993], [91.4775,   63.403], [94.2769, 63.403], [95.5863, 63.403], [96.9007,   63.403], [98.2013, 63.403], [99.523, 63.0067], [100.839,   63.0067], [102.147, 63.0067], [104.812, 63.0067], [106.131,   63.0067], [107.448, 62.6104], [108.789, 62.6104], [110.099,   62.6104], [111.414, 62.6104], [112.756, 62.6104], [115.404,   62.6104], [116.737, 62.6104], [118.059, 62.6104], [119.376,   62.2142], [120.7, 62.2142], [122.025, 62.2142], [123.367,   62.2142], [126.006, 62.2142], [127.352, 62.2142], [128.674,   62.2142], [129.999, 62.2142], [131.333, 62.2142], [132.652,   61.8179], [133.987, 61.8179], [136.652, 61.8179], [137.987,   61.8179], [139.324, 61.8179], [140.645, 61.8179], [141.989,   61.8179], [143.316, 61.8179], [144.646, 61.8179], [147.323,   61.8179], [148.651, 61.4216], [149.978, 61.4216], [151.303,   61.4216], [152.637, 61.4216], [153.949, 61.4216], [155.268,   61.4216], [157.929, 61.4216], [159.252, 61.4216], [160.581,   61.0254], [161.889, 61.0254], [163.219, 61.0254], [164.546,   61.0254], [165.873, 61.0254], [168.524, 61.0254], [169.852,   60.6291], [171.173, 60.6291], [172.494, 60.6291], [173.83,   60.6291], [175.16, 60.6291], [176.473, 60.6291], [179.122,   60.2328], [180.459, 60.2328], [181.78, 60.2328], [183.106,   59.8366], [184.439, 59.8366], [185.757, 59.8366], [187.079,   59.8366], [189.73, 59.4403], [191.06, 59.4403], [192.386,   59.044], [193.717, 59.044], [195.037, 59.044], [196.372,   58.6478], [197.693, 58.6478], [200.361, 58.2515], [201.689,   57.8552], [203.008, 57.8552], [204.335, 57.459], [205.651,   57.459], [206.996, 57.0627], [208.324, 57.0627], [210.977,   56.6664], [212.292, 56.2702], [213.602, 55.8739], [214.947,   55.8739], [216.259, 55.4776], [217.583, 55.0813], [218.895,   55.0813], [221.56, 54.2888], [222.865, 53.8925], [224.195,   53.4963], [225.505, 53.1], [226.824, 53.1], [228.192,   52.7037], [229.508, 52.3075], [232.191, 51.5149], [233.507,   51.1187], [234.836, 50.7224], [236.161, 50.3261], [237.479,   49.9299], [238.807, 49.5336], [240.13, 48.741], [242.774,   47.9485], [244.096, 47.5522], [245.417, 47.156], [246.74,   46.7597], [248.057, 46.3634], [249.382, 45.5709], [250.457,   46.2508], [251.117, 46.0365], [252.897, 45.4166], [253.112,   45.3033], [254.115, 44.9356], [255.112, 44.5707], [256.107,   44.2179], [257.936, 43.554], [259.077, 43.1208], [260.395,   42.633], [260.555, 42.5688], [261.881, 42.0749], [263.515,   41.4803], [264.389, 41.108], [264.559, 41.047], [264.724,   40.9939], [264.891, 40.9305], [265.058, 40.8548], [265.219,   40.807], [265.384, 40.7401], [265.549, 40.6799], [265.707,   40.6116], [265.872, 40.5381], [266.041, 40.4919], [266.218,   40.4131], [266.379, 40.365], [266.532, 40.2826], [266.705,   40.2386], [266.876, 40.1762], [267.038, 40.1109], [267.201,   40.0532], [267.373, 39.9972], [267.543, 39.9247], [267.71,   39.8555], [267.874, 39.7915], [268.035, 39.7336], [268.194,   39.6811], [268.367, 39.6056], [268.542, 39.5387], [268.705,   39.4762], [268.866, 39.3997], [269.029, 39.3423], [269.201,   39.2861], [269.374, 39.22], [269.539, 39.1538], [269.698,   39.0823], [269.87, 39.0231], [270.036, 38.9709], [270.194,   38.8925], [270.371, 38.8394], [270.536, 38.7638], [270.697,   38.7051], [270.868, 38.6425], [271.038, 38.574], [271.208,   38.5115], [271.367, 38.4455], [271.54, 38.3782], [271.711,   38.3156], [271.868, 38.2524], [272.03, 38.1917], [272.2,   38.1262], [272.361, 38.0706], [272.525, 38.0021], [274.156,   37.3821], [274.229, 37.3188], [274.362, 37.27], [274.524,   37.2048], [274.705, 37.1411], [274.877, 37.0827], [275.029,   37.0104], [275.196, 36.9607], [275.37, 36.8851], [275.526,   36.8172], [275.683, 36.7629], [275.855, 36.6982], [276.025,   36.6346], [276.185, 36.5695], [276.345, 36.5063], [276.506,   36.4347], [276.68, 36.3789], [276.845, 36.315], [277.012,   36.2632], [277.184, 36.1689], [277.35, 36.1295], [277.517,   36.0462], [277.676, 35.9975], [277.839, 35.9304], [278.009,   35.8579], [278.18, 35.8111], [278.341, 35.7453], [278.5,   35.6707], [278.67, 35.6269], [278.836, 35.5542], [279.005,   35.486], [279.176, 35.4266], [279.333, 35.3622], [279.502,   35.2938], [279.671, 35.2304], [279.833, 35.1688], [279.997,   35.1014], [280.169, 35.0437], [280.337, 34.9787], [280.498,   34.9181], [280.67, 34.8453], [280.837, 34.7763], [280.995,   34.7195], [281.164, 34.6492], [281.332, 34.5889], [281.5,   34.523], [281.664, 34.4628], [281.83, 34.3953], [282.004,   34.3438], [282.164, 34.2835], [282.324, 34.214], [282.493,   34.1503], [282.67, 34.087], [282.833, 34.0348], [282.996,   33.9549], [283.164, 33.9051], [284.809, 33.291], [284.864,   33.2491], [284.993, 33.1839], [285.16, 33.111], [285.323,   33.0541], [285.487, 32.9876], [285.655, 32.9227], [285.821,   32.8639], [285.991, 32.8064], [286.164, 32.7342], [286.329,   32.6806], [286.494, 32.606], [286.66, 32.562], [286.824,   32.4928], [286.985, 32.4283], [287.149, 32.3634], [287.321,   32.2873], [287.487, 32.2322], [287.644, 32.1765], [287.817,   32.1208], [287.988, 32.0552], [288.146, 31.9866], [288.306,   31.9345], [288.476, 31.8621], [288.646, 31.8045], [288.815,   31.7536], [288.985, 31.6971], [289.14, 31.6201], [289.307,   31.551], [289.481, 31.4893], [289.645, 31.444], [289.81,   31.3815], [289.975, 31.311], [290.141, 31.252], [290.306,   31.1829], [290.466, 31.1326], [290.636, 31.0702], [290.807,   31.011], [290.967, 30.9446], [291.139, 30.8823], [291.304,   30.8283], [291.47, 30.7563], [291.646, 30.706], [291.809,   30.6516], [291.974, 30.578], [292.141, 30.5155], [292.306,   30.4571], [292.47, 30.3986], [292.634, 30.3283], [292.805,   30.2715], [292.971, 30.212], [293.127, 30.1523], [293.299,   30.0891], [293.475, 30.0263], [293.634, 29.9634], [293.794,   29.9016], [295.41, 29.3469], [295.467, 29.2859], [295.612,   29.2319], [295.787, 29.1886], [295.954, 29.1198], [296.116,   29.0857], [296.277, 29.0178], [296.445, 28.9598], [296.617,   28.9062], [296.78, 28.8415], [296.946, 28.7719], [297.111,   28.7037], [297.277, 28.6442], [297.444, 28.5661], [297.603,   28.5065], [297.772, 28.438], [297.945, 28.3781], [298.11,   28.3183], [298.27, 28.2553], [298.436, 28.1991], [298.608,   28.1362], [298.775, 28.089], [298.934, 28.0215], [299.089,   27.9536], [299.26, 27.8898], [299.431, 27.8379], [299.596,   27.7736], [299.77, 27.7224], [299.939, 27.6707], [300.1,   27.6075], [300.268, 27.5584], [300.436, 27.4834], [300.595,   27.445], [300.761, 27.3805], [300.928, 27.3191], [301.085,   27.2599], [301.248, 27.2049], [301.416, 27.1372], [301.582,   27.0904], [301.742, 27.0381], [301.902, 26.9595], [302.077,   26.9182], [302.254, 26.8639], [302.409, 26.7965], [302.573,   26.751], [302.744, 26.6902], [302.904, 26.6289], [303.075,   26.5635], [303.24, 26.5116], [303.398, 26.4599], [303.571,   26.3952], [303.75, 26.3447], [303.915, 26.2879], [304.07,   26.2416], [304.236, 26.1701], [304.396, 26.1088]]

#print(exp[0][0])

def dchi(a, b, c, d):
    sd=0
    for i in range(len(exp)):
        CH_bound = functools.partial(CH, T=exp[i][0], a=a, b=b, c=c, d=d, H=315, TC=1000, C=2)
        ch1 = sp.optimize.root_scalar(f=CH_bound, x0=0)
        sd=sd+(ch1.root - exp[i][1]) ** 2
    return sd

#print(dchi(1,1,1,1))

first_guess = [1,1,1,1]
res = sp.optimize.minimize(dchi, first_guess)

#print(res)

错误：

Traceback (most recent call last):
  File "C:\Users\...\PycharmProjects\pythonProject6\main.py", line 23, in <module>
    res = sp.optimize.minimize(dchi, first_guess)
  File "C:\Users\...\PycharmProjects\pythonProject6\venv\lib\site-packages\scipy\optimize\_minimize.py", line 705, in minimize
    res = _minimize_bfgs(fun, x0, args, jac, callback, **options)
  File "C:\Users\...\PycharmProjects\pythonProject6\venv\lib\site-packages\scipy\optimize\_optimize.py", line 1419, in _minimize_bfgs
    sf = _prepare_scalar_function(fun, x0, jac, args=args, epsilon=eps,
  File "C:\Users\...\PycharmProjects\pythonProject6\venv\lib\site-packages\scipy\optimize\_optimize.py", line 383, in _prepare_scalar_function
    sf = ScalarFunction(fun, x0, args, grad, hess,
  File "C:\Users\...\PycharmProjects\pythonProject6\venv\lib\site-packages\scipy\optimize\_differentiable_functions.py", line 158, in __init__
    self._update_fun()
  File "C:\Users\...\PycharmProjects\pythonProject6\venv\lib\site-packages\scipy\optimize\_differentiable_functions.py", line 251, in _update_fun
    self._update_fun_impl()
  File "C:\Users\...\PycharmProjects\pythonProject6\venv\lib\site-packages\scipy\optimize\_differentiable_functions.py", line 155, in update_fun
    self.f = fun_wrapped(self.x)
  File "C:\Users\...\PycharmProjects\pythonProject6\venv\lib\site-packages\scipy\optimize\_differentiable_functions.py", line 137, in fun_wrapped
    fx = fun(np.copy(x), *args)
TypeError: dchi() missing 3 required positional arguments: 'b', 'c', and 'd'

Process finished with exit code 1

要解决至少几个问题，请使用一个复杂的，正确地解压缩参数，并检查根收敛：

import functools

import numpy as np
import scipy as sp

def CH(
    chi: float, T: float, a: float, b: float, c: float, d: float, H: float, TC: float, C: float,
) -> float:
    chi = chi + 0j
    err = (
        a*(T - TC)*chi**(1/c)
        + b * chi**(1/c) * chi**(1/d) * H**(1/d)
        + C/(T-TC)
        - 1
    )
    return np.abs(err)

exp = (
    (10.0052, 87.9716), (11.0433, 86.3866), (12.359, 84.8015),
    (169.852,   60.6291), (171.173, 60.6291), (172.494, 60.6291), (173.83,   60.6291),
)

def dchi(params: np.ndarray) -> float:
    a, b, c, d = params
    sd = 0
    x0 = 0
    for exp_a, exp_b in exp:
        CH_bound = functools.partial(CH, T=exp_a, a=a, b=b, c=c, d=d, H=315, TC=1000, C=2)
        ch1 = sp.optimize.root_scalar(
            f=CH_bound, x0=x0,
        )
        assert ch1.converged
        sd += (ch1.root - exp_b)**2
        x0 = ch1.root
    return sd

first_guess = (1,1,1,1)
res = sp.optimize.minimize(
    fun=dchi,
    x0=first_guess,
)
assert res.success, res.message
print(res.x)

这确实更进一步，但在上下文中，您的函数是有问题的，最终不会收敛。你所展示的变量有界限吗？你真的需要找出这些应该是什么。
总之，运行内部根查找是个坏主意。您应该运行一个具有向量化根约束的上层minimize;不要循环播放
你的气的功能

import numpy as np
from scipy.optimize import NonlinearConstraint, minimize

def CH(
    chi: np.ndarray, T: float, H: float, TC: float, C: float, a: float, b: float, c: float, d: float,
) -> np.ndarray:
    err = (
        a*(T - TC)*chi**(1/c)
        + b * chi**(1/c) * chi**(1/d) * H**(1/d)
        + C/(T-TC)
        - 1
    )
    return err

那么一个（缓慢但）成功的优化可能看起来像

def dchi(params: np.ndarray) -> float:
    chi = params[4:]
    chi_err = chi - exp_chi
    err = chi_err.dot(chi_err)
    return err

def dchi_jac(params: np.ndarray) -> np.ndarray:
    chi = params[4:]
    jac = np.zeros_like(params)
    jac[4:] = 2*(chi - exp_chi)
    return jac

def dchi_root(params: np.ndarray) -> np.ndarray:
    # for each value of exp_T, there is some unknown optimal value of chi such that
    # CH(chi, T) == 0
    abcd, chi = np.split(params, (4,))
    a, b, c, d = abcd
    root_err = CH(
        chi=chi, T=exp_T, a=a, b=b, c=c, d=d, H=315, TC=1000, C=2,
    )
    return root_err

exp_T, exp_chi = np.array(exp).T
first_guess = np.concatenate((
    (0.79, 0.00893, 0.428, 1.38),
    exp_chi,
))
result = minimize(
    fun=dchi, jac=dchi_jac,
    x0=first_guess,
    constraints=NonlinearConstraint(
        fun=dchi_root,
        lb=0, ub=0,
    ),
    options={'maxiter': 10, 'disp': True},
)
print(result.message)
print(f'Chi error: {result.fun:.3f}')
print(f'Root error: {np.abs(dchi_root(result.x)).sum():.3f}')
print(f'Parameters: {result.x[:4]}')

Iteration limit reached    (Exit mode 9)
            Current function value: 1543.01943809294
            Iterations: 10
            Function evaluations: 16
            Gradient evaluations: 10
Iteration limit reached
Chi error: 1543.019
Root error: 57.920
Parameters: [-1.61941337e-05 -7.62007792e-06  6.02274220e-01  1.32563818e+00]

如果你给它给予更多的迭代次数，特别是如果你为根写了第二个雅可比矩阵，这会做得更好。

scipy 如何通过多个参数最小化一个函数？(第一部分)

1条答案

相关问题

热门标签

最新问答