huazi071783 发表于 2011-5-12 22:09

用荣格-库塔法,如何实现外部力为环境激励问题

本帖最后由 huazi071783 于 2011-5-12 22:11 编辑

用荣格-库塔发解结构振动微分方程,求系统相应x,Mx''+cx'+kx=F(t),我设激振力F(t)为环境激励,就是用randn生成的随机数组。不知道能不能以我下面编的程序实现,就是逐个的解,虽然能计算,但计算速度特慢。如果激励为函数形势的如正弦函数就特别快。还有一下红字部分问题。下面是我的程序
clear;clc
global m k c f
m=0.5;
c=0.02;
k=80;
f=5;
x0=;
tspan=;
=ode45(@randvibration,tspan,x0);
figure(1)
set(gca,'FontSize',30)
plot(y(2000:10000,1),y(2000:10000,2),'r')
调用函数
function dy=randvibration(t,y)
global m k c f
load('randsignal.mat');
F=randsignal;
n=floor(t./0.01);
%n=t./0.01;这样会出错,提示n到一定大的时候就不是整数了,按理说不会出现小数的。为什么?
% dy= 这个计算特别块
dy(2)=-(c./m)*y(2)-(k./m)*y(1)+f*F(n+1);
dy(1)=y(2);
dy=dy';

结果不知道是否正确,请高手指点

huazi071783 发表于 2011-5-12 22:46

这是算出的结果画出的相轨迹图

meiyongyuandeze 发表于 2011-5-13 08:23

本帖最后由 meiyongyuandeze 于 2011-5-13 08:25 编辑

不知道楼主程序里面为什么都是用的点除呢?程序运行慢是因为没此迭代都要读取数据的原因,可以尝试在主程序中先.mat数据后用参数传递的方法进入状态方程,这样可能会节省一定的计算时间。没有你的.mat文件,所以也没法运行,随便说说看法!

huazi071783 发表于 2011-5-13 10:15

回复 3 # meiyongyuandeze 的帖子

请问怎么用参数传递呢?没有明白你的意思?指教!

huazi071783 发表于 2011-5-13 14:16

回复 3 # meiyongyuandeze 的帖子

谢谢,解决了,果然是你说的那样,改完后计算非常快了

meiyongyuandeze 发表于 2011-5-13 14:30

回复 5 # huazi071783 的帖子

呵呵,客气,一起学习!
页: [1]
查看完整版本: 用荣格-库塔法,如何实现外部力为环境激励问题