怎样用Matlab求复杂符号表达式的数值积分

2019-07-17 13:56发布

那位高手能帮我看看这个程序哪里出错了,最后一步积分积不出来,能帮我修改一下吗,非常感谢!
clc;
clear;
a=4.15.*(10).^-6;
l=1.6.*(10).^-6;
b=62.5.*(10).^-6;
n1=1.4504;
n2=1.4447;
n3=1.4786;
k=2.*pi./l;
V=k.*a.*((n1.^2-n2.^2).^(1./2));
cn=((n1.*k).^2-(2.27645.*V-0.306.*V.^2-0.992)./a.^2).^(1./2);
r=(cn.^2-n2.^2.*k.^2).^(1./2);
ka=(n1.^2.*k.^2-cn.^2).^(1./2);
ev=2;
R=8;
syms y;
X20=(R.*10.^-3./(2.*k.^2.*n2.^2)).^(2./3).*(cn.^2+y.^2-k.^2.*n2.^2);
X2a=(R.*10.^-3./(2.*k.^2.*n2.^2)).^(2./3).*(cn.^2+y.^2-k.^2.*n2.^2.*(1+2.*a./(R.*10.^-3)));
X2b=(R.*10.^-3./(2.*k.^2.*n2.^2)).^(2./3).*(cn.^2+y.^2-k.^2.*n2.^2.*(1+2.*b./(R.*10.^-3)));
X3b=(R.*10.^-3./(2.*k.^2.*n3.^2)).^(2./3).*(cn.^2+y.^2-k.^2.*n3.^2.*(1+2.*b./(R.*10.^-3)));
X2=(2.*k.^2.*n2.^2./(R.*10.^-3)).^(2./3).*(-X2b);
X3=(2.*k.^2.*n2.^2./(R.*10.^-3)).^(2./3).*(-X3b);
Oy=2./3.*(-X2b).^(3./2)+pi./4;
f=exp(-a.*(r.^2+y.^2).^(1/2))./(r.^2+y.^2).^(1/2).*1./pi.*(X20./3).^(1./2).*besseli(1./3,2./3.*X20.^(3./2))./(X2a./3).^(1./2)./(besselk(-1./3,2./3.*X2a.^(3./2))+besselk(1./3,2./3.*X2a.^(3./2))).*(X2).^(1/2).*(X3).^(1/2)./(X2.*(cos(Oy)).^2+X3.*(sin(Oy)).^2);
J=quadl(@(t)subs(f,y,t),0,inf)
友情提示: 此问题已得到解决,问题已经关闭,关闭后问题禁止继续编辑,回答。
该问题目前已经被作者或者管理员关闭, 无法添加新回复
4条回答
913688247
1楼-- · 2019-07-17 17:41
建议认真理清自己的思路   这样 看上去 你的代码太乱了。
z00
2楼-- · 2019-07-17 19:14
计算积分的方法:
  1. 方法1:
  2. function   y=fun(x)% %定义被积函数,生成一个m文件
  3. y=sin(x)+cos(x);   

  4. function  main()%生成另一个m文件
  5. T=quad(@fun,-1,1)           %对fun进行积分,下限为-1,上限为1

  6. T=
  7.     1.6829

  8. 方法2:
  9. fun=@(x)sin(x)+cos(x);   %定义被积函数
  10. T=quad(fun,-1,1)          %对fun进行积分,下限为-1,上限为1
复制代码
pxl0223
3楼-- · 2019-07-17 19:24
z00 发表于 2014-6-7 19:19
计算积分的方法:

谢谢你的回复~~~
pxl0223
4楼-- · 2019-07-17 19:34
913688247 发表于 2014-6-7 16:39
建议认真理清自己的思路   这样 看上去 你的代码太乱了。

是挺乱的,不怎么会编程,只是想计算出积分值

一周热门 更多>