Prev: Filtering
Next: cmyk tif bitmap imwrite fails
From: tlife on 29 Jun 2010 11:33
Dear all,

I recently started using Matlab and hope you can help me with the
following problem: I try to evaluate a double integral with a zeroth
order bessel function in the inner integral and a sum (n to infinity)
in the outer integral.

I tried nested trapz functions so far (I didn't succeed here, but I
guess it would be possible) and the dblquad function. My code seems to
work with the dblquad function, even though I'm not sure if Matlab
really interpretes the code the way I do (is there a way to set
breakpoints within the dblquad function?). When evaluating I get a lot
of warnings (Infinite or NaN function value detected... this seems to
be a problem of the bessel function ). So I hope you can give me some
advice here regarding the code, if its correct and if I somehow can
maybe improve it (computation time for a finer time discretization).

The formula loooks like this (simplified but I think you should get an
idea):

http://texify.com/?$out=\int_0^{tup} f(rd,tup)\cdot\left[\int_0^1
Bessel(rd,rdprime,td) \ drdprime\right]\cdot\left[1+\sum_{n=1}^{\infty}
f(n,tupprime)\right]dtupprime$

Thanks alot,
tlife


function doubletest
format long

maxtime=3000000;
n=10000;
for count=1:maxtime/n,
t=count*n;


a=7.9e-5;
b=3.62e-5;
c=6.38e-4;
d=0.115;
L=1200;
fh=194;
ft=0.05;
z=200;
r=0.;
total=200;
tup=(a*t)/(c*L^2);
rd=r/L;
hd=(tup/L)*sqrt(a/b);

result(count)= ((d*(L^2)*ft*1.3*10^-3)/(2*a*total))*dblquad(@integrnd,
0,1,0,tup,[],[],rd,tup,hd,fh,total,z,ft);

tplot(count)=t;

end

loglog((tplot),(result),'-');
grid off

function out = integrnd(rdprime,tupprime,rd,tup,hd,fh,total,z,ft)

summation=0;
i=1:100;
i=(exp(-((i.^2*(pi^2)*tup)/(hd^2)))).*cos(i.*pi*(fh/
total)).*cos(i.*pi*(z/total));

summation=sum(i);

extend=1+2*summation;

B=besselj(0,((rd.*rdprime)/(2.*tupprime)));

out = ((exp(-(rd^2/(4*tup))))/(2*tup)).*B.*exp(-((rdprime.^2)/
(4*tup))).*rdprime*extend;
 | 
Pages: 1
Prev: Filtering
Next: cmyk tif bitmap imwrite fails