24小时热门版块排行榜    

查看: 1628  |  回复: 0

caoxinchd

新虫 (初入文坛)

[求助] 用Newmark方法计算系统的动力学响应的matlab程序

请大家帮忙看看这个程序有什么问题?用Newmark方法计算系统的动力学响应,结果大的惊人。
function[Q,V,AA]=newmarkb
E=2.1e11;P=7850;D1=0.405;d1=0.375;D2=0.375;d2=0.335;D3=0.335;d3=0.285;D4=0.285;d4=0.225;D5=0.225;d5=0.150;
A=(pi*(D1^2-d1^2))/4;
I=(pi*(D1^4-d1^4))/64;
M1= Mass (P,A,I,0,0,13,0);
A=(pi*(D2^2-d2^2))/4;
I=(pi*(D2^4-d2^4))/64;
M2= Mass (P,A,I,13,0,26,0);
A=(pi*(D3^2-d3^2))/4;
I=(pi*(D3^4-d3^4))/64;
M3= Mass (P,A,I,26,0,39,0);
A=(pi*(D4^2-d4^2))/4;
I=(pi*(D4^4-d4^4))/64;
M4= Mass (P,A,I,39,0,52,0);
A=(pi*(D5^2-d5^2))/4;
I=(pi*(D5^4-d5^4))/64;
M5= Mass (P,A,I,52,0,65,0);
M=zeros(18,18);
M= MAssemble(M,M1,1,2);
M= MAssemble(M,M2,2,3);
M= MAssemble(M,M3,3,4);
M= MAssemble(M,M4,4,5);
M= MAssemble(M,M5,5,6);% 整体质量矩阵
A=(pi*(D1^2-d1^2))/4;
I=(pi*(D1^4-d1^4))/64;
K1= LStiffness1 (E,A,I,0,0,13,0);
A=(pi*(D2^2-d2^2))/4;
I=(pi*(D2^4-d2^4))/64;
K2= LStiffness2 (E,A,I,13,0,26,0);
A=(pi*(D3^2-d3^2))/4;
I=(pi*(D3^4-d3^4))/64;
K3= LStiffness3 (E,A,I,26,0,39,0);
A=(pi*(D4^2-d4^2))/4;
I=(pi*(D4^4-d4^4))/64;
K4= LStiffness4 (E,A,I,39,0,52,0);
A=(pi*(D5^2-d5^2))/4;
I=(pi*(D5^4-d5^4))/64;
K5= LStiffness5 (E,A,I,52,0,65,0);
K=zeros(18,18);
K= KAssemble(K,K1,1,2);
K= KAssemble(K,K2,2,3);
K= KAssemble(K,K3,3,4);
K= KAssemble(K,K4,4,5);
K= KAssemble(K,K5,5,6);% 整体刚度矩阵
C=0.00776*M+0.00398*K;          %瑞利阻尼矩阵
%newmark系数
dt=0.01;
nt=20;%计算相应时间
betae=0.25;
alfa=0.5;
a0=1/(betae*dt^2);
a1=alfa/(betae*dt);
a2=1/(betae*dt);
a3=1/(2*betae)-1;
a4=alfa/betae-1;
a5=dt/2*(alfa/betae-2);
a6=dt*(1-alfa);
a7=dt*alfa;
%系数定义完
KE=K+a0*M+a1*C;
[L,U]=lu(KE);
q=zeros(18,1);  
v=zeros(18,1);  
pp=[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,-100,0,100].';
a=M^(-1)*(pp-K*q-C*v);  
t=0;
Q(:,1)=q;
V(:,1)=v;
AA(:,1)=a;
PP(:,1)=pp;
for i=1:nt-1
  PP(:,i+1)=PP(:,i)+M*(a0*Q(:,i)+a2*V(:,i)+a3*AA(:,i))+C*(a1*Q(:,i)+a4*V(:,i)+a5*AA(:,i));
  ik=L\PP(:,i+1);
  ik=U\ik;
  Q(:,i+1)=L'\ik;
  AA(:,i+1)=a0*(Q(:,i+1)-Q(:,i))-a2*V(:,i)-a3*AA(:,i);
  V(:,i+1)=V(:,i)+a6*AA(:,i)+a7*AA(:,i+1);
end
回复此楼

» 猜你喜欢

已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

智能机器人

Robot (super robot)

我们都爱小木虫

相关版块跳转 我要订阅楼主 caoxinchd 的主题更新
最具人气热帖推荐 [查看全部] 作者 回/看 最后发表
[考研] 299求调剂 +4 shxchem 2026-03-20 6/300 2026-03-21 11:54 by 无际的草原
[考研] 336求调剂 +3 rmc8866 2026-03-21 3/150 2026-03-21 11:16 by 求调剂zz
[考研] 268求调剂 +7 简单点0 2026-03-17 7/350 2026-03-21 08:32 by JourneyLucky
[考研] 求调剂 +6 Mqqqqqq 2026-03-19 6/300 2026-03-21 08:04 by JourneyLucky
[考研] 301求调剂 +10 yy要上岸呀 2026-03-17 10/500 2026-03-21 03:14 by JourneyLucky
[考研] 299求调剂 +6 △小透明* 2026-03-17 6/300 2026-03-21 02:42 by JourneyLucky
[考研] 328求调剂,英语六级551,有科研经历 +4 生物工程调剂 2026-03-17 8/400 2026-03-21 02:12 by JourneyLucky
[考研] 求调剂 +3 Ma_xt 2026-03-17 3/150 2026-03-21 02:05 by JourneyLucky
[考研] 085700资源与环境308求调剂 +12 墨墨漠 2026-03-18 13/650 2026-03-21 01:42 by JourneyLucky
[考研] 南京大学化学376求调剂 +3 hisfailed 2026-03-19 6/300 2026-03-20 23:43 by hisfailed
[考研] 材料与化工 322求调剂 +4 然11 2026-03-19 4/200 2026-03-20 22:12 by luoyongfeng
[考研] 290求调剂 +7 ^O^乜 2026-03-19 7/350 2026-03-20 21:43 by JourneyLucky
[考研] 260求调剂 +3 朱芷琳 2026-03-20 3/150 2026-03-20 20:35 by 学员8dgXkO
[考研] 一志愿吉林大学材料学硕321求调剂 +11 Ymlll 2026-03-18 15/750 2026-03-20 19:40 by 丁丁*
[考研] 085410人工智能专硕317求调剂(0854都可以) +4 xbxudjdn 2026-03-18 4/200 2026-03-20 09:07 by 不168
[考研] 材料专硕306英一数二 +10 z1z2z3879 2026-03-16 13/650 2026-03-18 14:20 by 007_lilei
[考研] 277调剂 +5 自由煎饼果子 2026-03-16 6/300 2026-03-17 19:26 by 李leezz
[考研] 085601求调剂 +4 Du.11 2026-03-16 4/200 2026-03-17 17:08 by ruiyingmiao
[考研] 277材料科学与工程080500求调剂 +3 自由煎饼果子 2026-03-16 3/150 2026-03-16 14:10 by 运气yunqi
[考研] 本科南京大学一志愿川大药学327 +3 麦田耕者 2026-03-14 3/150 2026-03-14 20:04 by 外星文明
信息提示
请填处理意见