How can I use the Bisection method to converge on an arbitrary value?

7 ビュー (過去 30 日間)
Logan Parrish
Logan Parrish 2022 年 1 月 28 日
コメント済み: Logan Parrish 2022 年 1 月 29 日
I am working on coding a function to solve for the water level inside a water tank. The volume of the tank is described by:
V = pi * x^2 *((3*R-x)/3)
I have R = 4 m and I need to find the height (x) where V = 100 m^3. I chose to use the Bisection Method to find this height and use a while loop that should end when the height that gives out Volume equal to 100 m^3. Though I am having a hard time getting the loop to run more than just once and I don't know what could be going wrong. I'm new to Matlab so I might be missing something very obvious but any help would be appreciated.
Here is the code that I have managed to come up with so far:
% Define function
R = 4;
v = @(x) pi.*x.^2.*((3.*R-x)./3);
% Define variables
xL = 0;
xH = 2.*R;
vL = v(xL);
vH = v(xH);
% Define xr
xr = .5.*(xL+xH); % Bisection method
vR = v(xr);
% Define the loop's target value
tar = 100;
while v(xr) == tar
xr = .5.*(xL+xH);
% Update xH and xL
if vH>vR
xH = xr;
else
xL = xr;
end
end
root = xr;
Once again thanks for any help it is very appreciated

採用された回答

Voss
Voss 2022 年 1 月 28 日
In fact your loop runs 0 times (because v(xr) is not equal to 100).
% Define function
R = 4;
v = @(x) pi.*x.^2.*((3.*R-x)./3);
% Define variables
xL = 0;
xH = 2.*R;
vL = v(xL);
vH = v(xH);
% Define xr
xr = .5.*(xL+xH); % Bisection method
vR = v(xr);
% Define the loop's target value
tar = 100;
disp(v(xr));
134.0413
while v(xr) == tar
disp('looping');
xr = .5.*(xL+xH);
% Update xH and xL
if vH>vR
xH = xr;
else
xL = xr;
end
end
root = xr;
It should loop until v(xr) == tar, so: while v(xr) ~= tar, and in fact it probably should not be with exact (in)equality but using some small error, e.g., while abs(v(xr) - tar) > max_error. And you need to update your volumes as you go and use tar to tell which of xH or xL need to be updated each time.
% Define function
R = 4;
v = @(x) pi.*x.^2.*((3.*R-x)./3);
% Define variables
xL = 0;
xH = 2.*R;
vL = v(xL);
vH = v(xH);
% Define xr
xr = .5.*(xL+xH); % Bisection method
vR = v(xr);
% Define the loop's target value
tar = 100;
max_error = 1e-6;
while abs(v(xr) - tar) > max_error
xr = .5.*(xL+xH);
vR = v(xr); % update vR too
fprintf(1,'looping: xr = %12.8f, vR = %12.8f\n',xr,vR);
% Update xH and xL
if tar < vR % use tar to determine which side of the target you're on
xH = xr;
vH = v(xH); % update vH too
else
xL = xr;
vL = v(xL); % update vL too
end
end
looping: xr = 4.00000000, vR = 134.04128655 looping: xr = 2.00000000, vR = 41.88790205 looping: xr = 3.00000000, vR = 84.82300165 looping: xr = 3.50000000, vR = 109.03944502 looping: xr = 3.25000000, vR = 96.78396118 looping: xr = 3.37500000, vR = 102.88102348 looping: xr = 3.31250000, vR = 99.82405544 looping: xr = 3.34375000, vR = 101.35052611 looping: xr = 3.32812500, vR = 100.58677545 looping: xr = 3.32031250, vR = 100.20528511 looping: xr = 3.31640625, vR = 100.01463751 looping: xr = 3.31445312, vR = 99.91933825 looping: xr = 3.31542969, vR = 99.96698583 looping: xr = 3.31591797, vR = 99.99081115 looping: xr = 3.31616211, vR = 100.00272420 looping: xr = 3.31604004, vR = 99.99676765 looping: xr = 3.31610107, vR = 99.99974592 looping: xr = 3.31613159, vR = 100.00123506 looping: xr = 3.31611633, vR = 100.00049049 looping: xr = 3.31610870, vR = 100.00011820 looping: xr = 3.31610489, vR = 99.99993206 looping: xr = 3.31610680, vR = 100.00002513 looping: xr = 3.31610584, vR = 99.99997859 looping: xr = 3.31610632, vR = 100.00000186 looping: xr = 3.31610608, vR = 99.99999023 looping: xr = 3.31610620, vR = 99.99999604 looping: xr = 3.31610626, vR = 99.99999895 looping: xr = 3.31610629, vR = 100.00000041
format long;
root = xr
root =
3.316106289625168
v(root)
ans =
1.000000004077111e+02
  1 件のコメント
Logan Parrish
Logan Parrish 2022 年 1 月 29 日
Thanks a ton! I was really struggling trying to wrap my head around how to get my function to approach the target value and your explanation really understand where I went wrong.

サインインしてコメントする。

その他の回答 (0 件)

カテゴリ

Help Center および File ExchangeLoops and Conditional Statements についてさらに検索

製品


リリース

R2021b

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by