Not getting numerical answer of integration

height = 2.5; % in units of meter.
a=5; % considering square room of area 20*20 meter^2.
area=4*a^2;
psi_FOV = 30*pi/180;
psi_half = 60*pi/180; 
m = -log10(2)/log10(cos(psi_half));
alpha = (m+1)*height^(m+1)/(2*pi);
K = alpha^2;
N = (1e-20)/K;
r_fov_anlytical = height*tan(psi_FOV);
d_fov_anlytical = height/cos(psi_FOV);
beta = (m+3);
height = 2.5;
% lambdaaxis=  .05:.01:1;
% avgSINR=zeros(1,length(lambdaaxis));
% for lambdaIdx=1:length(lambdaaxis)
    syms x y z;
   
     ons= (vpaintegral(vpaintegral(exp(-2.*pi.*.05.*vpaintegral((z - z.*exp ((-y.*N.*(sin( psi_FOV )).^4)./((x.^2+height^2 ).^(-beta))) ...
         .*exp ((-y.*( z.^2+ height.^2) .^( - beta))./(( x.^2+height.^2).^(-beta)))),z,x,r_fov_anlytical)) ...
         .*exp(- .05.*pi.*(x).^2).*(2.*pi).*.05.*x,y,0,inf),x,0,r_fov_anlytical))   
    
%      lambdaIdx

Every time I run this code, i just see equation getting printed in the command window. It should rather give numerical result. When I use double, that too didn't worked. matlabFunction and integral2 combination also didn't worked. I have to calculate this equation for many values , Is there any way by which I can compute it in less time.

[madhan ravi]

Hhh... stop asking the same question billion times.

[HIMANSHU GAUTAM]

I didn't got answer even one time. I thought framing it in good way might help.

[Walter Roberson]

Is this different than https://www.mathworks.com/matlabcentral/answers/543824-how-to-solve-this-integral-equation-with-integral2-or-monte-carlo-methods ?

[HIMANSHU GAUTAM]

The equation mentioned there is slightly different than equation mentioned here. I implemented that equation (previous question) , here is the code for that:

height = 2.5; % in units of meter.
a=5; % considering square room of area 20*20 meter^2.
area=4*a^2;
psi_FOV = 30*pi/180;
psi_half = 60*pi/180; 
m = -log10(2)/log10(cos(psi_half));
alpha = (m+1)*height^(m+1)/(2*pi);
K = alpha^2;
N = (1e-20)/K;
r_fov_anlytical = height*tan(psi_FOV);
d_fov_anlytical = height/cos(psi_FOV);
beta = (m+3);
height = 2.5;
lambdaaxis=  .05:.01:1;
avgSINR=zeros(1,length(lambdaaxis));
for lambdaIdx=1:length(lambdaaxis)
    syms x y z;
   
     avgSINR(lambdaIdx)= vpaintegral(vpaintegral(exp(-2.*pi.*lambdaaxis(lambdaIdx) ...
         .*vpaintegral((z - z.*exp ((-y.*( z.^2+ height.^2) .^( - beta) )./(( x.^2+height.^2).^(-beta)))),z,x,r_fov_anlytical)) ...
         .* exp ((-y.*N.*(sin( psi_FOV )).^4)./((x.^2+height^2 ).^(-beta))).* exp(- lambdaaxis(lambdaIdx).*pi.*(x).^2).*(2.*pi).*lambdaaxis(lambdaIdx).*x,y,0,inf),x,0,r_fov_anlytical);   
    
     lambdaIdx
    end

I got the result for this. Since those results were not fully correct, so i tweaked the equation. But this time I didn't got numerical answer. I getting following in my command window.

     
ons =
 
vpaintegral(vpaintegral((x*pi*exp(-(pi*vpaintegral(z - z*exp(-(y*(x^2 + 25/4)^4)/(z^2 + 25/4)^4)*exp(-(3358452346175993*y*(x^2 + 25/4)^4)/21267647932558653966460912964485513216), z, x, (5*3^(1/2))/6))/10)*exp(-(pi*x^2)/20))/10, y, 0, Inf), x, 0, (5*3^(1/2))/6)
 
>>