sympy solve（）は暗黙的/間違った答えを返します

Question

私はsympyを使って16の方程式と16の未知数で方程式系を解こうとしていますが、それはうまく解けないようです。

[K] [d] = [f]ここで[K]は係数行列、[d]は未知数、[f]は定数です。私はいくつかの未知数 "d"といくつかの定数 "f"を知っているので、方程式と未知数の数は同じですが、これらの値を方程式に代入して解くと、 ""の結果は "dx8" 。行列の行列式を調べたところ、正の値をとるので、固有の答えを得なければなりません。ここで

はコードです：

import sympy as sp 
import numpy as np 

K = np.array([[560000000.0, 0.0, -480000000.0, 80000000.0, 0.0, 0.0, 0.0, 0.0,-80000000.0, 120000000.0, 0.0, -200000000.0, 0.0, 0.0, 0.0, 0.0], 
    [0.0, 393333333.3, 120000000.0, -180000000.0, 0.0, 0.0, 0.0, 0.0,80000000.0, -213333333.3, -200000000.0, 0.0, 0.0, 0.0, 0.0, 0.0], 
    [-480000000.0, 120000000.0, 1120000000.0, -200000000.0,-480000000.0, 80000000.0, 0.0, 0.0, 0.0, 0.0, -160000000.0,200000000.0, 0.0, -200000000.0, 0.0, 0.0], 
    [80000000.0, -180000000.0, -200000000.0, 786666666.7, 120000000.0,-180000000.0, 0.0, 0.0, 0.0, 0.0, 200000000.0, -426666666.7,-200000000.0, 0.0, 0.0, 0.0], 
    [0.0, 0.0, -480000000.0, 120000000.0, 1120000000.0, -200000000.0,-480000000.0, 80000000.0, 0.0, 0.0, 0.0, 0.0, -160000000.0,200000000.0, 0.0, -200000000.0], 
    [0.0, 0.0, 80000000.0, -180000000.0, -200000000.0, 786666666.7,120000000.0, -180000000.0, 0.0, 0.0, 0.0, 0.0, 200000000.0,-426666666.7, -200000000.0, 0.0], 
    [0.0, 0.0, 0.0, 0.0, -480000000.0, 120000000.0, 560000000.0,-200000000.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, -80000000.0, 80000000.0], 
    [0.0, 0.0, 0.0, 0.0, 80000000.0, -180000000.0, -200000000.0,393333333.3, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 120000000.0,-213333333.3], 
    [-80000000.0, 80000000.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 560000000.0,-200000000.0, -480000000.0, 120000000.0, 0.0, 0.0, 0.0, 0.0], 
    [120000000.0, -213333333.3, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,-200000000.0, 393333333.3, 80000000.0, -180000000.0, 0.0, 0.0, 0.0,0.0], 
    [0.0, -200000000.0, -160000000.0, 200000000.0, 0.0, 0.0, 0.0, 0.0,-480000000.0, 80000000.0, 1120000000.0, -200000000.0, -480000000.0,120000000.0, 0.0, 0.0], 
    [-200000000.0, 0.0, 200000000.0, -426666666.7, 0.0, 0.0, 0.0, 0.0,120000000.0, -180000000.0, -200000000.0, 786666666.7, 80000000.0,-180000000.0, 0.0, 0.0], 
    [0.0, 0.0, 0.0, -200000000.0, -160000000.0, 200000000.0, 0.0, 0.0,0.0, 0.0, -480000000.0, 80000000.0, 1120000000.0, -200000000.0,-480000000.0, 120000000.0], 
    [0.0, 0.0, -200000000.0, 0.0, 200000000.0, -426666666.7, 0.0, 0.0,0.0, 0.0, 120000000.0, -180000000.0, -200000000.0, 786666666.7,80000000.0, -180000000.0], 
    [0.0, 0.0, 0.0, 0.0, 0.0, -200000000.0, -80000000.0, 120000000.0,0.0, 0.0, 0.0, 0.0, -480000000.0, 80000000.0, 560000000.0, 0.0], 
    [0.0, 0.0, 0.0, 0.0, -200000000.0, 0.0, 80000000.0, -213333333.3,0.0, 0.0, 0.0, 0.0, 120000000.0, -180000000.0, 0.0, 393333333.3]]) 



x = [sp.var('dx'+ str(i+1)) for i in range(8)] 
y = [sp.var('dy'+ str(i+1)) for i in range(8)] 
fx = [sp.var('fx'+ str(i+1)) for i in range(8)] 
fy = [sp.var('fy'+ str(i+1)) for i in range(8)] 

xy = list(sum(zip(x, y),())) 
fxy = list(sum(zip(fx, fy),())) 

M = sp.Matrix(K)*sp.Matrix(xy) 
Ec = [sp.Eq(M[i], fxy[i]) for i in range(16)] 

#known values 
d_kwn = [(dy1, 0), (dy2, 0), (dy3, 0), (dy4, 0)] 

f_kwn = [(fx5, 0), (fy5, 0), (fx6, 0), (fy6, -3000), (fx7, 0), (fy7, -3000),(fx8, 0), (fy8, 0), (fx1, 0), (fx2, 0), (fx3, 0), (fx4, 0)] 

for var in d_kwn: 
    for i, eq in enumerate(Ec): 
     Ec[i] = eq.subs(var[0], var[1]) 

for var in f_kwn: 
    for i, eq in enumerate(Ec): 
     Ec[i] = eq.subs(var[0], var[1]) 

Sols = sp.solvers.solve(Ec) 
sp.Matrix(sorted(Sols.items(), key=str)) 

そして、これは私が取得しています出力されます：

{dx1: dx8−3.54468009860439⋅10−6, 
dx2: dx8−1.8414987360977⋅10−6, 
dx3: dx8−2.11496606381994⋅10−7, 
dx4: dx8+2.05943267588118⋅10−7, 
dx5: dx8−1.24937663359153⋅10−6, 
dx6: dx8−1.55655946713284⋅10−6, 
dx7: dx8−1.08797652070783⋅10−6, 
dy5: −2.10639657360695⋅10−6, 
dy6: −6.26959460018537⋅10−6, 
dy7: −6.32191585665888⋅10−6, 
dy8: −2.7105825114088⋅10−6, 
fy1: 439.746516706791, 
fy2: 2640.65618690176, 
fy3: 2399.44807607611, 
fy4: 520.14922031534} 

私はdx8ため、結果が届かない理由を私は知りません。理論的にはdx1 = dx4、dx2 = dx3、dx5 = dx8、dx6 = dx7などのように多くの方程式を追加しようとしました。しかしそれは私と空リストを与える。 助けていただければ幸いです。

Answer 1

あなたがSympyを使用する必要がある場合は、次のように動作することがあります。我々はすでに、次のように我々は方程式の縮小セットを解決するためにnumpyの配列のインデックスを使用することができます知られているfとdの要素を知っています。まず、未知数d値の場合にのみ、減少した方程式系を解くことができる。そして、dの値がすべてわかったら、未知数（）を未知数f（以下のコードでは実装されていません）に対してのみ実行して、未知のfの値を計算することができます。

import sympy as sp 
import numpy as np 

K = np.array([[560000000.0, 0.0, -480000000.0, 80000000.0, 0.0, 0.0, 0.0, 0.0,-80000000.0, 120000000.0, 0.0, -200000000.0, 0.0, 0.0, 0.0, 0.0], 
    [0.0, 393333333.3, 120000000.0, -180000000.0, 0.0, 0.0, 0.0, 0.0,80000000.0, -213333333.3, -200000000.0, 0.0, 0.0, 0.0, 0.0, 0.0], 
    [-480000000.0, 120000000.0, 1120000000.0, -200000000.0,-480000000.0, 80000000.0, 0.0, 0.0, 0.0, 0.0, -160000000.0,200000000.0, 0.0, -200000000.0, 0.0, 0.0], 
    [80000000.0, -180000000.0, -200000000.0, 786666666.7, 120000000.0,-180000000.0, 0.0, 0.0, 0.0, 0.0, 200000000.0, -426666666.7,-200000000.0, 0.0, 0.0, 0.0], 
    [0.0, 0.0, -480000000.0, 120000000.0, 1120000000.0, -200000000.0,-480000000.0, 80000000.0, 0.0, 0.0, 0.0, 0.0, -160000000.0,200000000.0, 0.0, -200000000.0], 
    [0.0, 0.0, 80000000.0, -180000000.0, -200000000.0, 786666666.7,120000000.0, -180000000.0, 0.0, 0.0, 0.0, 0.0, 200000000.0,-426666666.7, -200000000.0, 0.0], 
    [0.0, 0.0, 0.0, 0.0, -480000000.0, 120000000.0, 560000000.0,-200000000.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, -80000000.0, 80000000.0], 
    [0.0, 0.0, 0.0, 0.0, 80000000.0, -180000000.0, -200000000.0,393333333.3, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 120000000.0,-213333333.3], 
    [-80000000.0, 80000000.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 560000000.0,-200000000.0, -480000000.0, 120000000.0, 0.0, 0.0, 0.0, 0.0], 
    [120000000.0, -213333333.3, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,-200000000.0, 393333333.3, 80000000.0, -180000000.0, 0.0, 0.0, 0.0,0.0], 
    [0.0, -200000000.0, -160000000.0, 200000000.0, 0.0, 0.0, 0.0, 0.0,-480000000.0, 80000000.0, 1120000000.0, -200000000.0, -480000000.0,120000000.0, 0.0, 0.0], 
    [-200000000.0, 0.0, 200000000.0, -426666666.7, 0.0, 0.0, 0.0, 0.0,120000000.0, -180000000.0, -200000000.0, 786666666.7, 80000000.0,-180000000.0, 0.0, 0.0], 
    [0.0, 0.0, 0.0, -200000000.0, -160000000.0, 200000000.0, 0.0, 0.0,0.0, 0.0, -480000000.0, 80000000.0, 1120000000.0, -200000000.0,-480000000.0, 120000000.0], 
    [0.0, 0.0, -200000000.0, 0.0, 200000000.0, -426666666.7, 0.0, 0.0,0.0, 0.0, 120000000.0, -180000000.0, -200000000.0, 786666666.7,80000000.0, -180000000.0], 
    [0.0, 0.0, 0.0, 0.0, 0.0, -200000000.0, -80000000.0, 120000000.0,0.0, 0.0, 0.0, 0.0, -480000000.0, 80000000.0, 560000000.0, 0.0], 
    [0.0, 0.0, 0.0, 0.0, -200000000.0, 0.0, 80000000.0, -213333333.3,0.0, 0.0, 0.0, 0.0, 120000000.0, -180000000.0, 0.0, 393333333.3]]) 



x = [sp.var('dx'+ str(i+1)) for i in range(8)] 
y = [sp.var('dy'+ str(i+1)) for i in range(8)] 
fx = [sp.var('fx'+ str(i+1)) for i in range(8)] 
fy = [sp.var('fy'+ str(i+1)) for i in range(8)] 

xy = list(sum(zip(x, y),())) 
fxy = list(sum(zip(fx, fy),())) 

M = sp.Matrix(K)*sp.Matrix(xy) 
Ec = [sp.Eq(M[i], fxy[i]) for i in range(16)] 

#known values 
d_kwn = [(dy1, 0), (dy2, 0), (dy3, 0), (dy4, 0)] 

f_kwn = [(fx5, 0), (fy5, 0), (fx6, 0), (fy6, -3000), (fx7, 0), (fy7, -3000),(fx8, 0), (fy8, 0), (fx1, 0), (fx2, 0), (fx3, 0), (fx4, 0)] 

for var in d_kwn: 
    for i, eq in enumerate(Ec): 
     Ec[i] = eq.subs(var[0], var[1]) 

for var in f_kwn: 
    for i, eq in enumerate(Ec): 
     Ec[i] = eq.subs(var[0], var[1]) 

Ec_part = [] 
for i in [0,2,4,6,8,9,10,11,12,13,14,15]: 
    Ec_part.append(Ec[i]) 

unknwns = [*x, *y[4:8]] 

Sols = sp.linsolve(Ec_part,unknwns) 

Sols = next(iter(Sols)) 

#sp.Matrix(sorted(Sols.items(), key=str)) 

Answer 2

Numpy自体で線形方程式の系を解くのが便利です。あなたが解決しているシステムのタイプは、しばしば境界条件を持つ有限要素解析に現れます。私たちだけがナンシーを使用すればいいのですか？はいの場合、次のコードでジョブを実行します。

import numpy as np 

# The NxN Coefficients matrix 
K = np.array([[560000000.0, 0.0, -480000000.0, 80000000.0, 0.0, 0.0, 0.0, 0.0,-80000000.0, 120000000.0, 0.0, -200000000.0, 0.0, 0.0, 0.0, 0.0], 
    [0.0, 393333333.3, 120000000.0, -180000000.0, 0.0, 0.0, 0.0, 0.0,80000000.0, -213333333.3, -200000000.0, 0.0, 0.0, 0.0, 0.0, 0.0], 
    [-480000000.0, 120000000.0, 1120000000.0, -200000000.0,-480000000.0, 80000000.0, 0.0, 0.0, 0.0, 0.0, -160000000.0,200000000.0, 0.0, -200000000.0, 0.0, 0.0], 
    [80000000.0, -180000000.0, -200000000.0, 786666666.7, 120000000.0,-180000000.0, 0.0, 0.0, 0.0, 0.0, 200000000.0, -426666666.7,-200000000.0, 0.0, 0.0, 0.0], 
    [0.0, 0.0, -480000000.0, 120000000.0, 1120000000.0, -200000000.0,-480000000.0, 80000000.0, 0.0, 0.0, 0.0, 0.0, -160000000.0,200000000.0, 0.0, -200000000.0], 
    [0.0, 0.0, 80000000.0, -180000000.0, -200000000.0, 786666666.7,120000000.0, -180000000.0, 0.0, 0.0, 0.0, 0.0, 200000000.0,-426666666.7, -200000000.0, 0.0], 
    [0.0, 0.0, 0.0, 0.0, -480000000.0, 120000000.0, 560000000.0,-200000000.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, -80000000.0, 80000000.0], 
    [0.0, 0.0, 0.0, 0.0, 80000000.0, -180000000.0, -200000000.0,393333333.3, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 120000000.0,-213333333.3], 
    [-80000000.0, 80000000.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 560000000.0,-200000000.0, -480000000.0, 120000000.0, 0.0, 0.0, 0.0, 0.0], 
    [120000000.0, -213333333.3, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,-200000000.0, 393333333.3, 80000000.0, -180000000.0, 0.0, 0.0, 0.0,0.0], 
    [0.0, -200000000.0, -160000000.0, 200000000.0, 0.0, 0.0, 0.0, 0.0,-480000000.0, 80000000.0, 1120000000.0, -200000000.0, -480000000.0,120000000.0, 0.0, 0.0], 
    [-200000000.0, 0.0, 200000000.0, -426666666.7, 0.0, 0.0, 0.0, 0.0,120000000.0, -180000000.0, -200000000.0, 786666666.7, 80000000.0,-180000000.0, 0.0, 0.0], 
    [0.0, 0.0, 0.0, -200000000.0, -160000000.0, 200000000.0, 0.0, 0.0,0.0, 0.0, -480000000.0, 80000000.0, 1120000000.0, -200000000.0,-480000000.0, 120000000.0], 
    [0.0, 0.0, -200000000.0, 0.0, 200000000.0, -426666666.7, 0.0, 0.0,0.0, 0.0, 120000000.0, -180000000.0, -200000000.0, 786666666.7,80000000.0, -180000000.0], 
    [0.0, 0.0, 0.0, 0.0, 0.0, -200000000.0, -80000000.0, 120000000.0,0.0, 0.0, 0.0, 0.0, -480000000.0, 80000000.0, 560000000.0, 0.0], 
    [0.0, 0.0, 0.0, 0.0, -200000000.0, 0.0, 80000000.0, -213333333.3,0.0, 0.0, 0.0, 0.0, 120000000.0, -180000000.0, 0.0, 393333333.3]]) 

# A logical array for indexing 
N = K.shape[0] # The number of columns in K 
N_2 = int(N/2); 

# Prepare the 'f' 
fx = np.zeros(N_2); 
fy = np.zeros(N_2); 

fx[ [0,1,2,3,4,5,6,7] ] = np.array([0]*N_2) # Known values of fx 
fy[ [4,5,6,7] ] = np.array([0,-3000,-3000,0]) 

f = np.concatenate((fx,fy)) 

# Solve for the unknown equations only 
d = np.zeros(N) 
rows = np.array([0,1,2,3,4,5,6,7,12,13,14,15]) 
rows = rows[:, np.newaxis] 
columns = np.array([0,1,2,3,4,5,6,7,12,13,14,15]) 

d[ columns ] = np.linalg.solve(K[ rows, columns ], f[ columns ]) 

# Calculate unknown f values 
f[ [8,9,10,11] ] = K[ [8,9,10,11], [8,9,10,11] ]*d[[8,9,10,11]]