Problem on integral and ploting...

1 次查看(过去 30 天)
ZhiH
ZhiH 2012-12-17
Hello there!
I do have a problem with symbolic integration. The original code looks like this:
function f = fmbloss(phi)
n1 = 1.491;
n2 = 1.405;
ns1 = n1*n1;
ns2 = n2*n2;
lambda = 0.58e-6;
ra1 = 0.5e-3;
kw = 2*pi/lambda;
bend_r = 10e-3;
thetac = acos(n2/n1);
rad = @(theta) (bend_r+ra1)*n1./(n2*cos(theta));
n1_phi = @(theta) (sqrt(ns1-ns1*cos(theta).*cos(theta)));
n2_phi = @(theta) (sqrt(ns1*cos(theta).*cos(theta)-ns2));
p = @(theta) (rad(theta)/(bend_r+ra1));
ps = @(theta) (p(theta).*p(theta));
ps1 = @(theta) (sqrt(ps(theta)-1));
log_p = @(theta) (log(p(theta)+ps1(theta)));
abs_t = @(theta) (4*n1_phi(theta).*n2_phi(theta)/(ns1-ns2));
expf = @(theta) (-2*kw*ra1*n1*cos(theta).*(log_p(theta)-ps1(theta)./p(theta)));
bt = @(theta) (abs_t(theta).*expf(theta));
sub_sin = @(theta) (sqrt(sin(theta).*sin(theta)-sin(thetac)*sin(thetac)));
sub_sins = @(theta) (sin(theta)+sub_sin(theta).*sub_sin(theta));
bt1 = @(theta) (4*sin(theta).*sub_sin(theta)./(sub_sins(theta)));
exp_f1 = @(theta) exp(-bt(theta)*phi./(2*theta));
exp_f2 = @(theta) exp(-bt(theta)*bend_r*phi.*theta/(4*ra1));
exp_f3 = @(theta) exp(-bt1(theta)*bend_r*phi/(4*ra1));
fun1 = @(theta) (exp_f1(theta));
fun2 = @(theta) (exp_f2(theta));
fun3 = @(theta) (exp_f3(theta));
theta_max = sqrt(2*ra1/bend_r);
int1 = quadl(fun1,0,theta_max);
int2 = quadl(fun2,theta_max,thetac);
int3 = quadl(fun3,thetac,pi/2);
f = pi*(int1+int2+int3)/2;
and input fmbloss(pi/3) in the command window firstly,
and error message:
>> fmbloss(pi/3)
Warning: Infinite or Not-a-Number function value encountered.
> In quadl at 105
In fmbloss at 40
Warning: Infinite or Not-a-Number function value encountered.
> In quadl at 105
In fmbloss at 41
ans =
Inf
input fmbloss(0) in the command window secondly,
>> fmbloss(0) ans = 2.467401100272340
It is not working at all.I hope someone can give me a good advise, how to deal with that. Thank you!
  1 个评论
Walter Roberson
Walter Roberson 2012-12-17
Your routine contains no symbolic integration ?
You are hitting infinity in your integration. You should learn how to use the debugger to figure out why. In particular,
dbstop if naninf

请先登录,再进行评论。

回答(0 个)

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by