December 16, 2024

Solve a System of Nonlinear Equations in MATLAB


In this post, we will learn how to numerically solve systems of nonlinear equations using the MATLAB programming language. These types of problems appear in various engineering and physics disciplines. A video accompanying this post is given below.

Consider the following system of equations:

(1)   \begin{align*}2x_{1}-x_{2}-e^{-x_{1}}=0 \notag  \\-x_{1}+2x_{2}-e^{-x_{2}}=0 \end{align*}


where x_{1},x_{2}\in \mathbb{R} are uknown variables. Our task is simple: compute the solution of the above system of equations. This example is taken from the MATLAB explanation of the fsolve() function and can be found here.

In order to solve this system, we first need to define a MATLAB function that returns the value of the left-hand side of (1). The function is given below

function F=nonlinear_function(x)
    F=[2*x(1)-x(2)-exp(-x(1));
     -x(1)+2*x(2)-exp(-x(2))];
end

For a given value of the vector x=[x(1); x(2)], this function computes the left-hand side of the system (1). To use this function from another MATLAB script, we need to save it as a file. Save this function as the MATLAB entitled “nonlinear_function.m”. Furthermore, save this file in a folder and add this folder to a MATLAB path. This can be obtained by selecting “HOME->Set Path” in MATLAB.

After saving this function and adding it to the MATLAB path, we can call this function from any MATLAB script. Next, open a new MATLAB script, and type the following piece of code:

clear, pack, clc
% Possible algorithms 'trust-region-dogleg', 'trust-region', or 'levenberg-marquardt'
options = optimoptions(@fsolve,'Algorithm','trust-region','Display','iter','UseParallel',true,'OptimalityTolerance',1.0000e-8)

initial_guess=[0;0]

[solution] = fsolve(@(x)nonlinear_function(x),initial_guess,options);

F=[2*solution(1)-solution(2)-exp(-solution(1));
  -solution(1)+2*solution(2)-exp(-solution(2))];

A few comments are in order. We use the MATLAB function fsolve() to solve the nonlinear system of equations. On the code line 3 we set the solver options. We use the “trust-region-dogleg” algorithm. We set the ‘Display’ option to ‘iter’ since we want to monitor and display the solver progress. Since we want to speed up the computations, we set ‘UseParallel’ -> true (although for a small system that we are dealing with, the parallel option does not lead to any significant increase in the computation speed). Finally, with the option ‘OptimalityTolerance’ -> 1.0000e-8, we set the solution tolerance. On the code line 4, we choose an initial condition. Since MATLAB solves the nonlinear system using iterative methods, we need to initialize the solver with an initial_guess. Finally, on the code line 7, we solve the system. The code lines 9-10 are used to verify the solution. By executing this code, we obtain:


                                         Norm of      First-order 
 Iteration  Func-count     f(x)          step          optimality
     0          3               2                             2
     1          6       0.0226976       0.707107          0.171      
     2          9     3.40349e-06      0.0937779        0.00204      
     3         12     7.72066e-14     0.00117685       3.08e-07      
     4         15     3.19489e-29    1.77304e-07       6.26e-15      

Equation solved.

fsolve completed because the vector of function values is near zero
as measured by the selected value of the function tolerance, and
the problem appears regular as measured by the gradient.

<stopping criteria details>


solution =

    0.5671
    0.5671

The function fsolve(), returns the solution vector “solution”. solution(1) and solution(2) correspond to x_{1} and x_{2}, respectively. Finally, by executing “norm(F)”, we can obtain the Euclidian norm of the solution vector in order to check the accuracy. The result is 5.6523e-15, which means that the solution is computed with high accuracy.