英文:
Nonlinear equation system with constraints
问题
from scipy.optimize import least_squares
import numpy as np
def equations(x):
eq1 = x[7] - 1 / (np.exp((5.24677 - 1598.673 / (60.0 - -46.424)) * np.log(10.0)))
eq2 = x[4] - 1 / (np.exp((5.08354 - 1663.125 / (60.0 - -45.622)) * np.log(10.0)))
eq3 = 0.0 - (0.5 * 0.4) + x[2] * x[9] - x[14] * (0.02 * 0.1)
eq4 = 0.0 - (0.5 * 0.6) + x[2] * x[5] - x[12] * (0.02 * 0.1)
eq5 = 0.0 - (0.4 * 0.6) + x[3] * x[0] + x[14] * (0.02 * 0.1)
eq6 = 0.0 - (0.4 * 0.4) + x[3] * x[15] + x[12] * (0.02 * 0.1)
eq7 = 0.0 - (x[10] - x[7] * x[13])
eq8 = 0.0 - (x[8] - x[4] * x[11])
eq9 = x[14] - (1.0 * x[6] * (x[13] - x[9]) + x[9] * (x[14] + x[12]))
eq10 = x[12] - (1.0 * x[6] * (x[11] - x[5]) + x[5] * (x[14] + x[12]))
eq11 = x[14] - (0.01 * x[1] * (x[0] - x[10]) + x[0] * (x[14] + x[12]))
eq12 = x[12] - (0.01 * x[1] * (x[15] - x[8]) + x[15] * (x[14] + x[12]))
eq13 = 1.0 - (x[9] + x[5])
eq14 = 1.0 - (x[0] + x[15])
eq15 = 1.0 - (x[13] + x[11])
eq16 = 1.0 - (x[10] + x[8])
return [eq1, eq2, eq3, eq4, eq5, eq6, eq7, eq8, eq9, eq10, eq11, eq12, eq13, eq14, eq15, eq16]
# constraints
bounds = (np.zeros(16), np.array([1, np.inf, np.inf, np.inf, 1, np.inf, np.inf, 1, 1, 1, 1, 1, np.inf, 1, np.inf, 1]))
# inital
x0 = np.zeros(16)
# x0 = np.array([1, 1e-4, 1, 1, 1, 1, 1e-4, 1, 1, 1, 1, 1, 1, 1, 1, 1])
res = least_squares(equations, x0, bounds=bounds)
print(res.x)
英文:
Hello I want to solve a thermodynamic equation system. Unfortunately, I am quite a beginner in Python and need some help. I tried to solve it with the least_squares method but I always got my initial vector as a result.
I hope someone can help!
from scipy.optimize import least_squares
import numpy as np
def equations(x):
eq1 = x[7] - 1 / (np.exp((5.24677 - 1598.673 / (60.0 - -46.424)) * np.log(10.0)))
eq2 = x[4] - 1 / (np.exp((5.08354 - 1663.125 / (60.0 - -45.622)) * np.log(10.0)))
eq3 = 0.0 - (0.5 * 0.4) + x[2] * x[9] - x[14] * (0.02 * 0.1)
eq4 = 0.0 - (0.5 * 0.6) + x[2] * x[5] - x[12] * (0.02 * 0.1)
eq5 = 0.0 - (0.4 * 0.6) + x[3] * x[0] + x[14] * (0.02 * 0.1)
eq6 = 0.0 - (0.4 * 0.4) + x[3] * x[15] + x[12] * (0.02 * 0.1)
eq7 = 0.0 - (x[10] - x[7] * x[13])
eq8 = 0.0 - (x[8] - x[4] * x[11])
eq9 = x[14] - (1.0 * x[6] * (x[13] - x[9]) + x[9] * (x[14] + x[12]))
eq10 = x[12] - (1.0 * x[6] * (x[11] - x[5]) + x[5] * (x[14] + x[12]))
eq11 = x[14] - (0.01 * x[1] * (x[0] - x[10]) + x[0] * (x[14] + x[12]))
eq12 = x[12] - (0.01 * x[1] * (x[15] - x[8]) + x[15] * (x[14] + x[12]))
eq13 = 1.0 - (x[9] + x[5])
eq14 = 1.0 - (x[0] + x[15])
eq15 = 1.0 - (x[13] + x[11])
eq16 = 1.0 - (x[10] + x[8])
return [eq1, eq2, eq3, eq4, eq5, eq6, eq7, eq8, eq9, eq10, eq11, eq12, eq13, eq14, eq15, eq16]
# constraints
bounds = (np.zeros(16), np.array([1, np.inf, np.inf, np.inf, 1, np.inf, np.inf, 1, 1, 1, 1, 1, np.inf, 1, np.inf, 1]))
# inital
x0 = np.zeros(16)
# x0 = np.array([1, 1e-4, 1, 1, 1, 1, 1e-4, 1, 1, 1, 1, 1, 1, 1, 1, 1])
res = least_squares(equations, x0, bounds=bounds)
print(res.x)
答案1
得分: 1
以下是已经翻译好的部分:
"One probable problem is that the scale of the first and second RHS is enormous compared to the others. x4 and x7 appear to be fixed, but you do not reflect this in your initial values. And indeed, the upper bound for x4 and x7 is 1, which makes absolutely no sense given the first and second equations. So that needs to be fixed for your system of equations to be sensible.
I show below that disregarding the first two equations (reordered to be equations 4 and 7) produces a seemingly reasonable solution."
"这里可能存在的问题是,与其他方程相比,第一个和第二个右侧的比例巨大。x4和x7似乎是固定的,但您在初始值中没有反映这一点。而且,x4和x7的上限是1,这在考虑第一个和第二个方程时毫无意义。因此,必须修正这个问题,使您的方程系统变得合理。
我在下面展示,忽略了前两个方程(重新排列为方程4和7),得到了一个看似合理的解决方案。"
接下来的部分包括Python代码和输出,您可以根据需要研究它们。如果您需要更多翻译,请告诉我需要翻译的具体内容。
英文:
One probable problem is that the scale of the first and second RHS is enormous compared to the others. x4 and x7 appear to be fixed, but you do not reflect this in your initial values. And indeed, the upper bound for x4 and x7 is 1, which makes absolutely no sense given the first and second equations. So that needs to be fixed for your system of equations to be sensible.
I show below that disregarding the first two equations (reordered to be equations 4 and 7) produces a seemingly reasonable solution.
from scipy.optimize import least_squares
import numpy as np
def sides(x: np.ndarray) -> np.ndarray:
# x order as original, y reordered so that x indices are closer to their use
return np.array((
(x[14], 0.01*x[1]*(x[ 0] - x[10]) + x[ 0]*(x[14] + x[12])),
(x[12], 0.01*x[1]*(x[15] - x[ 8]) + x[15]*(x[14] + x[12])),
(0, -0.5*0.6 + x[2]* x[ 5] - x[12]*(0.02 * 0.1)),
(0, -0.5*0.4 + x[2]* x[ 9] - x[14]*(0.02 * 0.1)),
(x[4], 1 / np.exp((5.08354 - 1663.125/(60.0 - -45.622)) * np.log(10))),
(x[12], x[6]*(x[11] - x[5]) + x[5]*(x[14] + x[12])),
(x[14], x[6]*(x[13] - x[9]) + x[9]*(x[14] + x[12])),
(x[7], 1 / np.exp((5.24677 - 1598.673/(60.0 - -46.424)) * np.log(10))),
(1, x[10] + x[ 8]),
(1, x[ 9] + x[ 5]),
(0, x[10] - x[7]*x[13]),
(0, x[ 8] - x[4]*x[11]),
(0, -0.4*0.4 + x[3]* x[15] + x[12]*(0.02 * 0.1)),
(1, x[13] + x[11]),
(0, -0.4*0.6 + x[3]* x[ 0] + x[14]*(0.02 * 0.1)),
(1, x[ 0] + x[15]),
))
def equations(x: np.ndarray) -> np.ndarray:
left, right = sides(x)[[
0, 1, 2, 3,
5, 6,
8, 9, 10, 11, 12, 13, 14, 15,
], :].T
return left - right
# Ordered as original
x0 = np.array((
1, 1e-4, 1, 1,
10, # 1 / np.exp((5.08354 - 1663.125/(60.0 - -45.622)) * np.log(10)),
1, 1e-4,
10, # 1 / np.exp((5.24677 - 1598.673/(60.0 - -46.424)) * np.log(10)),
1, 1, 1, 1, 1, 1, 1, 1,
))
# Ordered as original, with x[4] and x[7] bounds changed to inf (not fixed)
bounds = np.array((
np.zeros(16),
(
1, np.inf, np.inf, np.inf,
np.inf, # 1,
np.inf, np.inf,
np.inf, # 1,
1, 1, 1, 1, np.inf, 1, np.inf, 1,
),
))
res = least_squares(fun=equations, x0=x0, bounds=bounds)
assert res.success, res.message
for i, (x, (left, right)) in enumerate(zip(res.x, sides(res.x))):
print(f'y{i:<2} = {left:8.3g} ~ {right:9.3g}, '
f'err={right-left:8.1e}, '
f'x{i:<2} = {x:.3g}')
y0 = 6.25 ~ 6.25, err= 5.3e-15, x0 = 0.591
y1 = 1.41 ~ 1.41, err= 0.0e+00, x1 = 361
y2 = 0 ~ -1.25e-13, err=-1.3e-13, x2 = 0.515
y3 = 0 ~ 1.78e-13, err= 1.8e-13, x3 = 0.385
y4 = 1.56 ~ 4.6e+10, err= 4.6e+10, x4 = 1.56
y5 = 1.41 ~ 1.41, err=-6.9e-15, x5 = 0.588
y6 = 6.25 ~ 6.25, err= 1.2e-14, x6 = 139
y7 = 0.266 ~ 5.96e+09, err= 6.0e+09, x7 = 0.266
y8 = 1 ~ 1, err=-2.2e-15, x8 = 0.884
y9 = 1 ~ 1, err= 1.9e-13, x9 = 0.412
y10 = 0 ~ 4.16e-17, err= 4.2e-17, x10 = 0.116
y11 = 0 ~ 2.22e-16, err= 2.2e-16, x11 = 0.565
y12 = 0 ~ 1.07e-13, err= 1.1e-13, x12 = 1.41
y13 = 1 ~ 1, err= 1.8e-13, x13 = 0.435
y14 = 0 ~ -7.43e-14, err=-7.4e-14, x14 = 6.25
y15 = 1 ~ 1, err=-2.2e-16, x15 = 0.409
This sidekick program will show you some options where values in the two problematic equations could be substituted to produce a reasonable solution in the original problem:
import numpy as np
from sympy import Symbol, Eq, solveset
for y, values in (
(1.690, (5.08354, 1663.125, 60, -45.622)),
(0.138, (5.24677, 1598.673, 60, -46.424)),
):
a = Symbol('a', real=True, infinite=False)
b = Symbol('b', real=True, infinite=False)
c = Symbol('c', real=True, infinite=False)
d = Symbol('d', real=True, infinite=False)
variables = a, b, c, d
eq = Eq(-np.log10(y), a - b/(c - d))
print('Do one of:')
for target, old_val in zip(variables, values):
seq = eq.subs({
s: v for s, v in zip(variables, values)
if s is not target
})
replacement, = solveset(seq, target)
print(f'Replace {target}={old_val} with {replacement:.1f}')
print()
Do one of:
Replace a=5.08354 with 15.5
Replace b=1663.125 with 561.0
Replace c=60 with 267.5
Replace d=-45.622 with -253.1
Do one of:
Replace a=5.24677 with 15.9
Replace b=1598.673 with 466.8
Replace c=60 with 318.0
Replace d=-46.424 with -304.4
通过集体智慧和协作来改善编程学习和解决问题的方式。致力于成为全球开发者共同参与的知识库,让每个人都能够通过互相帮助和分享经验来进步。
评论