MATLAB 方程式求根


前言

b站课程《MATLAB教程_台大郭彦甫(14课)》学习记录


一、Problem Statement

在这里插入图片描述

二、Symbolic Root Finding Approach 符号求根法

Performing mathematics on symbols, NOT numbers
在符号上进行数学运算,而不是数字
The symbols math are performed using “symbolic variables“
符号数学是使用“符号变量”来执行的
Use sym or syms to create symbolic variable

>> syms x
>> x + x + x
ans =
3*x
>> (x + x + x)/4
ans =
(3*x)/4
>> x = sym('x');
>> x + x + x
ans =
3*x
>> (x + x + x)/4 
ans =
(3*x)/4
>> y = x^2 - 2*x -8
y =
x^2 - 2*x - 8
y =
x^2 - 2*x - 8

三、Symbolic Root Finding: solve()

Function solve finds roots for equations
函数解是方程的解
𝑦 = 𝑥 ∙ sin(𝑥) − 𝑥 = 0

syms x
solve(x*sin(x)-x, x)

or

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

在这里插入图片描述

>> syms x
y = (cos(x))^2 - (sin(x))^2
solve(y,x)
y =
cos(x)^2 - sin(x)^2
ans =
pi/4

在这里插入图片描述

>> clear
syms x
y  = (cos(x))^2 + (sin(x))^2
solve(y,x)
y =
cos(x)^2 + sin(x)^2
ans =
Empty sym: 0-by-1

四、Solving Multiple Equations 解决多元方程

Solve this equation using symbolic approach:
用符号方法解此方程
在这里插入图片描述

>> clear
syms x y
eq1 = x - 2*y - 5;
eq2 = x + y - 6;
A = solve(eq1,eq2,x,y)
A = 
  包含以下字段的 struct:
    x: [1×1 sym]
    y: [1×1 sym]
>> A.x
ans =
17/3
>> A.y
ans =
1/3

五、Solving Equations Expressed in Symbols 解符号方程

What if we are given a function expressed in symbols?
在这里插入图片描述

>> clear
syms x a b
solve(a*x^2-b)
ans =
 b^(1/2)/a^(1/2)
-b^(1/2)/a^(1/2)

𝑥 is always the first choice to be solved
What if one wants to express 𝑏 in terms of 𝑎 and 𝑥?

>> clear
syms x a b
solve(a*x^2-b,b)
ans =
a*x^2

Exercise

Solve this equation for 𝑥 using symbolic approach
在这里插入图片描述

>> clear
syms x y a b r
t = (x - a)^2 + (y - b)^2 - r^2;
solve(t,x)
ans =
a + (b + r - y)^(1/2)*(r - b + y)^(1/2)
a - (b + r - y)^(1/2)*(r - b + y)^(1/2)

Find the matrix inverse using symbolic approach
在这里插入图片描述

>> clear
syms a b c d
A = [a b;c d];
inv(A)
ans =
[ d/(a*d - b*c), -b/(a*d - b*c)]
[-c/(a*d - b*c),  a/(a*d - b*c)]

六、Symbolic Differentiation: diff() 符号微分法

Calculate the derivative of a symbolic function:
计算一个符号函数的导数:
在这里插入图片描述

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

(直接将微分求了出来)

Exercise:

在这里插入图片描述

>> clear
syms x
f(x) = (exp(x^2))/(x^3 - x + 3);
fprime = diff(f(x))
fprime =
(2*x*exp(x^2))/(x^3 - x + 3) - (exp(x^2)*(3*x^2 - 1))/(x^3 - x + 3)^2
>> clear
syms x y
g(x) = (x^2 + x*y -1)/(y^3 + x +3);
gprime = diff(g(x),x)
gprime =
(2*x + y)/(y^3 + x + 3) - (x^2 + y*x - 1)/(y^3 + x + 3)^2

七、Symbolic Integration:

Calculate the integral of a symbolic function:
在这里插入图片描述

>> clear
syms x; 
y = x^2*exp(x);
z = int(y)
z =
exp(x)*(x^2 - 2*x + 2)

未规定常数时,会自动取一个常数

>> clear
syms x; 
y = x^2*exp(x);
z = int(y);
z = z - subs(z,x,0)
z =
exp(x)*(x^2 - 2*x + 2) - 2
>> clear
syms x
y = (x^2 - x + 1)/(x + 3);
z = int(y,0,10) 
z =
log(302875106592253/1594323) + 10
>> double(z)
ans =
   29.0624

八、Symbolic vs. Numeric 符号VS数字

在这里插入图片描述

九、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
先新建一个空白.m文件,先将下面的代码进行保存

function [y] = xy_plot(input,x)
% 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

再在另一个.m文件中进行运行

Try: xy_plot(@sin,0:0.01:2*pi);

在这里插入图片描述

十、Using Function Handles 使用函数处理

xy_plot(@sin,0:0.01:2*pi);

在这里插入图片描述

xy_plot(@atan,0:0.01:2*pi);

在这里插入图片描述

十一、fsolve()

A numeric root solver
一个数值根解计算器
For example, solve this equation:
𝑓(𝑥) = 1.2𝑥 + 0.3 + 𝑥 ∙ sin(𝑥)

>> f2 = @(x) (1.2*x+0.3+x*sin(x));
fsolve(f2,0)
ans =
   -0.3500

Exercise

在这里插入图片描述

>> clear
syms x
f = @(x) ([2*x(1) - x(2) - exp(-x(1));-x(1) + 2*x(2) - exp(-x(2))]);
t = fsolve(f,[-5,-5]);
x = t(1)
y = t(2)
x =
    0.5671
y =
    0.5671

十二、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)
f =
  包含以下值的 function_handle:
    @(x)x.^2
正在退出 fzero: 将终止搜索包含符号变化的区间
 因为在搜索期间遇到 NaN 或 Inf 函数值。
(-1.37296e+154 处的函数值为 Inf。)
请检查函数或使用其他起始值重试。
ans =
   NaN

(fzero只有在穿过x轴才可以解出来,相当于二分法)

>> f=@(x)x.^2
fsolve(f,0.1)
f =
  包含以下值的 function_handle:
    @(x)x.^2
Equation solved at initial point 
fsolve completed because the vector of function values at the initial point
is near zero as measured by the value of the function tolerance, and
the problem appears regular as measured by the gradient.
<stopping criteria details>
ans =
     0.0063

Options:

f=@(x)x.^2
options=optimset('MaxIter',1e3,'TolFun',1e-10);  【1e3为迭代次数1e-10为误差】
fsolve(f,0.1,options)
fzero(f,0.1,options)
>> f=@(x)x.^2
options=optimset('MaxIter',1e3,'TolFun',1e-10);
fsolve(f,0.1,options)
f =
  包含以下值的 function_handle:
    @(x)x.^2
Equation solved.
fsolve completed because the vector of function values is near zero
as measured by the value of the function tolerance, and
the problem appears regular as measured by the gradient.
<stopping criteria details>
ans =
   1.9532e-04

更接近0了

>> f=@(x)x.^2
options=optimset('MaxIter',1e3,'TolFun',1e-10);
%fsolve(f,0.1,options)
fzero(f,0.1,options)
f =
  包含以下值的 function_handle:
    @(x)x.^2
正在退出 fzero: 将终止搜索包含符号变化的区间
 因为在搜索期间遇到 NaN 或 Inf 函数值。
(-1.37296e+154 处的函数值为 Inf。)
请检查函数或使用其他起始值重试。
ans =
   NaN

十三、Finding Roots of Polynomials: roots()

Find the roots of this polynomial:

在这里插入图片描述

>> roots([1 -3.5 2.75 2.125 -3.875 1.25])
ans =
   2.0000 + 0.0000i
  -1.0000 + 0.0000i
   1.0000 + 0.5000i
   1.0000 - 0.5000i
   0.5000 + 0.0000i

roots only works for polynomials
根只对多项式有效
Find the roots of the polynomial:
在这里插入图片描述

>> roots([1 -6 -12 81])
ans =
   -3.5969
    5.5097
    4.0872

十四、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:

  1. Bracketing methods (e.g., bisection method) 交叉法
    Start with an interval that contains the root
    从包含根的间隔开始
  2. Open methods (e.g., Newton-Raphson method)开型法
    Start with one or more initial guess points
    从一个或多个初始猜测点开始

Roots are found iteratively until some criteria are satisfied:
迭代查找根,直到满足某些条件:

  1. Accuracy
  2. Number of iteration

十六、Bisection Method (Bracketing) 二分法

在这里插入图片描述
Algorithm: 运算法则
Loop 循环

  1. 𝑟 = (𝑙 + 𝑢)/2
  2. If 𝑓(𝑟) ∙ 𝑓(𝑢) < 0 then new interval [𝑟, 𝑢]
    If 𝑓(𝑙) ∙ 𝑓(𝑟) < 0 then new interval [𝑙, 𝑟]
    End

十七、Bisection Algorithm Flowchart 二分算法流程图

在这里插入图片描述

十八、Newton-Raphson Method (Open) 牛顿迭代法

在这里插入图片描述

十九、Newton-Raphson Algorithm Flowchart 牛顿迭代算法流程图

在这里插入图片描述

二十、Bisection vs. Newton-Raphson 二分法和牛顿迭代法

在这里插入图片描述

二十一、Recursive Functions 递归函数

在这里插入图片描述

二十二、Factorial Recursive Function 阶乘递归函数

The function includes a recursive case and a base case
该函数包括递归情况和基本情况
The function stops when it reaches the base case
函数在到达基本情况时停止

function output = fact(n)
% fact recursively finds n!
if n==1
    output = 1;
else
    output = n * fact(n-1);
end
end

在这里插入图片描述


总结

继续加油吧。

  • 4
    点赞
  • 55
    收藏
    觉得还不错? 一键收藏
  • 0
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值