24小时热门版块排行榜    

查看: 1875  |  回复: 4

明明那天

新虫 (初入文坛)

[求助] matlab 数值积分速度上不去

我是新手,在用matlab计算数值积分的程序,可能是存在奇点问题,运行速度非常慢。不知有没有高手可以帮忙指导一下。谢谢了!

clear all   
clc     
syms x y u1 u2 k1 k2 r1 r2 ss1 ss2 a1 c1 s1 cc1 c_c1
A11=1.58e8; A66=1.05e7;A22=7.99e7;

A12=5.625e5;
D1=A66*A22/(A12^2-A11*A22);
D2=A66*A11/(A12^2-A11*A22);
D3=-2*A66*A12/(A12^2-A11*A22)+1;
u=[1.875;4.694;7.855;10.996;14.173;17.279];%x
q=[4.730;7.853;10.996;14.137;17.279;20.378];%y
FUN5=0;


X=[];Y=[];int1=zeros(36,36);int2=zeros(36,36);int3=zeros(36,36);row=0;
for mode_num2=1:6  
      u_2=u(mode_num2);k_2=q(mode_num2);   
    r_2=(cos(u_2)-cosh(u_2))/(sinh(u_2)-sin(u_2));  
     s_2=(sinh(k_2)-sin(k_2))/(cos(k_2)+cosh(k_2));
    X=[X,cosh(u_2*x)+cos(u_2*x)+r_2*(sin(u_2*x)+sinh(u_2*x))];
    Y=[Y,cosh(k_2*y)-cos(k_2*y)+s_2*(sin(k_2*y)-sinh(k_2*y))];
end
for mode_num=1:6
       W5=0;F5=0;L=0;row=row+1;
       u_1=u(mode_num);k_1=q(mode_num);   
       r_1=(cos(u_1)-cosh(u_1))/(sinh(u_1)-sin(u_1));  
       s_1=(sinh(k_1)-sin(k_1))/(cos(k_1)+cosh(k_1));
  X5=[sinh(u1*x)^2*u1^2,-2*sinh(u1*x)*u1^2*sin(u1*x),2*sinh(u1*x)*u1^2*r1*cos(u1*x),2*sinh(u1*x)*u1^2*r1*cosh(u1*x),sin(u1*x)^2*u1^2,-2*sin(u1*x)*u1^2*r1*cos(u1*x),-2*sin(u1*x)*u1^2*r1*cosh(u1*x),r1^2*cos(u1*x)^2*u1^2,2*r1^2*cos(u1*x)*u1^2*cosh(u1*x),r1^2*cosh(u1*x)^2*u1^2];
  Y5=[sinh(k1*y)^2*k1^2,2*sinh(k1*y)*k1^2*sin(k1*y),2*sinh(k1*y)*k1^2*ss1*cos(k1*y),-2*sinh(k1*y)*k1^2*ss1*cosh(k1*y),sin(k1*y)^2*k1^2,2*sin(k1*y)*k1^2*ss1*cos(k1*y),-2*sin(k1*y)*k1^2*ss1*cosh(k1*y),ss1^2*cos(k1*y)^2*k1^2,-2*ss1^2*cos(k1*y)*k1^2*cosh(k1*y),ss1^2*cosh(k1*y)^2*k1^2];

  for i=1:10
    for  j=1:10
      W5=W5-u1^2*k1^2*X5(i)*Y5(j);
    end
  end


  syms c_c1 c_c101 c_c201 c_c301 c_c401 c_c501
c_c=[c_c1, c_c101, c_c201, c_c301, c_c401, c_c501];
   if mode_num==1
      str5=c_c1;
        for i=2:100
    syms(['c_c',num2str(i)])
    str5=[str5,['c_c',num2str(i)]];
        end   

   else
     i_1=2+(mode_num-1)*100;i_2=100+(mode_num-1)*100;
     str5=c_c(mode_num);
    for  i=i_1:i_2
    syms(['c_c',num2str(i)])
    str5=[str5,['c_c',num2str(i)]];
    end
   end
%str5
x_5=[cosh(u1*x),cos(u1*x),sin(u1*x),sinh(u1*x)];
  y_5=[cosh(k1*y),cos(k1*y),sin(k1*y),sinh(k1*y)];
  xm1=[cosh(u_1*x),cos(u_1*x),sin(u_1*x),sinh(u_1*x)];
  ym1=[cosh(k_1*y),cos(k_1*y),sin(k_1*y),sinh(k_1*y)];
  intmode1=zeros(100,36);intmode2=zeros(100,36);
  intmode3=zeros(100,36);
  for i=1:4
      for j=i:4  
          for r=1:4
              for k=r:4
          L=L+1;
          F5=F5+str5(L)*x_5(i)*x_5(j)*y_5(r)*y_5(k);
      for mu=1:6
          for s=1:6
          m=s+(mu-1)*6;
         
Ax1=diff(X(mu),x,2);
Ay1=Y(mu);  
Fx1=xm1(i)*xm1(j);
Fy1=diff(ym1(r)*ym1(k),y,2);
Sx1=vectorize(Fx1*X(s)*Ax1);
Sy1=vectorize(Fy1*Y(s)*Ay1);
intx1=quadl(Sx1,0,1) ;                                    
inty1=quadl(Sy1,0,1) ;
intmode1(L,m)=intx1*inty1;
         
Ax2=diff(X(mu),x);
Ay2=diff(Y(mu),y);  
Fx2=diff(xm1(i)*xm1(j),x);
Fy2=diff(ym1(r)*ym1(k),y);
Sx2=vectorize(Fx2*X(s)*Ax2);
Sy2=vectorize(Fy2*Y(s)*Ay2);
intx2=quadl(Sx2,0,1) ;                                    
inty2=quadl(Sy2,0,1) ;
intmode2(L,m)=intx2*inty2;

Ax3=X(mu);
Ay3=diff(Y(mu),y,2);  
Fx3=diff(xm1(i)*xm1(j),x,2);
Fy3=ym1(r)*ym1(k);
Sx3=vectorize(Fx3*X(s)*Ax3);
Sy3=vectorize(Fy3*Y(s)*Ay3);
intx3=quadl(Sx3,0,1) ;                                    
inty3=quadl(Sy3,0,1) ;
intmode3(L,m)=intx3*inty3;
          end
      end         
              end
          end
      end
  end
回复此楼

» 猜你喜欢

» 本主题相关价值贴推荐,对您同样有帮助:

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

Honbein

木虫 (正式写手)

【答案】应助回帖

★
感谢参与,应助指数 +1
jjdg: 金币+1, 感谢参与 2013-02-12 16:24:58
嵌套的for循环太多,多用向量表示能提高速度。MATLAB做循环是硬伤…

[ 发自手机版 http://muchong.com/3g ]
2楼2013-02-12 09:16:26
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

weiwei43

至尊木虫 (正式写手)

★
jjdg: 金币+1, 一针见血 2013-02-12 16:24:39
编程快是Matlab的优点,计算速度是他的软肋。要想速度提高,还是用别的语言吧。
3楼2013-02-12 10:14:19
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

nono2009

超级版主 (文学泰斗)

No gains, no pains.

优秀区长优秀区长优秀区长优秀区长优秀版主

【答案】应助回帖

感谢参与,应助指数 +1
xzhdty: 谢谢参与 2013-02-17 09:55:33
用向量计算代替循环计算;
尽量多用matlab提供的函数,包括数值积分。
4楼2013-02-12 19:42:18
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

lv123520

银虫 (初入文坛)

【答案】应助回帖

上面说的都有道理,我也用matlab,速度也是上不去的,但是它编程简单可以用进行符号运算。但是你要是求计算速度,还是fortran最快。如果你对速度又要求还是用别的语言好点。好运。我水平有限。
5楼2013-03-06 15:23:37
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
相关版块跳转 我要订阅楼主 明明那天 的主题更新
最具人气热帖推荐 [查看全部] 作者 回/看 最后发表
[教师之家] 某top大学教授说“能够在市场中兑现的能力才是真能力”无比同意! +3 zju2000 2026-09-26 3/150 2026-09-27 16:26 by lala1008
[硕博家园] 上海工程技术大学 激光智能制造课题组 2027级博士研究生招生 +5 oxidpl 2026-09-20 7/350 2026-09-27 10:43 by Equinoxhua
[博后之家] 售SCI文章,我:8O.5.5.1O.54,科目全,可十急 +3 ViqFlLSxHDGt 2026-09-26 3/150 2026-09-27 10:42 by 4IW0sJvtEqX8
[找工作] 售SCI一区T0P文章,我:8O.55.1.O.54,科目全,可伽急 +3 ViqFlLSxHDGt 2026-09-26 3/150 2026-09-27 10:26 by 4IW0sJvtEqX8
[博后之家] 售SCI一区T0P文章,我:8O.55.1.O.5.4,科目齐全,可+急 +3 ViqFlLSxHDGt 2026-09-26 3/150 2026-09-27 10:21 by 4IW0sJvtEqX8
[公派出国] 售SCI一区T0P文章,我:8.O.55.1.O54,科目全,可伽急 +3 ViqFlLSxHDGt 2026-09-26 3/150 2026-09-27 10:15 by 4IW0sJvtEqX8
[考博] 售SCI一区T0P文章,我:8.O.55.1.O.54,科目齐全,可+急 +3 oeyPlfmMUpOC 2026-09-26 3/150 2026-09-27 09:13 by 4IW0sJvtEqX8
[考博] 售SCI一区T0P文章,我:8.O.55.1.O.5.4,科目全,可+急 +3 GxESXJmptOtQ 2026-09-25 4/200 2026-09-27 06:42 by XgCC8uTcwILl
[找工作] 售SCI一区T0P文章,我:8.O.55.1.O54,科目全,可伽急 +3 GxESXJmptOtQ 2026-09-25 3/150 2026-09-27 06:32 by XgCC8uTcwILl
[找工作] 售SCI一区文章,我:8.O.55.1.O.54,科目齐全,可伽急 +3 vQZoDrrm7mUF 2026-09-25 3/150 2026-09-27 06:08 by XgCC8uTcwILl
[考研] 售一区SCI文章T0P,我:8O.551.O54,科目全,可十急 +3 vQZoDrrm7mUF 2026-09-25 3/150 2026-09-27 05:56 by XgCC8uTcwILl
[找工作] 售SCI一区文章,我:8.O.551.O.5.4,科目全,可伽急 +5 vQZoDrrm7mUF 2026-09-25 5/250 2026-09-27 05:53 by XgCC8uTcwILl
[博后之家] 售SCI一区T0P文章,我:8.O55.1.O.54,科目全,可十急 +4 vQZoDrrm7mUF 2026-09-25 4/200 2026-09-27 05:50 by XgCC8uTcwILl
[考研] 售SCI一区T0P文章,我:8O.55.1.O.54,科目全,可伽急 +3 vQZoDrrm7mUF 2026-09-25 5/250 2026-09-27 05:44 by XgCC8uTcwILl
[公派出国] 售SCI文章,我:8O5.5.1.O.54,科目齐全,可+急 +4 vQZoDrrm7mUF 2026-09-25 4/200 2026-09-27 05:44 by XgCC8uTcwILl
[考博] 售一区SCI文章T0P,我:8O.551.O54,科目全,可十急 +3 vQZoDrrm7mUF 2026-09-25 3/150 2026-09-27 05:44 by XgCC8uTcwILl
[论文投稿] 售SCI一区T0P文章,我:8O.55.1.O.54,科目全,可伽急 +3 vQZoDrrm7mUF 2026-09-25 6/300 2026-09-27 05:42 by XgCC8uTcwILl
[找工作] 售SCI一区T0P文章,我:8O.55.1.O.5.4,科目齐全,可+急 +4 PGTSwSC3F6nQ 2026-09-24 5/250 2026-09-25 11:19 by yhq9807
[公派出国] 售SCI-T0P文章,我:8O.5.5.1.O.54,科目齐全,可+急 +4 PGTSwSC3F6nQ 2026-09-24 7/350 2026-09-25 11:18 by yhq9807
[教师之家] 课题组招2027级博士 上海工程技术大学 激光智能制造方向 +4 oxidpl 2026-09-20 6/300 2026-09-24 20:37 by oxidpl
信息提示
请填处理意见