Integration HELP Mie Scattering

Dear friends, i am really in trable, i am trying to find the specific attenuation due to rain but i have problem in my code and i am asking you your help. My problem is in the last command which is an integration operation, and i am sure you will do your best for helping me. The link below is for the equations and it's followed by the code.

http://www.mediafire.com/view/u7ktl5tj5jy700u/nnnn.png

radius = 1;

nMax = 40; % maximum mode number

No=8*10^3; %m^ -4

R1=140; %rainfall rate=1mm/hr

AS1=8.2*R1^-0.21;

N1D=No*exp(-AS1*radius*1e-3);

% mode numbers

mode = 1:nMax;

frequency = 6e9;

% speed of light

c = 299792458.0;

lambda = c / ( frequency ) ;

for n=1:10;

n2 = (2*n+1);

% radian frequency

w = 2.0*pi*frequency;

% wavenumber

k = w/c;

% conversion factor between cartesian and spherical Bessel/Hankel function

s = sqrt(0.5*pi/(k*radius));

% compute spherical bessel, hankel functions

[J(mode)] = besselj(mode + 1/2, k*radius); J = J*s;

[H(mode)] = besselh(mode + 1/2, 2, k*radius); H = H*s;

[J2(mode)] = besselj(mode + 1/2 - 1, k*radius); J2 = J2*s;

[H2(mode)] = besselh(mode + 1/2 - 1, 2, k*radius); H2 = H2*s;

% derivatives of spherical bessel and hankel functions % Recurrence relationship, Abramowitz and Stegun Page 361

kaJ1P(mode) = (k*radius*J2 - mode .* J );

kaH1P(mode) = (k*radius*H2 - mode .* H );

% Ruck, et. al. (3.2-1)

An = -((1i).^mode) .* ( J ./ H ) .* (2*mode + 1) ./ (mode.*(mode + 1));

% Ruck, et. al. (3.2-2), using derivatives of bessel functions

Bn = ((1i).^(mode+1)) .* (kaJ1P ./ kaH1P) .* (2*mode + 1) ./ (mode.*(mode + 1));

Qt = (lambda^2/2*pi)*sum(n2).*sum(real(An + Bn));

end

Co1=4.343*No;

==============================================================

NOW SHOULD BE THE INTEGRATION OPERATION

Answer 1

[Mike Hosea]

Since I don't know much about this application area, I'm not in position to work on the code here. I would approach a problem like this in stages. The first stage is to write a function that you're integrating, with all parameters being represented as inputs to that

Step 1: Write a function that takes all your parameters as inputs and produces the value of the integrand for those inputs. That is to say, if I wanted to integrate f(x,a,b)*g(x,c) over x, I would write a function with the signature fun(x,a,b,c). In this case it will be fun(r,m,lambda).

Step 2: If possible make this function "vectorized" in the integration variable. Usually this means using .* instead of *, ./ instead of /, and .^ instead of ^. You will not need to vectorize it in all the variables, as these other variables will be fixed scalars during the integration. Sometimes vectorization is not practical, in which case you can skip this step, but you will need to add "'ArrayValued',true" to the argument list of INTEGRAL in step 3.

Step 3: Define the integral as an anonymous function. It looks like the integral you want is an antiderivative, meaning that you probably want to fix the left-hand limit of integration and to consider the right-hand limit to be r, so I'll write it this way. Just realize that you could easily make the limits variable. If the integrand is fun(r,m,lambda), then you would write

kscalar = @(r,m,lambda)integral(@(r1)fun(r1,m,lambda),a,r);

If you skipped step 2, then

kscalar = @(r,m,lambda)integral(@(r1)fun(r1,m,lambda),a,r,'ArrayValued',true);

Step 4: Vectorize the integral as needed for plotting, tabulating, or whatever.

k = @(r,m,lambda)arrayfun(kscalar,r,m,lambda);

With this definition, all inputs need to be the same size, so if you wanted scalar expansion on one input, e.g. r, you'd need to write something like k(r*ones(size(m)),m,lambda) where you use k.

Not sure if that will help or not, but the bottom line is at first ignore the integral. First focus on writing a MATLAB function that accepts all the variables and parameters (and at first, only scalar values of these), and then once you have that working, move towards integrating. There's no sense in trying to integrate anything until you have coded up a function that can correctly evaluate the integrand. -- Mike

[Mike Hosea]

Some additional explanation on

kscalar = @(r,m,lambda)integral(@(r1)fun(r1,m,lambda),a,r);

Here we are definining a function named kscalar. The @(r,m,lambda) defines the argument list for this function. The body of the function is to evaluate INTEGRAL. The first argument to INTEGRAL is an "anonymous function". This is just a function that has been constructed "on the fly". Since INTEGRAL is a 1-D integration function, it only accepts functions of 1 variable. The expression @(r1)fun(r1,m,lambda) defines a function of one variable, r1. The value of the function is obtained by evaluating the function "fun" at the arguments r1, m, and lambda. Since m and lambda do not appear in the argument list, their values are copied from the current values of these variables in the workspace at the moment the function is defined, which will be the instant INTEGRAL is called. Since r1 appears in the argument list, however, it is just a "dummy variable". That is to say, it doesn't matter what the variable is named. You can call it "r" instead of "r1". MATLAB will not be confused between this "r" and any other "r" variables floating around. The function "fun" is the one you were working on in the early stages. In the present context it is a function of 3 variables r, m, and lambda defined mathematically (per your png file) as fun(r,m,lambda) = Q_t(r,lambda,m)*n(r). Sorry for swapping the order of m and lambda, but I started out doing that above, and I don't want to add confusion by changing it at this point.

Once you have defined kscalar, you can evaluate it with scalar inputs. For example if r = 1, m = 2, and lambda = 3, then you can evaluate kscalar(1,2,3) to get the value of k corresponding to those values. Vectorizing this in the next step will help you if you want to plot it, i.e. so you can write k(1:10,2*ones(1,10),3*ones(1,10)) to evaluate k at r = 1:10 with m = 2 and lambda = 3.

[sese]

Thank you so much Mike for the detailed explanation BUT what is a and r??

[sese]

Sorry Mike for so many question but i have no choice because i have to present my progress report of this semester next week.

I have wrote in my m-file the following command

kscalar = @(radius,m,lambda)integral(@(radius1)(Qt.*N1D),a,r)

and this was the MATLAB answer

kscalar =

@(r,m,lambda)integral(@(r1)(Qt.*N1D),a,r)

[Mike Hosea]

Your reference defines k as an integral, but not a definite integral. This is the so-called Fundamental Theorem of Calculus. Here k is the anti-derivative of a function of r with m and lambda as parameters. So k is a function of r with m and lambda as parameters. I am doing numerical integration, not symbolic, so the only way I can represent an anti-derivative is as a definite integral with a variable upper bound, e.g. ln(t) = integral(@(x)1./x,1,t). For example, in MATLAB

>> lns = @(t)integral(@(x)1./x,1,t);
>> lns(3)
ans =
       1.0986
>> log(3)
ans =
       1.0986

But this function only works for scalar inputs

>> lns(1:10)
Error using integral (line 85)
A and B must be floating point scalars.
Error in @(t)integral(@(x)1./x,1,t)

To make my function more useful in MATLAB, I need to vectorize it:

>> ln = @(t)arrayfun(lns,t);
>> ln(1:4)
ans =
            0      0.69315       1.0986       1.3863

To use the kscalar or k functions, you provide the arguments to evaluate them, e.g., k(1,2,3) should return a number corresponding to the value of r=1, m=2, and lambda=3. If r should have been integrated out, then the antiderivative in your reference would need to be a definite integral with fixed limits.

[sese]

I tried to apply what you explained but Unfortunately couldn't get the answer................ Could you apply my function to pleaseeeeeeeeeeeeee