我在解一个矩阵形式的微分方程,每次我把4个小矩阵连接成一个大矩阵,打印出来的结果太奇怪了,我都不知道对不对。
import numpy as np
from sympy import *
from sympy import init_printing
# from scipy.integrate import odeint
from scipy.integrate import solve_ivp
init_printing()
params = {'m1': 1, 'm2': 2, 'k1': 1, 'k2': 2, 'c1': 0.1, 'c2': 0.1}
params.update({'N' : 3})
params['dimM'] = 4*params['N']
##---------------------------------------## M1 et M2 ##--------------------------------------------------##
M1 = np.zeros([2*params['N'],2*params['N']])
print('M1 = \n', M1)
M2 = np.eye(2*params['N'])
print('M2 = \n', M2)
# m1 pour les indices pair
# m2 pour les indices impair
##---------------------------------------## M3 ##--------------------------------------------------------##
v3 = np.zeros(2*params['N'])
v3[::2]= -(params['k1']/params['m2'] - params['k2']/params['m2'])
v3[1::2]= -(params['k1']/params['m1'] - params['k2']/params['m1'])
print("Vecteur v3:", v3, "\n")
v3_b = np.zeros(2*params['N']-1)
v3_b[::2] = params['k2']/params['m1']
v3_b[1::2] = params['k1']/params['m2']
print("Vecteur v3_b:", v3_b, "\n")
M3 = np.diag(v3, 0) + np.diag(v3_b, +1) + np.diag(v3_b, -1)
print('M3 = \n', M3)
##---------------------------------------## M4 ##--------------------------------------------------------##
v4 = np.zeros(2*params['N'])
v4[::2]= -(params['c1']/params['m2'] - params['c2']/params['m2'])
v4[1::2]= -(params['c1']/params['m1'] - params['c2']/params['m1'])
print("Vecteur v4:", v4, "\n")
v4_b = np.zeros(2*params['N']-1)
v4_b[::2] = params['c2']/params['m1']
v4_b[1::2] = params['c1']/params['m2']
print("Vecteur v4_b:", v4_b, "\n")
M4 = np.diag(v4, 0) + np.diag(v4_b, +1) + np.diag(v4_b, -1)
print('M4 = \n', M4)
##---------------------------------------## Matrice complete M ##----------------------------------------##
M = [[M1, M2], [M3,M4]]
print('M = \n', M)
如果你执行这段代码,你就可以看到所有的矩阵。M1,M2,M3,M4都运行得很好。但是M显示了这个奇怪的混乱:
M =
[[array([[0., 0., 0., 0., 0., 0.],
[0., 0., 0., 0., 0., 0.],
[0., 0., 0., 0., 0., 0.],
[0., 0., 0., 0., 0., 0.],
[0., 0., 0., 0., 0., 0.],
[0., 0., 0., 0., 0., 0.]]), array([[1., 0., 0., 0., 0., 0.],
[0., 1., 0., 0., 0., 0.],
[0., 0., 1., 0., 0., 0.],
[0., 0., 0., 1., 0., 0.],
[0., 0., 0., 0., 1., 0.],
[0., 0., 0., 0., 0., 1.]])], [array([[0.5, 2. , 0. , 0. , 0. , 0. ],
[2. , 1. , 0.5, 0. , 0. , 0. ],
[0. , 0.5, 0.5, 2. , 0. , 0. ],
[0. , 0. , 2. , 1. , 0.5, 0. ],
[0. , 0. , 0. , 0.5, 0.5, 2. ],
[0. , 0. , 0. , 0. , 2. , 1. ]]), array([[0. , 0.1 , 0. , 0. , 0. , 0. ],
[0.1 , 0. , 0.05, 0. , 0. , 0. ],
[0. , 0.05, 0. , 0.1 , 0. , 0. ],
[0. , 0. , 0.1 , 0. , 0.05, 0. ],
[0. , 0. , 0. , 0.05, 0. , 0.1 ],
[0. , 0. , 0. , 0. , 0.1 , 0. ]])]]
2条答案
按热度按时间ldioqlga1#
经过简短的谷歌搜索,Numpy有一个连接功能。
因此,这个问题的答案很简单,只需使用numpy的concatenate函数,如下所示:
piv4azn72#
您可以使用
np.block
: