用MATLAB模擬彈簧擺運(yùn)動(dòng)

1、用MATLAB模擬彈簧擺運(yùn)動(dòng)引言：物體在一定位置附近所作的往復(fù)運(yùn)動(dòng)稱為機(jī)械振動(dòng)。這種振動(dòng)現(xiàn)象在自然界是廣泛存在的。例如，擺的運(yùn)動(dòng)，一切發(fā)聲體的運(yùn)動(dòng)，機(jī)器開動(dòng)時(shí)各部分的微小顫動(dòng)等都是機(jī)械振動(dòng)。在不同的振動(dòng)現(xiàn)象中，最基本最簡單的振動(dòng)是簡諧振動(dòng)。一切復(fù)雜的振動(dòng)都可以分解為若干個(gè)簡諧振動(dòng)，彈簧振子和小角度單擺就是簡單的簡諧振動(dòng)。計(jì)算機(jī)模擬可以直觀表現(xiàn)物理系統(tǒng)的運(yùn)動(dòng)狀態(tài)，對(duì)于理解抽象的物理公式和研究系統(tǒng)的運(yùn)動(dòng)帶來很多方便，所以在科研中有廣泛的應(yīng)用。這里我們借助于MATLAB,在研究物理系統(tǒng)的運(yùn)動(dòng)的同時(shí)，也用動(dòng)畫模擬他們的運(yùn)動(dòng)。由于這些運(yùn)動(dòng)模擬式采用真實(shí)的從求解運(yùn)動(dòng)方程得到的數(shù)據(jù)。所以利用matlab的數(shù)

2、據(jù)擬合可以畫出精確的圖形。題目：設(shè)質(zhì)量為m的擺錘掛在勁度系數(shù)為k,原長為l0的輕彈簧上，彈簧的另一端懸掛于固定點(diǎn)O,系統(tǒng)靜止自然下垂時(shí)彈簧長度為，系統(tǒng)可在過O點(diǎn)的豎直平面內(nèi)自由擺動(dòng)，試研究擺錘的運(yùn)動(dòng)。原理：系統(tǒng)自由度為2。以O(shè)為極點(diǎn)，豎直向下的Ox軸為極軸，建立極坐標(biāo)系，如圖所示，r為質(zhì)點(diǎn)m到O點(diǎn)距離，為Ox軸與彈簧間的夾角，則系統(tǒng)的拉格朗日函數(shù)為L=T-W=12mdrdt2+r2ddt2-mgrcos+12kr-l02由拉格朗日方程可求出系統(tǒng)的運(yùn)動(dòng)微分方程為 上式未作小擺角近似，因此可用以研究彈簧擺的大擺角運(yùn)動(dòng)，但是難以求出解析解。由數(shù)值計(jì)算結(jié)果所畫出彈簧擺的大擺角運(yùn)動(dòng)軌跡，可發(fā)現(xiàn)它的運(yùn)動(dòng)情

3、況很復(fù)雜，不作數(shù)值計(jì)算是無法想象的。令，則成為 源程序：function thbglobal L m k gtheta0=pi/10; %初始角度,可由讀者設(shè)不同的值m=1;k=80;g=9.8; L0=1;L=L0+m*g/k; %L0為彈簧原來長度，L為彈簧靜止時(shí)長度t,u1=ode45(thbfun,0:0.005:15,L0 0 theta0 0);y1,x1=pol2cart(u1(:,3),u1(:,1);y1 = -y1;%將極坐標(biāo)換為直角坐標(biāo) figureymax=max(abs(y1);axis(-1.2 1.2 -1.2*ymax 0.2);%設(shè)置坐標(biāo)范圍axis offt

4、itle(彈簧擺,fontsize,14)hold on; R =0.055 ; %設(shè)置彈簧半徑y(tǒng)y = -L0 : 0.01 : 0;xx = R*sin(yy./L0*30*pi);%用正弦曲線表示彈簧a,r = cart2pol(xx,yy); %用坐標(biāo)變換來畫初始位置的彈簧a = a + theta0;xx,yy = pol2cart(a,r);%彈簧的數(shù)據(jù) line(-1 1,0 0,color,r,linewidth,2)%畫橫桿，球與彈簧ball = line(xx(1),yy(1),color,r,marker,.,markersize,70,erasemode,xor);ba

5、ll2 = line(xx(1),yy(1),color,0.5 0.51 0.6,linestyle,-,linewidth,1.3,erasemode,none);spring = line(xx,yy,color,r,linewidth,2,erasemode, xor);pause(0.5) for i = 1 : length(t) %畫出每一步的彈簧位置 yy = -u1(i,1) : 0.01 : 0; xx = R*sin(yy./u1(i,1)*30*pi); a,r = cart2pol(xx,yy); a = a + u1(i,3); xx,yy = pol2cart(a

6、,r); set(ball,XData,x1(i),YData,y1(i); set(ball2,XData,x1(i),YData,y1(i); set(spring,XData,xx,YData,yy); drawnow;end function F = thbfun(t,u)global L m k gF = u(2); u(1)*u(4)2 + g*cos(u(3) - k/m*(u(1) - L + m*g/k); u(4); -2*u(2)*u(4)/u(1) - g*sin(u(3)/u(1);程序說明：程序中解微分方程的過程比較簡單，主要的技巧是在做模擬動(dòng)畫上。有兩點(diǎn)值得注意，一是模擬彈簧的運(yùn)動(dòng)，另一點(diǎn)是在模擬彈簧運(yùn)動(dòng)的同時(shí)畫出軌跡圖。主程序的文件名是thb，子程序文件為thbfun，兩個(gè)程序共用的變量L，m，k，g用指令global來傳遞。程序運(yùn)行截圖：總結(jié)：這次小論文的編寫，雖然程序還是相對(duì)簡單的，但是體現(xiàn)了MATLAB真正