看了一篇柱塞幫浦離散化動力學建模的文章，感覺還挺有意思，於是嘗試做一下
二、matlab下的動力學方程總結
斜盤式軸向柱塞幫浦是一類常見的柱塞幫浦，本文以ｐｃｙ－２５型斜盤式軸向柱塞幫浦為研究物件，研究幫浦內機械振動的傳遞問題。由於該幫浦傳動軸與缸體之間為過盈配合，且柱塞滑靴元件位於缸體的柱塞腔內，因此，將傳動軸、缸體和柱塞滑靴元件視為乙個剛性體，該旋轉體即為轉子系統。
１．轉子系統是機械振動產生的主振源。轉子系統轉動時，傳動軸、缸體 和柱 塞滑靴元件的自激振
動，以及這些零部件與相接觸部件之間的相互作用，都是機械振動的激勵源。
２．幫浦殼是機械振動最終受體。該幫浦採用高強螺栓將前殼體與鐘形罩緊固連線，如認為幫浦的安裝為
完全固支，幫浦殼則是最終振動受體，但是三個殼體受到的振動會相互疊加作用。
３．後殼體振動最為複雜且劇烈。由安裝方式可以看出，整幫浦振動結構為懸臂梁結構，其主振型為垂
直於軸向的上下左右擺動，此外，其他兩個殼體承受的振動也會疊加作用到後殼體上。
４．該振 動 垂 直 於 幫浦 軸 由 內 向 外 傳 遞。幫浦 工 作時，旋轉元件圍繞幫浦軸高速轉動，在離心力及主振型作用下，振動將沿著幫浦軸的方向由內而外傳遞。
５．振動路徑貢獻度取決於接觸零部件之間的徑向等效剛度和阻尼。
各引數數值
該模型的振動傳遞路徑拉格朗日微分方程如下：
其中，激振力 ｆ ＝f0*sin(πｎ/60)*ｔ，ｎ為軸向柱塞幫浦轉速，這裡設定的n=1460r/min，當幫浦的轉速不同時，激振力頻率也會不同。由於軸向柱塞幫浦激振力和軸承的剛度與阻尼是時變的，因此，該模型屬於非線性時變方程組，因此，在理論上不存在解析解，只能得到數值解。
求解這個動力學方程使用的是變步長求解器ode45。主要的核心點就是這個6元的二階微分方程的方程建立，讓計算機能夠讀懂你的方程。在這個過程中需要使用到換元，來代替中間的一階項。
function dx=
test_fun
(t,x)
mtotal=
6.8;
mvp =
0.25
;mf =
8.1;
mm =
8.6;
mb =
8.5;
msp =
2.1;
kfe1 =
2.5*10^
8;kfe2 =3*
10^8;
kb1 =
5.11*10
^8+1.3*(
10^7)
*sin
(438
*pi*t)
;kb2 =
6.04*10
^8+1.55*(
10^7)
*sin
(1216
*pi*t)
;k1=
5.2*10^
8;k2=6.3*10
^8;kvp=
4.71*10
^8;kvm=10^
9;f0=250;n=
1460
;cvp=
1100
;csb=
370;
cb1=
(1900
+100
*sin
(438
*pi*t));
cb2=
(2500
+500
*sin
(1216
*pi*t));
dx=zeros(12
,1);
dx(1)
=x(7
);dx(
2)=x
(8);
dx(3)
=x(9
);dx(
4)=x
(10);
dx(5)
=x(11
);dx(
6)=x
(12);
dx(7)
=(-cb1*(dx
(1)-
dx(3)
)-cb2*(dx
(1)-
dx(4)
)-cvp*(dx
(1)-
dx(2)
)-csb*(dx
(1)-
dx(6)
)-kfe1*x(
1)-kb1*(x
(1)-
x(3)
)-kb2*(x
(1)-
x(4)
)+f0*
sin((2
*pi*n/60)
*t))
/mtotal;dx(
8)=(cvp*(dx
(1)-
dx(2)
)-kvp*(x
(2)-
x(4)
))/mvp;dx(
9)=(cb1*(dx
(1)-
dx(3)
)-kfe2*x(
3)+kb1*(x
(1)-
x(3)
)-k1*(x
(3)-
x(4)
))/mf;dx(
10)=(cb2*(dx
(1)-
dx(4)
)+kb2*(x
(1)-
x(4)
)+kvp*(x
(2)-
x(4)
)+k1*(x
(3)-
x(4)
)-k2*(x
(4)-
x(5)
))/mm;dx(
11)=(k2*(x
(4)-
x(5)
)-kvm*(x
(5)-
x(6)
))/mb;dx(
12)=(csb*(dx
(1)-
dx(6)
)+kvm*(x
(5)-
x(6)
))/msp;
end
[t,x]
=ode45
(@test_fun,[0
4],[
0000
0011
1111]);
利用Matlab求解Laplace方程
問題描述 求解泊松方程 u 1 求解區域為單位圓盤，邊界條件為在圓盤邊界上u 0 下面求它的數值解，編寫程式如下 1 問題定義 g circleg 單位圓 b circleb1 邊界上為零條件 c 1 a 0 f 1 2 產生初始的三角形網格 p,e,t initmesh g 3 迭代直至得到誤差允...
matlab求解平面方程的原理
已知三點p1,p2,p3,求其平面方程 syms x y z p1,p2,p3的座標由自己定義。p1 x1,y1,z1 p2 x2,y2,z2 p3 x3,y3,z3 那麼求解下面矩陣q行列式就是了 q ones 4,1 x,y,z p1 p2 p3 detb det q 最後令 q 0 這裡的求解...
MATLAB 求解方程（組）
eg.解方程x 2 x 2 0 1.roots p 函式 此 matlab 函式 以列向量的形式返回 p 表示的多項式的根。輸入 p 是乙個包含 n 1 多項式係數的向量，以 xn 係數開頭。0係數表示方程中不存在的中間冪。p 1 1,2 x roots p 2.solve函式 利用solve函式求...