学习笔记：郭彥甫 (Yan-Fu Kuo), 台大生機系 ，MATLAB教學 - 11方程式求根

Problem Statement

• Suppose you have a mathematical function $f(x)$ and you want to find $x_0$such that $f(x_o)=0$, e.g.

  $f(x) = x^2 - 2x - 8 = 0$

• Solution1:Analytical Solution

• Solution2:Graphical Illustration

• Solution3:Numerical Solution

Symbolic Root Finding Approach

• Performing mathematics on symbols, NOT numbers

• The symbols math are performed using “symbolic variables”

• Use sym or symsto create symbolic variables(符号变量).

  syms x
  % or x = sym('x');
  x + x + x
  (x + x + x)/4

• Define : $y = x^2-2x-8$ （y is a symbolic variable）

Symbolic Root Finding:solve()

• Function solvefinds roots for equations(方程)

  $$y = x\cdot{sin(x)} - x = 0$$

  syms x
  y = x*sin(x)-x;
  solve(y,x)
  %solve('x*sin(x)-x',x)

• Find the roots for:

  $cos(x)^2-sin(x)^2 = 0$ and $cos(x)^2 + sin(x)^2 = 0$

Solving Multiple Equations

• Solve this equation using symbolic approach:

  $$\begin{cases} x-2y= 5\\ x+y= 6 \end{cases}$$

syms x y
eq1 = x-2*y-5;
eq2 = x+y-6;
A = solve(eq1,eq2,x,y)

Solving Equations Expressed in Symbols

• What if we are given a function expressed in symbols?

  $ax^2-b=0$

  syms x a b
  solve('a*x^2-b')

• $x$ is always the first choice to be solved

• What if one wants to express $b$ in tems of $a$ and $x$ ?

  syms x a b
  solve('a*x^2-b','b')

Exercise

• Solve this equation for $x$ using symbolic approach

  $(x-1)^2 + (y-b)^2 = r^2$

  syms x y a b r
  solve((x-a)^2+(y-b)^2 - r^2)

• Find the matrix inverse using symbolic approach

syms a b c d
A = [a b; c d];
B = inv(A)

Symbolic Differentiation : diff()

• Caculate the derivative of a symbolic function:

  $y = 4x^5$

  syms x
  y = 4*x^5;
  yprime = diff(y)

• Exercise:

$f(x) = \frac{e^{x^2}}{x^3-x+3} , \frac{df}{dx}=?$

syms x
f = exp(x^2)/(x^3-x+3);
fprime = diff(f)

$f(x) = \frac{x^2+xy-1}{y^3+x+3} , \frac{\partial f}{\partial x}=?$

syms x y
f = (x^2+x*y-1)/(y^3+x+3);
fprime = diff(f)

Symbolic Integration : int()

• Calculate the integral of a symbolic function:

  $z = \int ydx = \int x^2e^xdx , z(0)=0$

  syms x;
  y = x^2*exp(x);
  z = int(y);
  z = z-subs(z,x,0)

Exercise:

syms x;
y = (x^2-x+1)/(x+3);
z = int(y,0,10)

Symbolic vs Numeric

Review of Function Handles(@)

• A handle is a pointer to a function

• Can be used to pass functions to other functions

• For example, the input of the following function is another function:

  function [y] = xy_plot(input,x)  % input is a function variable
  % xy_plot receives the handle of a function and plots that
  % function of x
  y = input(x);
  plot(x,y,'r--');
  xlabel('x');
  ylabel('function(x)');
  end

  • When we need to use the function xy_plot ,use xy_plot(@sin,0:0.01:2*pi);

fsolve()

• A numeric root solver

• For example, solve this equation:

  $f(x) = 1.2x+0.3+x\cdot sin(x)$

  f2 = @(x) (1.2*x+0.3+x*sin(x));  % inline function
  fsolve(f2,0)  % f2 is a function handle ; 0 is a initial guess

Exercise

• Find the root for this equation:

$f(x,y)=\begin{cases}2x-y-e^{-x} \\ -x+2y-e^{-y}\end{cases}$

• using initial value(x,y) = (-5,5)

fzero()

• Another numeric root solver

• Find the zero if and only if the function crosses the x-axis

  f = @(x) x.^2
  fzero(f,0.1)  % 0.1 is a initial guess
  fsolve(f,0)

• Options:

  f = @(x) x.^2
  options = optimset('MaxIter',le3,'TolFun',1e-10);  % 1e3 is 'Number of iterations'; 1e-10 is Tolerance(误差)
  fsolve(f,0.1,options)
  fzeros(f,0.1,options)

Finding Roots of Polynomials : roots()

• Find the roots of this polynomial:

  $f(x) =x^5-3.5x^4+2.75x^3+2.125x^2-3.875x+1.25$

  roots([1 -3.5 2.75 2.125 -3.875 1.25])

• roots only works for polynomials

How Do These Solvers Find the Roots?

• Now we are going to introduce more details of some numeric methods

Numeric Root Finding Methods

• Two major types:
  
  • Bracketing methods(3.g.,bisection method)
  • Start with an interval that contains the root
  • Open methods(3.g.,Newton-Raphson method)
  • Start with one or more initial guess points

• Roots are found iteratively until some criteria are satisfied:

  • Accuracy
  • Number of iteration

• Bisection vs Newton