Scaling dependence using equilibrate function

I want to use the function equilibrate, but I don't understand the following written in Rescaling to Solve a Linear System.

"equilibrate is most useful when the scales of the b and x vectors in the original system x = A\b are irrelevant. However, if the scales of b and x are relevant, then using equilibrate to rescale A only for the linear system solve is not recommended. The obtained solution does not generally yield a small residual for the original system, even if expressed using the original basis vectors."

Could someone explain the meaning of this with an example?

[Umar]

Hi @Nicholas ,

Addressing your query regarding, “I don't understand the following written in Rescaling to Solve a Linear System, equilibrate is most useful when the scales of the b and x vectors in the original system x = A\b are irrelevant. However, if the scales of b and x are relevant, then using equilibrate to rescale A only for the linear system solve is not recommended. The obtained solution does not generally yield a small residual for the original system, even if expressed using the original basis vectors." Could someone explain the meaning of this with an example?

Meaning

In layman’s terms, it helps improve the accuracy and efficiency of solving linear systems by rescaling the basis vectors in the matrix A. This can be useful when the scales of the vectors b and x in the original system x = A\b are not important.

An example to illustrate it

So, I created a matrix A with significantly different scales along its diagonal elements, then solve the original system x = A\b and equilibrate matrix A using the equilibrate function. Finally, solve the equilibrated system and compare the original solution x with the equilibrated solution xeq. So, you can observe now, how the use of equilibrate impacts the solution of the linear system and understand the implications of rescaling when the scales of b and x are relevant. Please see attached results.

Hope this helps, please let me know if you have any further questions.

[Nicholas]

Thanks for the explanation!

So with this in mind, when would we have a system of equations to solve where the scale is irrelevant? Maybe I have a hard time thinking about this because I come from an engineering background where units are almost always included.

[Umar]

Hi @Nicholas,

That is a very good question. Let me consider an example where you create a matrix A with varying scales in its diagonal elements. By solving the original system x = A\b and then equilibrating matrix A using the equilibrate function, you can compare the solutions x and xeq to observe the impact of rescaling on the linear system's solution elucidating the implications of rescaling when the scales of b and x are relevant. Please bear in mind that rescaling in linear systems, as mentioned in the context of equilibration, becomes crucial when the scales of the vectors b and x in the original system x = A\b are insignificant. Hope this helps.

Answer 1

[praguna manvi]

Hi @Nicholas,

Equilibration is a technique used to rescale the matrix A to improve numerical stability and conditioning while solving a Linear system. However, if the scales of b and x are important (i.e., they have physical or meaningful units that must be preserved), then rescaling (A) alone can lead to inaccurate results. This is because the rescaled system may not yield a solution that fits well with the original scales of b and x.

Below is an example demonstrating how the residuals can increase for a larger values of b using “equilibration”:

    n = 10; % Size of the matrix
    A = diag(linspace(1e-5, 1e5, n)); % Diagonal matrix with values ranging from 1e-5 to 1e5
    b = linspace(1e7, 1e15, n)'; % Vector with values ranging from 1e7 to 1e15
    % Use the equilibrate function to get P, R, and C
    [P, R, C] = equilibrate(A);
    % Rescale the matrix A and vector b
    A_eq = R * P * A * C;
    b_eq = R * P * b;
    % Solve the rescaled system
    x_eq = A_eq \ b_eq;
    % Transform the solution back to the original scale
    x_original_scale = C * x_eq;
    % Compute the residual
    residual = norm(A * x_original_scale - b);
    x = A \ b;
    disp("Residual Without equilibrate: " + norm(A * x - b))

Residual Without equilibrate: 0.0625

    disp("Residual With equilibrate: " + residual);

Residual With equilibrate: 0.28641