Plot integral of a function

Question

Elinor Ginzburg 2020-12-27

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/703407-plot-integral-of-a-function

评论： Star Strider 2020-12-28

Hello,

I'm tryign to plot the integral of the function rho in the following code, but I'm having trouble with plotting it. Perhaps anyone can help me with the code for plotting the integral of this rho function in the different areas? This is what I've tried:

q = 1.5*10^-16;
a = 1;
Na = 10^16;
xn = 3;
xp = -2;
rho_righ = @(x) a.*x;
rho_left = @(x) -q*Na;
nop = @(x) 0;
figure(1);
fplot(rho_righ, [0,xn]);
hold on;
fplot(rho_left, [xp,0]);
hold on;
fplot(nop, [-5,xp]);
hold on;
fplot(nop, [xn,5]);
grid on;
legend('{\rho}=qax','{\rho}=-qNa');
hold off;
upper_limit = linspace(-xp, 0);
upper_limit2 = linspace(0, xn);
xval = arrayfun(@(uplim) integral(rho, 0, uplim, 'ArrayValued',true), upper_limit);
xval2 = arrayfun(@(uplim) integral(rho, 0, uplim, 'ArrayValued',true), upper_limit2);
figure(2)
plot(upper_limit, xval);
hold on;
plot(upper_limit2, xval2);
grid on;

Thank you very much!

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Star Strider 2020-12-27

1
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/703407-plot-integral-of-a-function#answer_585332

在 MATLAB Online 中打开

The problem appears to be that ‘rho’ itself does not exist.

I have no idea what you want ‘rho’ to be, however defining it as:

rho = @(x) rho_left(x) .* (x <= 0) + rho_righ(x) .* (x > 0);

will likely be an appropriate place to begin, and that produces figure(2) (and myriad Warning messages).. Whatever you intend it to be, it needs to be defined.

4 个评论
显示 2更早的评论隐藏 2更早的评论

Elinor Ginzburg 2020-12-28

在 MATLAB Online 中打开

Thank you very much for your answer, I totally missed that. I managed to plot the first integral of rho, and now I'm also trying to plot the second integral of this, but for some reason it starts to run and doesn't finish even after few minutes. perhaps you have any idea what could I do wrong? here's my new code:

q = 3;
a = 1/2;
Na = 1;
xn = 2;
xp = -1;
rho_righ = @(x) q*a.*x;
rho_left = @(x) -q*Na;
rho = @(x) rho_left(x) .* (x <= 0 && x > xp) + rho_righ(x) .* (x > 0 && x <= xn);
figure(1);
fplot(rho);
grid on;
xlabel('x [{cm^{-3}}]');
ylabel('{\rho} [{cb/cm^{-3}}]');
title('Charge Density Distribution');
upper_limit = linspace(-10, 10);
EF = @(uplim) integral(rho, 0, uplim, 'ArrayValued',true) - 3;
Electrical_field = arrayfun(EF, upper_limit);
figure(2)
plot(upper_limit, Electrical_field);
grid on;
xlabel('x [{cm^{-3}}]');
ylabel('E [V/cm]');
title('Electrical Field');
V = @(uplim) integral(EF, 0, uplim, 'ArrayValued',true);
voltage = arrayfun(V, upper_limit);
figure(3)
plot(upper_limit, voltage);
grid on;
xlabel('x [{cm^{-3}}]');
ylabel('V [V]');
title('Voltage');

Thank you very much for your help!

Star Strider 2020-12-28

编辑：Star Strider 2020-12-28

在 MATLAB Online 中打开

My pleasure!

First, arrayfun has its uses, however it’s simply a for loop in disguise, so it’s more efficient to just use a for loop. Also, the ‘&&’ is a ‘short-circuit’ operator for strictly logical values. For vectors, use a single ‘&’ operator.

Second, I’m not certain what you’re doing in the second integration. It looks as though you’re doing a double integration over the same limits, where ‘V’ is integrating ‘EF’ is integrating ‘rho’. There might be ways to make that more tractable if I understand what you’re doing.

Meanwhile, while I wait for you to clarify what the the ‘voltage’ integration is (possibly to make it faster and more efficient), your (slightly-revised) code is:

q = 3;
a = 1/2;
Na = 1;
xn = 2;
xp = -1;
rho_righ = @(x) q*a.*x;
rho_left = @(x) -q*Na;
rho = @(x) rho_left(x) .* ((x <= 0) & (x > xp)) + rho_righ(x) .* ((x > 0) & (x <= xn));
figure(1);
fplot(rho);
grid on;
xlabel('x [{cm^{-3}}]');
ylabel('{\rho} [{cb/cm^{-3}}]');
title('Charge Density Distribution');
upper_limit = linspace(-10, 10);
EF = @(uplim) integral(rho, 0, uplim, 'ArrayValued',true) - 3;
for k = 1:numel(upper_limit)
    Electrical_field(k) = EF(upper_limit(k));
end
figure(2)
plot(upper_limit, Electrical_field);
grid on;
xlabel('x [{cm^{-3}}]');
ylabel('E [V/cm]');
title('Electrical Field');
V = @(uplim) integral(EF, 0, uplim, 'ArrayValued',true);
for k = 1:numel(upper_limit)
    voltage(k) = V(upper_limit(k));
end
figure(3)
plot(upper_limit, voltage);
grid on;
xlabel('x [{cm^{-3}}]');
ylabel('V [V]');
title('Voltage');

EDIT — (28 Dec 2020 at 13:28)

It only ended up taking 1647.612439 seconds (27.46 minutes) for all 100 loop iterations. In the end, it produced:

To save you the trouble of integrating it, the ‘voltage’ vector is:

voltage = [1.49999729339940 1.49999696906513 1.50015522691315 1.50000236066725 1.49999973680205 1.50000048145849 1.49999931842135 1.49999972895611 1.49999906076886 1.49999918636686 1.49999895100276 1.49999918002933 1.49999971120408 1.49999952794135 1.49999962815786 1.50000030206563 1.49999982037269 1.50000168699648 1.50000016177641 1.49999957361846 1.49999976226557 1.49999946542346 1.50000107234303 1.49999866049175 1.50000007607924 1.49999936777078 1.49999950297102 1.50000078110427 1.49999943944442 1.50000006989847 1.50000061270948 1.49999949173887 1.49999992644981 1.49999988813672 1.49999975368376 1.49999941182790 1.49999970577751 1.50000030183291 1.49999923594608 1.49999958590485 1.49999927394898 1.49999985898082 1.50000059556414 1.49999973974386 1.50000074795001 1.48760330578512 1.37128864401592 1.13253749617386 0.771349862258953 0.287725742271197 -0.302772650492271 -0.902134290564043 -1.48294494789750 -2.03283730066711 -2.53944402704733 -2.99039780521262 -3.37333131333744 -3.67587722959624 -3.88566823216349 -3.99033699921365 -4.00000064537387 -4.00000095524024 -3.99999924718646 -4.00001441974673 -3.99999995888804 -4.00000125703370 -4.00000118209333 -4.00000010402610 -3.99999897834425 -4.00000301314983 -3.99999985682291 -3.99999997325554 -4.00000049834742 -3.99999985595181 -4.00000027641440 -4.00000226576712 -4.00000150165796 -4.00000056720757 -3.99999992990422 -4.00000016104336 -3.99999930707085 -3.99999881231945 -4.00000009175432 -3.99999780633590 -3.99999761823918 -4.00000016397481 -3.99999902244021 -3.99999924718711 -3.99999907700015 -4.00000099591624 -3.99999769551857 -4.00000009842234 -3.99999631419898 -3.99999899595210 -3.99999286143871 -3.99999971820702 -3.99999724890175 -3.99999932916977 -3.99999697121052 -3.99999797250258];

I was able to recover it from my Workspace window and copy it to here.

Elinor Ginzburg 2020-12-28

Thank you so much!

yeah that was exactly what I ment, and the results match my calculations.

Thank you very much for your help and advices!

Star Strider 2020-12-28

As always, my pleasure!

请先登录，再进行评论。

Plot integral of a function

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

采纳的回答

4 个评论
显示 2更早的评论隐藏 2更早的评论

更多回答（0 个）

另请参阅

类别

标签

产品

版本

Community Treasure Hunt

Plot integral of a function

0 个评论 显示 -2更早的评论隐藏 -2更早的评论

采纳的回答

4 个评论 显示 2更早的评论隐藏 2更早的评论

更多回答（0 个）

另请参阅

类别

标签

产品

版本

Community Treasure Hunt

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

4 个评论
显示 2更早的评论隐藏 2更早的评论