4.17如圖4.17所示，均質(zhì)矩形薄板質(zhì)量為m，長(zhǎng)為a，寬為b，可繞豎直邊轉(zhuǎn)動(dòng)，初始角速度為ω0。轉(zhuǎn)動(dòng)時(shí)受到空氣阻力，其方向垂直于板面，每一小面積所受的阻力的大小與其面積和速度的n(n≥0)次方的乘積成正比，比例系數(shù)為kn(kn>0)。求薄板轉(zhuǎn)動(dòng)的角速度和角度。當(dāng)n取不同的正數(shù)時(shí)，角速度和角度與時(shí)間的關(guān)系曲線有什么特點(diǎn)？[解析]在薄板上取一長(zhǎng)為dr，寬為b的矩形，其面積為ds=bdr；矩形到軸的距離為r，當(dāng)薄板以角速度ω轉(zhuǎn)動(dòng)時(shí)，其速度為v=rω。矩形受到的阻力大小為df=knvnds=knbωnrndr，取角速度ω的方向?yàn)檎较?，矩形所受的阻力矩為dM=-rdf=-knbωnrn+1dr，負(fù)號(hào)表示力矩的方向與角速度的方向相反。薄板受到的總力矩為。(4.1)薄板繞軸的轉(zhuǎn)動(dòng)慣量為J=ma2/3，根據(jù)轉(zhuǎn)動(dòng)定理M=Jd2θ/dt2可列方程題4.17圖題4.17圖rωbadr即，(4.2)由于dθ/dt=ω，分離變量可得，(4.3)當(dāng)n=1時(shí)，微分方程簡(jiǎn)化為，積分可得，即。(4.4)當(dāng)n≠1時(shí)，(4.3)式積分可得，角速度為。設(shè)，角速度可表示為。(4.5)角度的微分為，設(shè)，則，因此，(4.6)當(dāng)n=2時(shí)，上式可化為，積分可得。當(dāng)t=0時(shí)，y=1時(shí)，θ=0，所以C=0。因此。(4.7)當(dāng)n≠2時(shí)，(4.6)式積分可得，因此角度為。(4.8)[討論]①設(shè)，當(dāng)n→1時(shí)，則ε→0。根據(jù)定義，由(4.5)式可得，這就是(4.4)式?？梢?jiàn)：(4.5)式包含(4.4)式。②設(shè)，當(dāng)n→2時(shí)，根據(jù)羅必塔法則，當(dāng)ε→0時(shí)，。取ε=n-2，，由(4.8)式可得，這就是(4.7)式?？梢?jiàn)：(4.8)式包含(4.7)式。③如果n=0，阻力是常數(shù)，角速度為，可知：薄板作勻變速直線運(yùn)動(dòng)。薄板旋轉(zhuǎn)的角度為。當(dāng)t=2t0/3時(shí)，薄板就會(huì)靜止，旋轉(zhuǎn)的角度為ω0t0/3。④如果0<n<1，當(dāng)薄板靜止時(shí)，經(jīng)過(guò)的時(shí)間為，薄板旋轉(zhuǎn)的角度為。⑤如果0<n<2，當(dāng)t→∞時(shí)薄板也會(huì)靜止，薄板旋轉(zhuǎn)的角度由上式?jīng)Q定。⑥如果2≤n，當(dāng)t→∞時(shí)薄板雖然會(huì)靜止，薄板旋轉(zhuǎn)的角度隨時(shí)間不斷增加。[算法]方法一：用解析式。取tn為時(shí)間單位，ω0為角速度單位，角速度可表示為ω=ω0exp(-t*)，(n=1)；(4.4*)，(n≠1)。(4.5*)其中，t*=t/tn。取θn=ω0tn為角度單位，角度可表示為，(n=2)；(4.7*)，(n≠2)。(4.8*)取指數(shù)n向量，再取時(shí)間向量，化為矩陣即可求角速度和角度。當(dāng)n為1或2時(shí)，需要加一個(gè)小量才能進(jìn)行極限運(yùn)算。(4.4*)式和(4.7*)式可驗(yàn)證極限運(yùn)算的正確性。[程序]P4_17.m如下。%薄板所受的阻力與速率的n次方成正比的運(yùn)動(dòng)(用解析式)clear%清除變量tm=5;%最大時(shí)間t=0:0.01:tm;%時(shí)間向量n=0:0.5:3;%指數(shù)向量n(n==1)=1+eps^(1/2);%為1的指數(shù)加小量n(n==2)=2+eps^(1/2);%為2的指數(shù)也加小量m=length(n);%指數(shù)個(gè)數(shù)[N,T]=meshgrid(n,t);%化為矩陣OMEGA=(1+3*(N-1)./(N+2).*T).^(1./(1-N));%求角速度THETA=(N+2)/3./(N-2).*((1+3*(N-1)./(N+2).*T).^((N-2)./(N-1))-1);%求角度OMEGA(N<1&T>(N+2)/3./(1-N))=nan;%將不合理的角速度改為非數(shù)THETA(N<1&T>(N+2)/3./(1-N))=nan;%將不合理的角度改為非數(shù)%----------------------------------------------------------figure%創(chuàng)建圖形窗口plot(t,OMEGA)%畫角速度曲線族gridon%加網(wǎng)格fs=16;%字體大小title('薄板所受阻力與速率的\itn\rm次方成正比的角速度與時(shí)間','fontsize',fs)%顯示標(biāo)題xlabel('時(shí)間\itt/t_n','fontsize',fs)%顯示橫坐標(biāo)ylabel('角速度\it\omega/\omega\rm_0','fontsize',fs)%顯示縱坐標(biāo)legend([repmat('\itn\rm=',m,1),num2str(n')])%圖例holdon%保持圖像plot(t,exp(-t),'.')%驗(yàn)證n=1時(shí)的角速度txt='\itt_n\rm=\itm/k_na^nb\omega\rm_0\it^n\rm^-^1';%時(shí)間文本text(0,0.1,txt,'fontsize',fs)%標(biāo)記時(shí)間文本figure%創(chuàng)建圖形窗口plot(t,THETA)%畫角度曲線族gridon%加網(wǎng)格title('薄板所受阻力與速率的\itn\rm次方成正比的角度與時(shí)間','fontsize',fs)%顯示標(biāo)題xlabel('時(shí)間\itt/t_n','fontsize',fs)%顯示橫坐標(biāo)ylabel('角度\it\theta/\omega\rm_0\itt_n','fontsize',fs)%顯示縱坐標(biāo)legend([repmat('\itn\rm=',m,1),num2str(n')],2)%圖例holdon%保持圖像plot(t,4/3*log(1+3/4*t),'.')%驗(yàn)證n=2時(shí)的角度[圖示](1)薄板轉(zhuǎn)動(dòng)的角速度如P4_17a圖所示，角速度隨時(shí)間的增加而減小。如果n=0，阻力是常數(shù)，薄板作勻變速直線運(yùn)動(dòng)，當(dāng)t=2t0/3時(shí)，薄板就會(huì)靜止。如果n=0.5，當(dāng)t=5t0.5/3時(shí)，薄板才會(huì)靜止。如果n=1，角速度按指數(shù)規(guī)律減少。(2)薄板轉(zhuǎn)動(dòng)的角度如P4_17b圖所示，角度隨時(shí)間增加而增加。如果n=0，薄板旋轉(zhuǎn)的最大角度為θ=ω0t0/3。如果n=0.5，薄板旋轉(zhuǎn)的最大角度為θ=5ω0t0.5/18。如果n=2，角度按對(duì)數(shù)規(guī)律增加。P4_17a圖P4_17b圖方法二：用一個(gè)一階常微分方程的數(shù)值解求解。利用時(shí)間單位，薄板的運(yùn)動(dòng)方程(4.2)式可化為，取約化時(shí)間為t*=t/tn，取約化角速度ω*=ω/ω0，可得，(4.2*)初始條件為ω*=1。薄板旋轉(zhuǎn)的角度為，求出角速度的數(shù)值解，通過(guò)數(shù)值積分就能求出薄板旋轉(zhuǎn)的角度，角度的單位是ω0tn。[程序]P4_17_.m計(jì)算部分如下。%薄板所受的阻力與速率的n次方成正比的運(yùn)動(dòng)(求一階微分方程的數(shù)值解)clear%清除變量tm=5;%最大時(shí)間dt=0.01;%時(shí)間間隔t=0:dt:tm;%時(shí)間向量n=0:0.5:3;%指數(shù)向量n(n==1)=1+eps^(1/2);%為1的指數(shù)加小量n(n==2)=2+eps^(1/2);%為2的指數(shù)也加小量m=length(n);%指數(shù)個(gè)數(shù)OMEGA=[];%角速度矩陣置空THETA=[];%角度矩陣置空f(shuō)ori=1:m%按指數(shù)循環(huán)s=num2str(n(i));%取指數(shù)并化為字符串f=inline(['-3/(2+'s')*w.^'s],'t','w');%被積內(nèi)線函數(shù)[tt,w]=ode45(f,t,1);%求微分方程數(shù)值解w(w<0)=nan;%負(fù)數(shù)改為非數(shù)w(imag(w)~=0)=nan;%復(fù)數(shù)改為非數(shù)th=cumsum(w)*dt;%累積角度OMEGA=[OMEGA,w];%連接角速度矩陣THETA=[THETA,th];%連接角度矩陣end%結(jié)束循環(huán)%----------------------------------------------------------(其他指令與上一程序的相應(yīng)部分是相同的。)[說(shuō)明]用常微分方程指令求微分方程的數(shù)值解特別簡(jiǎn)單。方法三：用二個(gè)一階常微分方程的數(shù)值解求解。角速度公式dθ/dt=ω可化為，取θn=ω0tn為角度單位，則得。其中，θ*=θ/θn，t*=t/tn，ω*=ω/ω0。設(shè)θ(1)=θ*，θ(2)=dθ*/dt*，可得，。當(dāng)t=0時(shí)，θ=0，ω=ω0。所以初始條件可表示為θ(1)=0，θ(2)=1。[程序]P4_17__.m計(jì)算部分如下。%薄板所受的阻力與速率的n次方成正比的運(yùn)動(dòng)(求二階微分方程數(shù)值解)clear%清除變量tm=5;%最大時(shí)間dt=0.01;%時(shí)間間隔t=0:dt:tm;%時(shí)間向量n=0:0.5:3;%指數(shù)向量n(n==1)=1+eps^(1/2);%為1的指數(shù)加小量n(n==2)=2+eps^(1/2);%為2的指數(shù)也加小量m=length(n);%指數(shù)個(gè)數(shù)OMEGA=[];%角速度矩陣置空THETA=[];%角度矩陣置空f(shuō)ori=1:m%按指數(shù)循環(huán)[tt,W]=ode45('P4_17__fun',t,[0,1],[],n(i));%求微分方程數(shù)值解th=W(:,1);%取角度向量w=W(:,2);%取角速度向量th(w<0|imag(w)~=0)=nan;%角速度小于零或者是復(fù)數(shù)時(shí)角度改為非數(shù)w(w<0|imag(w)~=0)=nan;%角速度小于零或者是復(fù)數(shù)時(shí)角速度改為非數(shù)OMEGA=[OMEGA,w];%連接角速度矩陣THETA=[THETA,th];%連接角度矩陣end%結(jié)束循環(huán)%----------------------------------------------------------(其他指令與上面兩個(gè)程序的相應(yīng)部分是相同的。)程序調(diào)用P4_17__fun.m求微分方程的數(shù)值解。%薄板所受的阻力與速率的n次方成正比的運(yùn)動(dòng)的函數(shù)functionf=fun(t,th,flag,n)f=[th(2);%速度-3/(n+2)*th(2)^n];