最近遇到了一些需要求函数的雅可比矩阵的问题,例如上次发布的 blog:Matlab 最速下降法 实列及代码实现,里面也需要用到求解雅可比矩阵,那篇blog中的雅可比也是自己编写的代码,之前搜过自带的函数 (官方函数叫:jacobian(f,v)
1),由于官方的自带的没有自己比较想要的一些功能,那就自己动手这样写出来的函数可以达到自己想要的样子。
本篇介绍的雅可比矩阵代码是上次那篇的改进升级版,更加的全面。
Matlab 版本 2022b
- 句柄函数和符号表达式 (将在第二章内详细讲解)
- 求导函数 diff()
- symvar()函数
- 进阶版(可以选择性了解):varargout,varargin,nargin
对于因变量为一元函数情况下,即 y(x1,x2,x3,…,xm)y(x_1,x_2,x_3,\dots,x_m)y(x1,x2,x3,…,xm) ,它的雅可比矩阵就是各各个自变量的偏导数组成了导数向量 JF=(∂y∂x1,∂y∂x2,∂y∂x3,…,∂y∂xm)J_F=(\frac{\partial y}{\partial x_1},\frac{\partial y}{\partial x_2},\frac{\partial y}{\partial x_3},\dots,\frac{\partial y}{\partial x_m})JF=(∂x1∂y,∂x2∂y,∂x3∂y,…,∂xm∂y)
对于一函数由 mmm 维自变量 X=(x1,x2,…,xm)TX=(x_1,x_2,\dots,x_m)^TX=(x1,x2,…,xm)T 映射到 nnn 维的应变量上,即为 F(X)={f1(X),f2(X),…,fn(X)}TF(X)=\{f_1(X),f_2(X),\dots,f_n(X)\}^TF(X)={f1(X),f2(X),…,fn(X)}T ,它的雅可比矩阵为JF=[∂f1∂x1∂f1∂x2⋯∂f1∂xm∂f2∂x1∂f2∂x2⋯∂f2∂xm⋮⋮⋱⋮∂fn∂x1∂fn∂x2⋯∂fn∂xm]J_F=\left[\begin{array}{cccc} \frac{\partial f_{1}}{\partial x_{1}} &\frac{\partial f_{1}}{\partial x_{2}} &\cdots & \frac{\partial f_{1}}{\partial x_{m}} \\ \frac{\partial f_{2}}{\partial x_{1}} &\frac{\partial f_{2}}{\partial x_{2}} &\cdots & \frac{\partial f_{2}}{\partial x_{m}} \\ \vdots& \vdots & \ddots & \vdots \\ \frac{\partial f_{n}}{\partial x_{1}} & \frac{\partial f_{n}}{\partial x_{2}} &\cdots & \frac{\partial f_{n}}{\partial x_{m}} \end{array}\right]JF=∂x1∂f1∂x1∂f2⋮∂x1∂fn∂x2∂f1∂x2∂f2⋮∂x2∂fn⋯⋯⋱⋯∂xm∂f1∂xm∂f2⋮∂xm∂fn
只是简单的回顾一下雅可比矩阵的定义,如果想要深入的了解就可以去翻阅高数书
在平时书写函数表达式的时候是这样的如:y=sinxy=\sin xy=sinx,f(x,y)=yexf(x,y)=ye^xf(x,y)=yex,F(x1,x2)=[x12−x2x22−x1]F(x_1,x_2)=\begin{bmatrix} x_1^2-x_2\\ x_2^2-x_1 \end{bmatrix}F(x1,x2)=[x12−x2x22−x1]等等一些表达式。那如何在matlab中把它们表示出来,接下来讲两种表达方式:(1)符号表达式,(2)句柄函数。
这一部分细讲需要花一篇来写,这里就是先普及了解一下matlab中有这么个东西就可以了,之后有时间就详细讲这一块内容。
在数学运算中,运算的结果如果是一个数值,可以称这类运算为数值运算;如果运算结果为表达式,在MATLAB中称为符号运算,符号计算是对未赋值的符号对象 (可以是常数、变量、表达式)进行运算和处理。
简单来说吧(个人看法),就是定义这个变量是个符号字母,你可以把它当做数字进行各种运算,只是一种不知道数值得一个量。
举得例子是用matlab中的实时在线脚本编写
首先通过syms 变量1 变量2 ... 变量n
定义变量,举个例子
syms x y z
可以看到这样的结果
接下来编写一些简单得函数表达式 y=sinx+zy=\sin x +zy=sinx+z
syms x y z
y=sin(x)+z
看到结果是这样的
在Matlab 2022b中会有这么一个比较智能的功能,本人是直接从2020b升级到2022b不知道2021的版本有无这些功能,大家可以尝试这些功能非常有趣的。
令x=1x=1x=1代入上述式子应该是y=sin1+zy=\sin 1+zy=sin1+z,用到sub()函数
实现,代码实现:
syms x y z
y=sin(x)+z
y=subs(y,x,1)
结果为:
函数句柄是一种表示函数的 MATLAB® 数据类型。函数句柄的典型用法是将函数传递给另一个函数。例如,您可以将函数句柄用作基于某个值范围计算数学表达式的函数的输入参数2。
函数句柄可以表示命名函数或匿名函数。要创建函数句柄,需要使用 @
运算符。
例如,创建用于计算表达式 x2 – y2 的匿名函数的句柄:
f = @(x,y) (x.^2 - y.^2);
相对于上述符号函数只需要一行代码还是相当简单的,但是函数句柄就没有符号函数上面那些功能。
接下来令x=1,y=2x=1, y=2x=1,y=2 结果是f(1,2)=−3f(1,2)=-3f(1,2)=−3,只需要天加一句代码就能计算函数值,非常简单。
f = @(x,y) (x.^2 - y.^2);
f(1,2)
可以看到这个写法非常符合平时学数学将变量值代入得写法,感觉这点matlab做的非常的人性化
函数句柄与符号表达式都有自己的优缺点,这个需要根据自己体会使用才知道两者中那个会比较适合自己的在编写代码的过程带来便捷,直观。对我个人而言我程序会使用两者,比如在带入点坐标时会使用句柄函数,比如要展现表达式的时候我会更倾向于这个符号表达式。
sym()
函数即可。举个例子如下:f = @(x,y) (x.^2 - y.^2)
f1 = sym(f)
matlabFunction()
函数即可。举个例子如下:f = @(x,y) (x.^2 - y.^2)
f1 = sym(f)
f2 = matlabFunction(f1)
str2sym()
函数即可。举个例子如下:f = 'x.^2 - y.^2'
f1 = str2sym(f)
注意:
字符串不能直接转化为句柄函数的,步骤需要:字符串 ⟶\longrightarrow⟶ 符号表达式 ⟶\longrightarrow⟶ 函数句柄函数 | 作用 |
---|---|
str2sym() | 将字符串转化为符号表达式 |
matlabFunction() | 将符号表达式转化为函数句柄形式 |
subs() | 使用新变量替换符号表达式中的某些旧变量 |
char() | 将符号表达式转化为字符串 |
double() | 将sym格式转化为浮点型(比如对一个sym格式的式子求解然后求出来的结果是个sym格式,可以通过此函数转为浮点型)好用 |
isa(变量,‘sym’) | 判断变量是否为sym格式 |
symvar() | 确定表达式中的符号变量 |
diff() | 求导函数 |
[J,Jf,var,jx] = Jacobi(f,_)
: 此函数为雅可比函数,在编写函数过程使用了varargout
和varargin
,前者可以不固定函数输出个数,后者为不固定函数输入参数。
代码链接放到了文章最后,存放到了个人的GitHub中
输入:
f
可以是字符串形式 / 符号函数形式;f
,其中_
处可以有输入参数也可以没有参数。 [J,Jf,var] = Jacobi(f)
;x0
,函数使用为[J,Jf,var,jx] = Jacobi(f,x0)
。输出:
J
函数f
雅可比矩阵(符号表达式形式展现);Jf
函数f
雅可比矩阵(句柄函数形式展现);var
为函数f
中的变量(字符串形式展现,以元胞数组(cell)形式存放中,n × 1维);jx
为在某点x0x_0x0处的雅可比矩阵(double),当有输入x0
的话才会有输出jx
。本函数常用形式:
J = Jacobi(f)
[J,Jf] = Jacobi(f)
[J,Jf,var] = Jacobi(f)
[J,Jf,var,jx] = Jacobi(f,x0)
[~,Jf,var,jx] = Jacobi(f,x0)
,当不想某些参数出现时可以使用~
符号。[~,~,~,jx] = Jacobi(f,x0)
,如果只想得到最后jx
可以这样写。function varargout=Jacobi(f,varargin)
% 求函数表达式雅可比矩阵 J
% 输出:
% 输出1:符号表达式
% 输出2:句柄函数
% 输出3:变量
% 输出4:若输入点x0, 输出带入点后的值
% @Author
% Copyright© 2022.10.18 CSDN name: cugautozp[x,f]=fx(f);n=nargin(f); % 找到输入参数个数df=[]; for i =1:ndf1 = diff(f,x(i));df = [df,df1];end
% J=matlabFunction(df); varargout{1}=df; % 输出为符号表达式varargout{2}=matlabFunction(df); % 输出为句柄函数for i=1:length(x)s{i}=char(x(i));endvarargout{3}=s; % 输出变量if ~isempty(varargin)varargout{4}=Jx(df,s,varargin{1}); % 输出代入点后的值 end
endfunction [x,f]=fx(f)
% 将用字符串写的函数表达式转化为句柄函数 if ~isa(f,'sym') % 判断f是否为符号函数格式。if iscolumn(f)f=str2sym(f);elsef=str2sym(f');end end x=symvar(f); % 搜寻函数中的符号变量f=matlabFunction(f);
endfunction Jk=Jx(J,x,x0)
% 将点 x0 代入雅可比矩阵 J 中求值
% 输出格式为:矩阵值n=nargin(matlabFunction(J));if n==0Jk=double(J);elsea=symvar(J); % 找雅可比矩阵中的符号变量for i=1:length(a)s=char(a(i));idx(i) = find(strcmp(x,s));endJk = subs(J,a,x0(idx));Jk = double(Jk);end
end
举几个示例
f={'x^2+y^2+z^2','y^z*cos(x)','x^z*sin(y)','x^(y^z)'};
[J,Jf,var] = Jacobi(f)
输出结果:
2. 求以下函数形式的雅可比矩阵和在点(1,2,3)处的结果
F(X)=(zexyxyz)F(X) = \begin{pmatrix} ze^{x^y}\\ x\\ y\\ z \end{pmatrix}F(X)=zexyxyz
使用了符号函数形式写出 f
syms x y z
f=[z*exp(x^y);x;z;y]
x0=[1,2,3];
[J,Jf,var,jx] = Jacobi(f,x0)
输出结果:
官方的雅可比函数jacobian(f,v)
只能求出雅可比矩阵,无法求出某点处的雅可比矩阵值。其中V
为变量向量
功能比较少,使得使用时候有些功能不是很全,不能满足个人想使用的要求。
举个例子
syms x y z
f=[z*exp(x^y);x;z;y]
J = Jacobi(f)
j1 = jacobian(f,[x,y,z])
结果展示
相关代码文件:
GitHub:code
Jacobian matrix - MATLAB jacobian - MathWorks 中国 ↩︎
函数句柄 - MATLAB & Simulink - MathWorks 中国 ↩︎