数值积分的问题

chenyf03

a=0:0.1:5;
num=exp(0.5);
mul=1./(2.*pi);
F=1:1:length(a);
R=1:1:length(a);
result=1:1:length(a);
for i=1:1:length(a)
    R(i)=a(i)*exp((-a(i)^2)./2);
    F(i)=quad(@(y)y.*exp(-y.^2./2).*(log(1+besl(R(i).*y)./num)+(besl(R(i).*y)./num).*log(1+num./besl(R(i).*y))),0,100);
    result(i)=log(2)-0.5.*F(i);
end
figure;
plot(a,result,'k-');

function r=besl(x)
  syms x;
  r=quad(@(z)(exp(x.*cos(z))),0,2.*pi);

老是出错，请问各位大虾，应该怎样改正，或是怎样调用下面的函数？
小弟感激不尽啊！！！

[ 本帖最后由 eight 于 2007-5-18 18:46 编辑 ]

sffei

quad 是数值积分，不能用于符号计算

chenyf03

那怎么用,符号积分?谢谢

sffei

你的F函数太复杂，没时间调试。改了一个简单的，能运行。你自己改吧
a=0:0.1:5;
num=exp(0.5);
mul=1./(2.*pi);
F=1:1:length(a);
R=1:1:length(a);
besel=@(x)quad(@(z)(exp(x.*cos(z))),0,2.*pi);
result=1:1:length(a);
for i=1:1:length(a)
    R(i)=a(i)*exp((-a(i)^2)./2);
    F(i)=quad(@(y)besel(R(i))*y,0,100);
    result(i)=log(2)-0.5.*F(i);
end
figure;
plot(a,result,'k-');

eight

原帖由 chenyf03 于 2007-5-18 15:33 发表

                               
登录/注册后可看大图

那怎么用,符号积分?谢谢

本版这类问题讨论得很多，如果仍有问题，请搜索一下，并多动手试试

chenyf03

F(i)=quad(@(y)besel(R(i))*y,0,100);
应该是将y放在besel()里,如:besel(R(i)*y)
在里边就出不来了,请sffei和eight帮忙解答一下,小弟先谢过了!
:@)

sffei

这个函数exp(x.*cos(z))对z进行积分根本积不出来，只能用数值积分，而数值积分想求这个besel(y)的表达式是做不到的

[ 本帖最后由 eight 于 2007-5-18 23:24 编辑 ]