Message Boards Message Boards

0
|
4785 Views
|
2 Replies
|
2 Total Likes
View groups...
Share
Share this post:

Convert the matlab bisection code into Wolfram Language?

Posted 6 years ago

(Convert matlab code into mathematica code)

(I'm fairly new to mathematica)

(bisect function matlab code:)

 function [root,fx,ea,iter]=bisect(func,xl,xu,es,maxit,varargin)

(% bisect: root location zeroes)

(% [root,fx,ea,iter]=bisect(func,xl,xu,es,maxit,p1,p2,...):)

(% uses bisection method to find the root of func)

(% input:)

(func = name of function)

(% xl, xu = lower and upper guesses)

(% es = desired relative error (default = 0.0001%))

(% maxit = maximum allowable iterations (default = 50))

(% p1,p2,... = additional parameters used by func)

(% output:)

(% root = real root)

(% fx = function value at root)

(% ea = approximate relative error (%))

(% iter = number of iterations)

 if nargin<3,error('at least 3 input arguments required'),end

 test = func(xl,varargin{:})*func(xu,varargin{:});

 if test>0,error('no sign change'),end

 if nargin<4|isempty(es), es=0.0001;end

 if nargin<5|isempty(maxit), maxit=50;end

 iter = 0; xr = xl; ea = 100;

while (1)

xrold = xr;

 xr = (xl + xu)/2;

 iter = iter + 1;

 if xr ~= 0,ea = abs((xr - xrold)/xr) * 100;end

 test = func(xl,varargin{:})*func(xr,varargin{:});

 if test < 0

 xu = xr;

elseif test > 0

 xl = xr;

 else

 ea = 0;

 end

 if ea <= es | iter >= maxit,break,end

 end

 root = xr; fx = func(xr, varargin{:});

(Problem code)

 clc

 clear all

 close all

 r1 = 0.5;

 r2 = 1;

 h = 1;

 pf= 200;

 pw = 1000;

 fh1 = @(h1)pf*pi*h/(3)*(r1^2 + r2^2 + r1*r2)-...

 pw *pi*(h - h1)/(3).*((r1 + (r2 - r1)/h*h1).^2 + r2^2 + (r1 + (r2 - 
 r1)/h*h1)*r2);

  h1 = (0:h/20:h);

  fhh =fh1(h1);

   plot (h1, fhh),grid

  [height f ea iter]= bisect(fh1, 0, 1);

(here's the picture of the exact problem and equations used to tackle the matlab problem.)

(I'm struggling with root function and bisection portion in mathematica)

https://imgur.com/a/LWaHEuM

(mathematica code i attempted )

func[root, fx, ea, iter] = Bisection[a0_, b0_, m_] :=
 Module[{}, a = N[a0]; b = N[b0]; c = (a + b)/2; k = 0; output = {{k, 
 a, c, b, f[c]}}; While[k < m, If[Sign[f[b]] == Sign[f[c]],
 b = c, a = c;]; c = a + b/2; i = i + 1;
 output = Append[output, {i, a, c, b, f[c]}];];
 Print[NumberForm[TableForm[output,
 TableHeadings \[RightArrow] {None, {"k", "ak", "ck", "bk",
 "f[ck]"}}], 16]];
 Print[" c = ", NumberForm[c, 16]];
 Print[" \[Laplacian]c = \[PlusMinus]", b \[Minus] a/2 ];
 Print["f[c] = ", NumberForm[f[c], 16] ];]

 r1 = 0.5;
 r2 = 1;
 h = 1;
 pf = 200;
 pw = 1000;
 fh1 = function[(h1)] pf*\[Pi]*h/(3)*(r1^2 + r2^2 + r1*r2) -
 pw *\[Pi]*(h - h1)/(3)*((r1 + (r2 - r1)/h*h1)^2 =
 r2^2 + (r1 + (r2 - r1)/h*h1)*r2)
 Plot[OutageProb, {h1, 0, h/20, h}, PlotRange -> {h1, fhh}, ,
 grid {h1, fhh}, Frame -> True]
 fhh = fh1 (h1)
 [height f ea iter] = bisection (fh1, 0, 1)
POSTED BY: Katya Claros
2 Replies

The main issue with func is a simple bug: a new variable i gets incremented each iteration instead of k. But starting with the line fh1 = ... the remaining code is gibberish in terms of Mathematica syntax, undefined functions, and the like. It really is not close enough to viability to offer much at this point-- it needs to be gone over with an instructor or TA.

POSTED BY: Daniel Lichtblau
Reply to this discussion
Community posts can be styled and formatted using the Markdown syntax.
Reply Preview
Attachments
Remove
or Discard

Group Abstract Group Abstract