matlab - 将 Newton Raphson 和 Gauss Seidel 用于非线性系统

问题描述

我需要帮助弄清楚如何结合 Newton Raphson 和 Gauss Seidel 方法来求解 Matlab 中的非线性方程组。这是我的逻辑：

通过求方程的导数将系统线性化并放入矩阵。转移系统，所以我们有[偏导数矩阵][H] = [原始系统]，其中 H 是我要解决的问题。H 是步长。
首先将我对 x 的初始猜测插入到第一个和最后一个矩阵中，然后每次求解 H 时，我都会将该步长添加到 x 的先前猜测中以创建一个新的 x。然后我将插入更新后的 x，找到另一个 h，然后继续该过程。
每次迭代在新旧猜测之间都有一个错误。当误差小于容差时，报告值。

我之前曾单独编写过 Newton Raphson 和 Gauss Seidel 方法。我的主要困惑来自如何准确地将其设置在一起并使用两个不同的变量。首先，我想知道我的逻辑是否正确，并且此设置的伪代码将非常有帮助。

我尝试对此进行编码，并且我相信矩阵是正确的，并且我认为 Gauss Seidel 函数是正确的，但是我对如何正确更新猜测感到非常困惑。现在，一旦我摆脱了错误，它就会输出 0 的迭代。我被告知对线性方程进行 [0,0] 的猜测，对非线性方程进行 [2,2] 的猜测。基本上，我不确定如何将 Newton Raphson 融入其中。我之前已经对 NR 方法进行了编码，但只使用了 1 个变量，所以我不明白如何将此逻辑转移到此代码中。这是我单独的 Newton Raphson 代码：

标签： matlab

您走在正确的轨道上，但我们可以在您的代码中改进以下几点：

首先，请注意 Newton-Raphson 方法和 Gauss-Seidel 方法都是求解方程组的数值方法。不惜一切代价避免使用符号变量。它不仅速度慢得多，而且不能解决所有的数值问题。
使用和滥用 MATLAB 语法来处理矩阵和向量。

您的问题是耦合 Newton-Raphson 和 Gauss-Seidel 求解器。我不会解决整个问题，正如您所说，您的 Gauss-Seidel 似乎正在工作。

我还将稍微更改您提供的一些符号，以便您可以完全了解 MATLAB 易于使用的语法。

非线性方程组：

所提供的非线性方程组可以用向量形式写成：

function fx = F(x)
    fx = zeros(2,1); % pre-allocate to create a column-vector
    fx(1) = x(1)-x(2)+1;
    fx(2) = x(1)^2+x(2)^2-4;
end

这将创建列向量函数F作为变量x1(x) 和x2(y) 的函数。

Newton-Raphson 算法需要对 Jacobian 进行评估。您需要计算这些导数。功能很简单，您可以手动完成。创建一个函数来计算雅可比矩阵作为变量的函数。

function J = jacob(x)
    J = zeros(2,2);
    J(1,1) = 1; % df/dx
    J(1,2) = -1; % df/dy
    J(2,1) = 2*x(1); % dg/dx
    J(2,2) = 2*x(2); % dgdy
end

您对 Newton-Raphson 的逻辑是正确的。让我们用它来编写我们的 Newton-Raphson：

function x = newton(fun,jacobian,x0)
    % This function uses the Newton-Raphson method to find a root for the
    % nonlinear system of equations F, with jacobian J, given an initial
    % estimate x0.

    % Here we take advantage of MATLAB syntax to easily implement our
    % routine. fun and jacobian are function handles to the equations and 
    % its derivatives. So it is computationally much more efficient to use
    % it instead of working with symbolic variables.

    % Following the steps you provided and the Newton-Raphson for one variable:
    TOL = 0.00001;
    MAXITER = 100;

    err = 1;
    iter = 0;
    x = x0;
    while err>TOL && iter<MAXITER % repeat until error is smaller than tolerance
        % step 1: linearize the system
        fx = fun(x);
        J = jacobian(x);
        % notice that the system is already in matrix form.

        % ********************************************* %
        % step 1.2: solve for stepsize H              * %
        H = -J\fx; % SOLVE LINEAR SYSTEM OF EQUATIONS * %
        % ********************************************* %

        % step 2: add stepsize to the previous guess of x
        x = x+H;

        % step 3: calculate error
        err = norm(H);
    end
end

上述函数使用 Newton-Raphson 方法求解非线性方程组的零点。让我们检查它是否在我们的主程序上运行：

clc; clear; close all force;

x = newton(@F,@jacob,[2;2]) % implemented Newton-Raphson
x0 = fsolve(@F,[2;2]) % Matlab built-in function

norm(x-x0) % compare both solutions

ans =

   8.9902e-10

我们的函数与 MATLAB 的内置函数之间的差异非常小。看来该功能有效。

请注意上面突出显示的步骤 1.2。这是重要的一步。在这一步中，我使用 MATLAB 反冲算子来求解线性系统 Ax=b。以下语句具有相同的功能（求解线性方程组）：

x = A\B
x = mldivide(A,B)

如果您必须使用 Gauss-Seidel 方法来求解线性方程组，我将把这些修改留给您去做。在步骤 1.2 中，您将更改 Gauss-Seidel 函数的反冲算子：

H = -J\fx;
H = gaussSeidel(-J,fx);

并正确编码您的gaussSeidel函数（在函数中转换您的代码）：

function x = gaussSeidel(A,b)
    % Solves the system of linear equations Ax=b by the Gauss-Seidel
    % method.
    %
    % ...
    %
end

这应该够了吧。

matlab - 将 Newton Raphson 和 Gauss Seidel 用于非线性系统

问题描述

解决方案

推荐阅读