声振论坛

 找回密码
 我要加入

QQ登录

只需一步,快速开始

查看: 1826|回复: 8

[线性振动] 求解一个非线性微分方程组

[复制链接]
发表于 2012-10-4 06:39 | 显示全部楼层 |阅读模式

马上注册,结交更多好友,享用更多功能,让你轻松玩转社区。

您需要 登录 才可以下载或查看,没有账号?我要加入

x
求解一个微分方程组 Mx''+Cx'+Kx=F
主程序算出了M C K 都为5×5的矩阵, F的表达式随时间t的变量,用ODE求解
主程序
[t,X]=ode45(@(t,X)thermal1(t,X,M,K,C,tt,F1,f1,f2,f3,f4),[0:120],[0 0 0 0 0 0 0 0 0 0]);
F1,f1,f2,f3,f4 为F的表达式
tt = linspace(0,120,120);
dT1=0;
dT2=0;
for ii=1:1:100
    dT1=dT1+((-1)^ii/ii^2)*(-ii^2*pi^2*k2/h^2)*exp(-ii^2*pi^2*k2*tt/h^2)*(cos(ii*pi+ii*pi*h1/2/h)-cos(-ii*pi*h1/2/h));
end
MT1=Efc*af2*(w*h1*h/2)*(qs*h/k1)*(-2/pi^2)*dT1;
F1=rho*A*f11*MT1/2/EI;
% test'error'
for ii=1:1:100
    dT2=dT2+((-1)^ii/ii^2)*(-ii^2*pi^2*k2/h^2)^2.*exp(-ii^2*pi^2*k2*tt/h^2)*(cos(ii*pi+ii*pi*h1/2/h)-cos(-ii*pi*h1/2/h));
end
MT2=Efc*af2*(w*h1*h/2)*(qs*h/k1)*(-2/pi^2)*dT2;
MTT2=-MT2/2/EI;
MTT1=-MT1/2/EI;
f1=-quadl(@(x)f_fs1(x,1),0,L)*MTT2*rho*A-quadl(@(x)f_fs1(x,1),0,L)*MTT1*Cdamp;
f2=-quadl(@(x)f_fs1(x,2),0,L)*MTT2*rho*A-quadl(@(x)f_fs1(x,2),0,L)*MTT1*Cdamp;
f3=-quadl(@(x)f_fs1(x,3),0,L)*MTT2*rho*A-quadl(@(x)f_fs1(x,3),0,L)*MTT1*Cdamp;
f4=-quadl(@(x)f_fs1(x,4),0,L)*MTT2*rho*A-quadl(@(x)f_fs1(x,4),0,L)*MTT1*Cdamp;
就是在主程序里面先将F求出来了,然后子程序
function dX=thermal1(t,X,M,K,C,tt,F1,f1,f2,f3,f4)
F1= interp1(tt,F1,t);
f1= interp1(tt,f1,t);
f2= interp1(tt,f2,t);
f3= interp1(tt,f3,t);
f4= interp1(tt,f4,t);
dX=[M*X(6:10); F-C*X(6:10)-K*X(1:5)];

这样子求可以不可以 为什么我求出的结果是错误的呢!
F=[F1;f1;f2;f3;f4;];

回复
分享到:

使用道具 举报

发表于 2012-10-4 18:04 | 显示全部楼层
你是否看过matlab 中的ode45例子,二阶常微分方程应该写成状态方程
发表于 2012-12-8 15:45 | 显示全部楼层
我想请问下,这种类型的方程,有很多,怎么转化为一阶微分方程求解
7.jpg
 楼主| 发表于 2012-12-9 11:14 | 显示全部楼层
看一下高数最后一章
发表于 2012-12-10 09:27 | 显示全部楼层
我是说用matlab编程解非线性方程组,mx''+kx=0,因为我的m为非对角阵,所以用微分方程解比较困难
 楼主| 发表于 2012-12-12 09:33 | 显示全部楼层
先把二阶的编程一阶的, 然后求解方程组, 但是我现在方程组怎么解都不对,化成向量就也不对
 楼主| 发表于 2012-12-12 09:34 | 显示全部楼层
我决跟你的M对焦不对焦关系不大吧 你的是几阶的呢
发表于 2013-1-3 17:13 | 显示全部楼层

参看线性系统的书。引入变量y=x的一阶导,写成状态向量形式就行!
发表于 2013-1-4 15:25 | 显示全部楼层
我的不能写成那样,因为列成微分方程后很难解!
您需要登录后才可以回帖 登录 | 我要加入

本版积分规则

QQ|小黑屋|Archiver|手机版|联系我们|声振论坛

GMT+8, 2024-5-9 23:36 , Processed in 0.076031 second(s), 21 queries , Gzip On.

Powered by Discuz! X3.4

Copyright © 2001-2021, Tencent Cloud.

快速回复 返回顶部 返回列表