将 sympy 表达式转换为向量以找到线性独立的子集
Turning a sympy expression into a vector to find linearly independent subset
我有一个表达式列表,例如 4.0*x[0] + 5.0*x[10] + 1 = 0
我想根据 [4.0, 0, 0, ..., 5.0, ... , 1]
等系数将它们转换为向量。原因是我的一些方程可能是线性相关的,我想从 numpy 库中 运行 QR 这样我就可以找到一个线性独立的子集。
我可以通过 expr.replace(x[i], 0)
和 i
通配符索引来获得常数项。我还可以通过 expr.atoms(Mul)
得到大部分其他术语,这给了我集合 4.0*x[0], 5.0*x[10]
然后对于这些表达式中的每一个,我可以做 expr.atoms(Indexed).pop()
和 expr.atoms(Float).pop()
来拆分部分.
问题是当我有一个像 x[0] + 5.0*x[10] + 1 = 0
这样的表达式时,其中第一个变量出现的隐式系数为 1
。该术语不再被识别为 Mul
对象。
无论如何,我认为可能有更好的方法来实现我的目标?
如果你给你的符号一个特定的顺序,就像下面的代码一样,你可以将表达式转换为多项式并得到它的系数:
>>> from sympy import *
>>> x, y, z, t = symbols('x y z t')
>>> a1, a2, a3, a4 = symbols('a[1], a[2], a[3], a[4]')
>>> used_symbols = (a1, a2, a3, a4)
>>> replacements = [(n, x**(enu+1)) for enu,n in enumerate(used_symbols)]
>>> expr = 5 + a1 + 4*a4
>>> Poly(expr.subs(replacements)).all_coeffs()
[4, 0, 0, 1, 5]
如果事先不知道使用以下递归函数,您也可以检索所用符号的列表:
def retrieve_used_symbols(expr):
"""Return the symbols used in the `expr` in a list."""
used_symbols = []
for term in expr.args:
if term.is_Atom and term.is_Symbol:
used_symbols.append(term)
else:
used_symbols.extend(retrieve_used_symbols(term))
return used_symbols
当你有混合符号时,后者会派上用场:
>>> crazy_expr = expr + 10*y-2*z
>>> crazy_expr
a[1] + 4*a[4] + 10*y - 2*z + 5
>>> used_symbols = retrieve_used_symbols(crazy_expr)
>>> replacements = [(n, x**(enu+1)) for enu,n in enumerate(used_symbols)]
>>> Poly(crazy_expr.subs(replacements)).all_coeffs()
[4, -2, 1, 10, 5]
>>> list(reversed(used_symbols))
[a[4], z, a[1], y]
对于IndexedBase对象,更简单:
coeffs = [expr.coeff(x[i]) for i in range(10)]
但是您仍然需要添加常数项,如您所说,您可以通过通配符替换获得该常数项:
ind = Wild('i')
constant_term = expr.replace(x[ind], 0)
{根据@(Oliver W.) 的要求}
给出
>>> x = IndexedBase('x')
>>> eqs = 4*x[0] + 5*x[5] + 1, x[1] - x[2]
>>> v = list(ordered(Tuple(*eqs).atoms(Indexed)))
可以这样做
>>> [[eq.coeff(vi) for vi in v] + [eq.as_coeff_Add()[0]] for eq in eqs]
[[4, 0, 0, 5, 1], [0, 1, -1, 0, 0]]
但是其中大部分可以通过矩阵方法获得jacobian
。但是要使用它你必须用符号替换 x[i] (因为 diff 只适用于函数是符号,IIRC):
>>> d = [Dummy() for vi in v]
>>> z = dict(zip(d, [0]*len(d)))
>>> m = Matrix([eq.xreplace(dict(zip(v, d))) for eq in eqs])
>>> m.jacobian(d)
Matrix([
[4, 0, 0, 5],
[0, 1, -1, 0]])
>>> m.subs(z)
Matrix([
[1],
[0]])
我有一个表达式列表,例如 4.0*x[0] + 5.0*x[10] + 1 = 0
我想根据 [4.0, 0, 0, ..., 5.0, ... , 1]
等系数将它们转换为向量。原因是我的一些方程可能是线性相关的,我想从 numpy 库中 运行 QR 这样我就可以找到一个线性独立的子集。
我可以通过 expr.replace(x[i], 0)
和 i
通配符索引来获得常数项。我还可以通过 expr.atoms(Mul)
得到大部分其他术语,这给了我集合 4.0*x[0], 5.0*x[10]
然后对于这些表达式中的每一个,我可以做 expr.atoms(Indexed).pop()
和 expr.atoms(Float).pop()
来拆分部分.
问题是当我有一个像 x[0] + 5.0*x[10] + 1 = 0
这样的表达式时,其中第一个变量出现的隐式系数为 1
。该术语不再被识别为 Mul
对象。
无论如何,我认为可能有更好的方法来实现我的目标?
如果你给你的符号一个特定的顺序,就像下面的代码一样,你可以将表达式转换为多项式并得到它的系数:
>>> from sympy import *
>>> x, y, z, t = symbols('x y z t')
>>> a1, a2, a3, a4 = symbols('a[1], a[2], a[3], a[4]')
>>> used_symbols = (a1, a2, a3, a4)
>>> replacements = [(n, x**(enu+1)) for enu,n in enumerate(used_symbols)]
>>> expr = 5 + a1 + 4*a4
>>> Poly(expr.subs(replacements)).all_coeffs()
[4, 0, 0, 1, 5]
如果事先不知道使用以下递归函数,您也可以检索所用符号的列表:
def retrieve_used_symbols(expr):
"""Return the symbols used in the `expr` in a list."""
used_symbols = []
for term in expr.args:
if term.is_Atom and term.is_Symbol:
used_symbols.append(term)
else:
used_symbols.extend(retrieve_used_symbols(term))
return used_symbols
当你有混合符号时,后者会派上用场:
>>> crazy_expr = expr + 10*y-2*z
>>> crazy_expr
a[1] + 4*a[4] + 10*y - 2*z + 5
>>> used_symbols = retrieve_used_symbols(crazy_expr)
>>> replacements = [(n, x**(enu+1)) for enu,n in enumerate(used_symbols)]
>>> Poly(crazy_expr.subs(replacements)).all_coeffs()
[4, -2, 1, 10, 5]
>>> list(reversed(used_symbols))
[a[4], z, a[1], y]
对于IndexedBase对象,更简单:
coeffs = [expr.coeff(x[i]) for i in range(10)]
但是您仍然需要添加常数项,如您所说,您可以通过通配符替换获得该常数项:
ind = Wild('i')
constant_term = expr.replace(x[ind], 0)
{根据@(Oliver W.) 的要求}
给出
>>> x = IndexedBase('x')
>>> eqs = 4*x[0] + 5*x[5] + 1, x[1] - x[2]
>>> v = list(ordered(Tuple(*eqs).atoms(Indexed)))
可以这样做
>>> [[eq.coeff(vi) for vi in v] + [eq.as_coeff_Add()[0]] for eq in eqs]
[[4, 0, 0, 5, 1], [0, 1, -1, 0, 0]]
但是其中大部分可以通过矩阵方法获得jacobian
。但是要使用它你必须用符号替换 x[i] (因为 diff 只适用于函数是符号,IIRC):
>>> d = [Dummy() for vi in v]
>>> z = dict(zip(d, [0]*len(d)))
>>> m = Matrix([eq.xreplace(dict(zip(v, d))) for eq in eqs])
>>> m.jacobian(d)
Matrix([
[4, 0, 0, 5],
[0, 1, -1, 0]])
>>> m.subs(z)
Matrix([
[1],
[0]])