Why do I get an error when I ran the integral?
Show older comments
Hi everyone, I got an error when I run following the code, I think because of the symbolic piecewise functions or the nested integrals, or maybe something else. Could anyone help here, please?
clear all; clear; clc;
hU = 100;
NU = 10;
x0 = 150;
a1 = 11.95;
b1 = 0.135;
drm = sqrt(hU^2 + x0^2);
DML = 0.57 * drm^1.6;
% AL = zeros(size(50:10:300));
j = 1;
Rs = 50:10:300;
for i = 50:10:300
xl = sqrt(i^2 - x0^2 + hU^2);
xh = sqrt(i^2 + x0^2 + hU^2);
syms x r
f(x) = piecewise(x >= hU & x <= xl, (2 * x) / i^2, x > xl & x <= xh, (2 * x) / (pi * i^2) * acos(((x^2) - hU^2 + x0^2 - i^2) / (2 * x0 * sqrt(x^2 - hU^2))));
f(r) = piecewise(r >= hU & r <= xl, (2 * r) / i^2, r > xl & r <= xh, (2 * r) / (pi * i^2) * acos(((r^2) - hU^2 + x0^2 - i^2) / (2 * x0 * sqrt(r^2 - hU^2))));
DLM(r) = r^1.6;
PL(r) = 1 / (1 + a1 * exp(-b1 * ((180/pi)*asin(hU / r) - a1)));
PN(r) = 1 - PL(r);
PL(x) = 1 / (1 + a1 * exp(-b1 * ((180/pi)*asin(hU / x) - a1)));
PN(x) = 1 - PL(x);
AL(j) = NU * vpa(int(PL(r) * f(r) * (int(PL(x) * f(x), r, xh) + int(PN(x) * f(x), DLM(r), xh))^((NU)-1), hU, DML), 6);
j = j + 1;
end
plot(Rs, AL, 'b--o')
grid on
xlabel('Avg. Cell radius in meters')
ylabel('Association Probability')
legend('LOS Association Prob.')
5 Comments
Walter Roberson
on 7 Feb 2024
Edited: Walter Roberson
on 7 Feb 2024
clear all; clear; clc;
hU = 100;
NU = 10;
x0 = 150;
a1 = 11.95;
b1 = 0.135;
drm = sqrt(hU^2 + x0^2);
DML = 0.57 * drm^1.6;
% AL = zeros(size(50:10:300));
j = 1;
%Rs = 50:10:300;
Rs = 50:10:100;
for i = Rs
xl = sqrt(i^2 - x0^2 + hU^2);
xh = sqrt(i^2 + x0^2 + hU^2);
syms x r
f(x) = piecewise(x >= hU & x <= xl, (2 * x) / i^2, x > xl & x <= xh, (2 * x) / (pi * i^2) * acos(((x^2) - hU^2 + x0^2 - i^2) / (2 * x0 * sqrt(x^2 - hU^2))));
f(r) = piecewise(r >= hU & r <= xl, (2 * r) / i^2, r > xl & r <= xh, (2 * r) / (pi * i^2) * acos(((r^2) - hU^2 + x0^2 - i^2) / (2 * x0 * sqrt(r^2 - hU^2))));
DLM(r) = r^1.6;
PL(r) = 1 / (1 + a1 * exp(-b1 * ((180/pi)*asin(hU / r) - a1)));
PN(r) = 1 - PL(r);
PL(x) = 1 / (1 + a1 * exp(-b1 * ((180/pi)*asin(hU / x) - a1)));
PN(x) = 1 - PL(x);
AL(j) = NU * vpa(int(PL(r) * f(r) * (int(PL(x) * f(x), r, xl, xh) + int(PN(x) * f(x), DLM(r), xh))^((NU)-1), hU, DML), 6);
AL(j)
j = j + 1;
end
dAL = double(AL);
plot(Rs, dAL, 'b--o')
grid on
xlabel('Avg. Cell radius in meters')
ylabel('Association Probability')
legend('LOS Association Prob.')
Pooja Kumari
on 7 Feb 2024
What is the error you are getting?
Hadeel Obaid
on 7 Feb 2024
Hadeel Obaid
on 7 Feb 2024
Walter Roberson
on 8 Feb 2024
Several of your int() do not indicate which variable to integrate over, and are guessing incorrectly.
Answers (1)
Walter Roberson
on 7 Feb 2024
f(x) = piecewise(x >= hU & x <= xl, (2 * x) / i^2, x > xl & x <= xh, (2 * x) / (pi * i^2) * acos(((x^2) - hU^2 + x0^2 - i^2) / (2 * x0 * sqrt(x^2 - hU^2))));
There is no clause for x > xh or x < hU
5 Comments
Hadeel Obaid
on 7 Feb 2024
Walter Roberson
on 7 Feb 2024
piecewise() returns NaN if none of the conditions are met.
You do not have conditions covering the cases x < hU or x > xh
Hadeel Obaid
on 8 Feb 2024
Wafaa
on 23 Mar 2024
Edited: Walter Roberson
on 23 Mar 2024
Warning: Unable to check whether the integrand
exists everywhere on the integration interval.
> In symengine
In sym/int (line 162)
In pp (line 31)
Warning: Unable to check whether the integrand
exists everywhere on the integration interval.
R(x,y)=piecewise(y <= x & x <= 1 & 0 <= y, y + 1, x < y & 0 <= x & y <= 1, x + 1);
%define xi
xi=zeros(1,m);
for i=1:m
xi(i)=(i-1)/(m-1);
end
%define psi(i)=L(R(x,y)) where y=xi
Lpsi=sym(zeros(1,m));
psi=sym(zeros(1,m));
for i=1:m
% psi(i)=ff(xi(i))-(R(x,xi(i))^3)*(x/2)+diff(R(x,xi(i))^3)*((x/3)-(x^2/2))+diff(R(x,xi(i))^3,2)*((x/4)-(2*x^2/4)+(x^3/3));%int(K(x,y)*(R(x,y))^3,y,0,xi(i))+int(K(x,y)*(R(x,y))^3,y,xi(i),1);
psi(i)=ff(xi(i))-int(K(x,y)*(R(x,y))^3,y,0,1);
Lpsi(i)=ff(x)-int(K(x,y)*subs(psi(i),x,y)^3,y,0,1);
% Lpsi(i)=psi(i)^3*(x/2)+diff(psi(i)^3)*((x/3)-(x^2/2))+diff(psi(i)^3,2)*((x/4)-(2*x^2/4)+(x^3/3));
end
please help me
your xl is mostly complex-valued, which is outside the valid piecewise() range, so you mostly get NaN
hU = 100;
NU = 10;
x0 = 150;
a1 = 11.95;
b1 = 0.135;
drm = sqrt(hU^2 + x0^2);
DML = 0.57 * drm^1.6;
% AL = zeros(size(50:10:300));
j = 1;
%Rs = 50:10:300;
Rs = 50:10:100;
for i = Rs
xl = sqrt(i^2 - x0^2 + hU^2)
xh = sqrt(i^2 + x0^2 + hU^2)
syms x r
f(x) = piecewise(x >= hU & x <= xl, (2 * x) / i^2, x > xl & x <= xh, (2 * x) / (pi * i^2) * acos(((x^2) - hU^2 + x0^2 - i^2) / (2 * x0 * sqrt(x^2 - hU^2))))
%f(r) = piecewise(r >= hU & r <= xl, (2 * r) / i^2, r > xl & r <= xh, (2 * r) / (pi * i^2) * acos(((r^2) - hU^2 + x0^2 - i^2) / (2 * x0 * sqrt(r^2 - hU^2))));
DLM(r) = r^1.6;
%PL(r) = 1 / (1 + a1 * exp(-b1 * ((180/pi)*asin(hU / r) - a1)));
%PN(r) = 1 - PL(r);
PL(x) = 1 / (1 + a1 * exp(-b1 * ((180/pi)*asin(hU / x) - a1)));
PN(x) = 1 - PL(x);
inner3 = PL(x) * f(x)
inner1 = int(inner3, x, xl, xh)
inner4 = PN(x) * f(x)
inner2 = int(inner4, x, DLM(r), xh)
outer = int(PL(r) * f(r) * (inner1 + inner2)^((NU)-1), r, hU, DML)
AL(j) = NU * vpa(outer, 6);
AL(j)
j = j + 1;
end
dAL = double(AL);
plot(Rs, dAL, 'b--o')
grid on
xlabel('Avg. Cell radius in meters')
ylabel('Association Probability')
legend('LOS Association Prob.')
Categories
Find more on Creating and Concatenating Matrices in Help Center and File Exchange
on 7 Feb 2024
on 23 Mar 2024
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!
