如何使用 sympy 从标准的作业行序列中本地切换到使用 "set of equations"?
How to use sympy to locally switch to use of "set of equations" from standard sequence of lines of assignments?
以下 Python 函数在计算 f(x) 的步骤顺序中存在缺陷。
def f(x):
f = np.zeros((3))
f[0] = x[0]
f[1] = f[2] + x[1]
f[2] = x[2]
return f
输入 x= [1,1,1] 给出 [1,1,1] 但应该给出 [1,2,1] 如果我们将代码从数学上解释为方程组。如何使用 sympy 整理评估顺序并获得“正确”答案?
请注意,我根本不是在寻找符号式求解方程组的方法,只是希望求值顺序正确,避免自己进行映射。
您从根本上误解了 =
在 Python 中的作用。这不是数学上的等式(双向都有效),而是 assignment。当你说 a = b
时,你说“忘记 a
是什么,从现在开始就是 b
”。
所以一步一步我们得到:
f[0] = x[0] # f[0] becomes x[0] = 1
f[1] = f[2] + x[1] # f[1] becomes f[2] + x[1] = 0 + 1 = 1
f[2] = x[2] # f[2] becomes x[2] = 1
Note, I am not at all looking for solving the equations system symbolically, but just want the evaluation order correct and avoid work with mapping out that myself.
您要么需要自己符号式求解方程,要么使用计算机代数系统来求解。常规 Python 代码是命令式的,它们都不会(开箱即用)。
参见例如kabooya 对使用 Sympy 建立这样一个方程组的回答。
你可以让sympy计算f0、f1、f2的公式:
import sympy as sym
f0, f1, f2 = sym.symbols('f0,f1,f2')
x0, x1, x2 = sym.symbols('x0,x1,x2')
eq1 = sym.Eq(f0, x0)
eq2 = sym.Eq(f1, f2 + x1)
eq3 = sym.Eq(f2, x2)
result = sym.solve([eq1,eq2, eq3],(f0, f1, f2))
print(result)
>>> {f0: x0, f1: x1 + x2, f2: x2}
使用 sympy,您可以将函数编写为:
import numpy as np
def f(x):
from sympy import symbols, Eq, solve
f = symbols('f0:3')
sol = solve([Eq(f[0], x[0]), Eq(f[1], f[2] + x[1]), Eq(f[2], x[2])], f)
return np.array([sol[fi].evalf() for fi in f], dtype=float)
f(np.array([1, 1, 1]))
这假设您只使用了一个简单的例子,而实际上方程式要复杂得多。请注意,如果您的方程式有多个有效解,则该方法可能会涉及更多。
lambdify()
可以将表达式转换为 numpy 格式。为了在这里工作,函数 f
应该获取一些符号输入参数并将它们转换为输出表达式。
这是它的样子。
import numpy as np
from sympy import symbols, Eq, solve, lambdify
def f(x):
y = symbols('y0:3')
sol = solve([Eq(y[0], x[0]), Eq(y[1], y[2] + x[1]), Eq(y[2], x[2])], y)
return [sol[yi].simplify() for yi in y]
x0, x1, x2 = symbols('x0:3')
f_np = lambdify((x0, x1, x2), f((x0, x1, x2)))
# calling the function on a 1D numpy array
f_np(*np.array([1, 1, 1]))
# calling the function on a 2D numpy array and convert the result back to a numpy array
a = np.array([[1, 1, 1], [2, 2, 2]])
np.array(f_np(a[:, 0], a[:, 1], a[:, 2])).T # array([[1, 2, 1], [2, 4, 2]])
这里,f_np
是一个 numpy 函数,它获得三个输入(*
unpacks 数组)并有一个 3 的列表作为输出。目前似乎无法直接使用数组作为输入和输出。
在交互式环境中,您可以输入f_np?
(或f_np??
)来查看函数来源:
>> f_np?
Signature: f_np(x0, x1, x2)
Docstring:
Created with lambdify. Signature:
func(x0, x1, x2)
Expression:
(x0, x1 + x2, x2)
Source code:
def _lambdifygenerated(x0, x1, x2):
return ((x0, x1 + x2, x2))
f_np
然后可以在 scipy 优化器中使用,甚至可以编译(Cython 或 numba)。这种方法假设方程式保持不变。
以下 Python 函数在计算 f(x) 的步骤顺序中存在缺陷。
def f(x):
f = np.zeros((3))
f[0] = x[0]
f[1] = f[2] + x[1]
f[2] = x[2]
return f
输入 x= [1,1,1] 给出 [1,1,1] 但应该给出 [1,2,1] 如果我们将代码从数学上解释为方程组。如何使用 sympy 整理评估顺序并获得“正确”答案?
请注意,我根本不是在寻找符号式求解方程组的方法,只是希望求值顺序正确,避免自己进行映射。
您从根本上误解了 =
在 Python 中的作用。这不是数学上的等式(双向都有效),而是 assignment。当你说 a = b
时,你说“忘记 a
是什么,从现在开始就是 b
”。
所以一步一步我们得到:
f[0] = x[0] # f[0] becomes x[0] = 1
f[1] = f[2] + x[1] # f[1] becomes f[2] + x[1] = 0 + 1 = 1
f[2] = x[2] # f[2] becomes x[2] = 1
Note, I am not at all looking for solving the equations system symbolically, but just want the evaluation order correct and avoid work with mapping out that myself.
您要么需要自己符号式求解方程,要么使用计算机代数系统来求解。常规 Python 代码是命令式的,它们都不会(开箱即用)。
参见例如kabooya 对使用 Sympy 建立这样一个方程组的回答。
你可以让sympy计算f0、f1、f2的公式:
import sympy as sym
f0, f1, f2 = sym.symbols('f0,f1,f2')
x0, x1, x2 = sym.symbols('x0,x1,x2')
eq1 = sym.Eq(f0, x0)
eq2 = sym.Eq(f1, f2 + x1)
eq3 = sym.Eq(f2, x2)
result = sym.solve([eq1,eq2, eq3],(f0, f1, f2))
print(result)
>>> {f0: x0, f1: x1 + x2, f2: x2}
使用 sympy,您可以将函数编写为:
import numpy as np
def f(x):
from sympy import symbols, Eq, solve
f = symbols('f0:3')
sol = solve([Eq(f[0], x[0]), Eq(f[1], f[2] + x[1]), Eq(f[2], x[2])], f)
return np.array([sol[fi].evalf() for fi in f], dtype=float)
f(np.array([1, 1, 1]))
这假设您只使用了一个简单的例子,而实际上方程式要复杂得多。请注意,如果您的方程式有多个有效解,则该方法可能会涉及更多。
lambdify()
可以将表达式转换为 numpy 格式。为了在这里工作,函数 f
应该获取一些符号输入参数并将它们转换为输出表达式。
这是它的样子。
import numpy as np
from sympy import symbols, Eq, solve, lambdify
def f(x):
y = symbols('y0:3')
sol = solve([Eq(y[0], x[0]), Eq(y[1], y[2] + x[1]), Eq(y[2], x[2])], y)
return [sol[yi].simplify() for yi in y]
x0, x1, x2 = symbols('x0:3')
f_np = lambdify((x0, x1, x2), f((x0, x1, x2)))
# calling the function on a 1D numpy array
f_np(*np.array([1, 1, 1]))
# calling the function on a 2D numpy array and convert the result back to a numpy array
a = np.array([[1, 1, 1], [2, 2, 2]])
np.array(f_np(a[:, 0], a[:, 1], a[:, 2])).T # array([[1, 2, 1], [2, 4, 2]])
这里,f_np
是一个 numpy 函数,它获得三个输入(*
unpacks 数组)并有一个 3 的列表作为输出。目前似乎无法直接使用数组作为输入和输出。
在交互式环境中,您可以输入f_np?
(或f_np??
)来查看函数来源:
>> f_np?
Signature: f_np(x0, x1, x2)
Docstring:
Created with lambdify. Signature:
func(x0, x1, x2)
Expression:
(x0, x1 + x2, x2)
Source code:
def _lambdifygenerated(x0, x1, x2):
return ((x0, x1 + x2, x2))
f_np
然后可以在 scipy 优化器中使用,甚至可以编译(Cython 或 numba)。这种方法假设方程式保持不变。