Matlab 自编雅可比矩阵 (jacobi) 函数与官方的Jacobian matrix（雅可比矩阵）函数对比及创新

0 引言

  最近遇到了一些需要求函数的雅可比矩阵的问题，例如上次发布的 blog：Matlab 最速下降法 实列及代码实现，里面也需要用到求解雅可比矩阵，那篇blog中的雅可比也是自己编写的代码，之前搜过自带的函数 (官方函数叫：jacobian(f,v)¹)，由于官方的自带的没有自己比较想要的一些功能，那就自己动手这样写出来的函数可以达到自己想要的样子。
  本篇介绍的雅可比矩阵代码是上次那篇的改进升级版，更加的全面。

Matlab 版本 2022b

1. 句柄函数和符号表达式 (将在第二章内详细讲解)
2. 求导函数 diff()
3. symvar()函数
4. 进阶版(可以选择性了解)：varargout，varargin，nargin

1 雅可比矩阵

  对于因变量为一元函数情况下，即 $y(x_1,x_2,x_3,\dots ,x_m)$ ，它的雅可比矩阵就是各各个自变量的偏导数组成了导数向量 $$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})$$
  对于一函数由 $m$ 维自变量 $X=(x_1,x_2,\dots ,x_m{)}^T$ 映射到 $n$ 维的应变量上，即为 F ( X ) = { f 1 ( X ) , f 2 ( X ) , … , f n ( X ) } T F(X)=\{f_1(X),f_2(X),\dots,f_n(X)\}^T F(X)={f1​(X),f2​(X),…,fn​(X)}T ，它的雅可比矩阵为 J F = [ ∂ f 1 ∂ x 1 ∂ f 1 ∂ x 2 ⋯ ∂ f 1 ∂ x m ∂ f 2 ∂ x 1 ∂ f 2 ∂ x 2 ⋯ ∂ f 2 ∂ x m ⋮ ⋮ ⋱ ⋮ ∂ f n ∂ x 1 ∂ f n ∂ x 2 ⋯ ∂ f n ∂ x m ] 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​​​ ​

只是简单的回顾一下雅可比矩阵的定义，如果想要深入的了解就可以去翻阅高数书

2 matlab中函数表达式两种方法

  在平时书写函数表达式的时候是这样的如：$y=\sin x$，$f(x,y)=y{e}^x$，$F(x_1,x_2)=\left[\begin{array}{c}{x}_1^2-x_2\\ x_2^2-x_1\end{array}\right]$等等一些表达式。那如何在matlab中把它们表示出来，接下来讲两种表达方式：（1）符号表达式，（2）句柄函数。

这一部分细讲需要花一篇来写，这里就是先普及了解一下matlab中有这么个东西就可以了，之后有时间就详细讲这一块内容。

2.1 符号表达式

  在数学运算中，运算的结果如果是一个数值，可以称这类运算为数值运算；如果运算结果为表达式，在MATLAB中称为符号运算，符号计算是对未赋值的符号对象 (可以是常数、变量、表达式)进行运算和处理。

简单来说吧（个人看法），就是定义这个变量是个符号字母，你可以把它当做数字进行各种运算，只是一种不知道数值得一个量。
举得例子是用matlab中的实时在线脚本编写

首先通过syms 变量1 变量2 ... 变量n定义变量，举个例子

syms x y z

可以看到这样的结果

接下来编写一些简单得函数表达式 $y=\sin x+z$

syms x y z
y=sin(x)+z

看到结果是这样的

在Matlab 2022b中会有这么一个比较智能的功能，本人是直接从2020b升级到2022b不知道2021的版本有无这些功能，大家可以尝试这些功能非常有趣的。

令$x=1$代入上述式子应该是$y=\sin 1+z$，用到sub()函数实现，代码实现：

syms x y z
y=sin(x)+z
y=subs(y,x,1)

结果为：

2.2 函数句柄

  函数句柄是一种表示函数的 MATLAB® 数据类型。函数句柄的典型用法是将函数传递给另一个函数。例如，您可以将函数句柄用作基于某个值范围计算数学表达式的函数的输入参数²。

函数句柄可以表示命名函数或匿名函数。要创建函数句柄，需要使用 @ 运算符。
例如，创建用于计算表达式 x2 – y2 的匿名函数的句柄：

f = @(x,y) (x.^2 - y.^2);

相对于上述符号函数只需要一行代码还是相当简单的，但是函数句柄就没有符号函数上面那些功能。

接下来令$x=1,y=2$ 结果是$f(1,2)=-3$，只需要天加一句代码就能计算函数值，非常简单。

f = @(x,y) (x.^2 - y.^2);
f(1，2)

可以看到这个写法非常符合平时学数学将变量值代入得写法，感觉这点matlab做的非常的人性化

2.3 函数句柄与符号表达式相互转化

  函数句柄与符号表达式都有自己的优缺点，这个需要根据自己体会使用才知道两者中那个会比较适合自己的在编写代码的过程带来便捷，直观。对我个人而言我程序会使用两者，比如在带入点坐标时会使用句柄函数，比如要展现表达式的时候我会更倾向于这个符号表达式。

• 函数句柄 $⟶$ 符号表达式：直接使用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)

• 注意： 字符串不能直接转化为句柄函数的，步骤需要：字符串 $⟶$ 符号表达式 $⟶$ 函数句柄

2.4 常会用到的一些函数

3 自编代码

[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为在某点$x_0$处的雅可比矩阵（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:n
        df1 = 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));
    end
    varargout{3}=s;                        % 输出变量
    if ~isempty(varargin)
        varargout{4}=Jx(df,s,varargin{1}); % 输出代入点后的值 
    end
end

function [x,f]=fx(f)
% 将用字符串写的函数表达式转化为句柄函数    
    if  ~isa(f,'sym')            % 判断f是否为符号函数格式。
        if iscolumn(f)
            f=str2sym(f);
        else
            f=str2sym(f');
        end                 
    end 
    x=symvar(f);             % 搜寻函数中的符号变量
    f=matlabFunction(f);
end

function Jk=Jx(J,x,x0)
% 将点 x0 代入雅可比矩阵 J 中求值
% 输出格式为：矩阵值
    n=nargin(matlabFunction(J));
    if n==0
        Jk=double(J);
    else
        a=symvar(J);  % 找雅可比矩阵中的符号变量
        for i=1:length(a)
            s=char(a(i));
            idx(i) = find(strcmp(x,s));
        end
        Jk = subs(J,a,x0(idx));
        Jk = double(Jk);
    end   
end

举几个示例

1. 求以下函数形式的雅可比矩阵
   F ( x , y , z ) = ( x 2 + y 2 + z 2 y z c o s ( x ) x z s i n ( y ) x y z ) F(x,y,z) = \begin{pmatrix} x^2+y^2+z^2\\ y^zcos(x)\\ x^zsin(y)\\ x^{y^z} \end{pmatrix} F(x,y,z)= ​x2+y2+z2yzcos(x)xzsin(y)xyz​ ​
   使用字符串形式作为函数 f 输入

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 ) = ( z e x y x y z ) 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)

输出结果：

4 官方函数

官方的雅可比函数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])

结果展示

5 参考文献

相关代码文件：
GitHub：code

到了这里，关于Matlab 自编雅可比矩阵 (jacobi) 函数与官方的Jacobian matrix（雅可比矩阵）函数对比及创新的文章就介绍完了。如果您还想了解更多内容，请在右上角搜索TOY模板网以前的文章或继续浏览下面的相关文章，希望大家以后多多支持TOY模板网！