How to fix **Warning: Matrix is close to singular or badly scaled. Results may be inaccurate. RCOND = 1.053110e-20.** ?

720 次查看（过去 30 天）
kaps
kaps 2022-5-18
评论： John D'Errico 2023-9-28
I am trying to solve the following PDE by using finite difference
For a uniform spacing $h$, I got the following equation,
for i =1,2,......,Nx+2 and j=1,2,...,Ny+2. After the implementation of the boundary conditions, I converted the system into the system of linear equations Au=f.
Now, I am trying to solve this system of linear equations, for certain value of $Nx$, I got the following warning;
**Warning: Matrix is close to singular or badly scaled. Results may be inaccurate. RCOND = 1.053110e-20.**
I believe this message is appearing because the matrix $A$ become singular after the implementation of the boundary conditions. The following a sub-code of my main code to solve this problem.
function [v1] = new_v1(w,Nx,Ny,dx,dy)
% -------------------------------------------------------
Iint = 1:Nx+1; % Indices of interior point in x direction
Jint = 1:Ny+2; % Indices of interior point in y direction
%---------------------------------------------
%assembly of the tridiagonal matrix(LHS)
sx = 1/(dx^2);
sy = 1/(dy^2);
e=ones(Nx+1,1);
T=spdiags([sx.*e,((-2*sx)+(-2*sy)).*e,sx.*e],-1:1,Nx+1,Nx+1);
T(1,Nx+1)= sx;
T(Nx+1,1)= sx;
D=spdiags(sy.*e,0,Nx+1,Nx+1);
A=blktridiag(T,D,D,Ny+2);
for i=1:Nx+1
for j=1:Nx+1
A(i,j)=(1/2).*A(i,j);
A((Nx+1)*(Ny+1)+i,(Nx+1)*(Ny+1)+j)=(1/2).*A((Nx+1)*(Ny+1)+i,(Nx+1)*(Ny+1)+j);
end
end
%---------------------------------------------------------------
%Solve the linear system
rhs = w ;
for i=1:Nx+1
rhs(i,1)=(1/2).*rhs(i,1);
rhs(i,Ny+2)=(1/2).*rhs(i,Ny+2);
end
%convert the rhs into column vector
F = reshape(rhs, (Nx+1)*(Ny+2),1);
uvec = A\F;
v1(Iint, Jint)= reshape(uvec, Nx+1,Ny+2);
end
  1 个评论
Sharmin Kibria
Sharmin Kibria 2022-5-23
Hi kaps,
I believe the warning is coming from the matrix inversion uvec = A\F;. These warnings are generally caused when the singular value of the matrix A becomes very close to zero. In that case the determinant of the matrix becomes very close to zero.
Sometimes these errors are fixed by modifying the singular values to ensure they do not becomes close to numerically zero.
[U,S,V] = svd(A) performs a singular value decomposition of matrix A, such that A = U*S*V'. Setting the singular values to (S+e) (e is a very small number) instead of S should resolve the issue.
Thanks
Sharmin

请先登录，再进行评论。

请先登录，再回答此问题。

回答（1 个）

Christine Tobler
Christine Tobler 2022-5-25
Yes, the matrix A becomes singular after applying the last two for-loops (which I think are the boundary conditions).
You can verify that by calling null(full(A)) on an example matrix (I used Nx = Ny = 10, dx = dy = 0.1). This showed that there is a null space of dimension one, and the vector in that null space had all elements of equal value.
That means that if you insert all values of x as being the same number, then A*x = 0 is always true. That means there is a degree of freedom left in this linear system.
Thinking about your initial equations, this makes sense: None of the three equations above tie the value of u to any specific number: The derivatives are specified, but we just know that u(0, y) = u(1, y) is required, which could be true for any function u that is shifted by a constant.
I would suggest just setting one of the values of x to an arbitrary value, and then rearranging the rest of this linear system based on that value.
  3 个评论
显示 1更早的评论
Christine Tobler
Christine Tobler 2022-5-27
Yes, when there is no unique solution, the linear system is singular, a warning is given, and the results are less reliable.
You can make the solution unique by adding an additional equation which sets any one of the values in the solution to zero (or any other value, basically this fixes the c). To keep the linear system square (which is cheaper to solve typically), you could remove one of the existing equations, since the system being singular means that one of these equations is redundant.
Or you could use the lsqminnorm function, which computes the solution x which has minimal 2-norm among all solutions. But keep in mind that this can get expensive as the size of the system increases.
John D'Errico
John D'Errico 2023-9-28
As Christine has said, that is exactly the source of your problem. You admit that you can essentially translate the problem by any constant value, and the solution is as good. This means you need to choose any arbitrary point, and force the solution to pass through that point. If you don't, then the problem becomes singular. And backslash won't solve a singular problem. (At least not in a way that will make you happy.)
Also as Christine has suggested, one solution is to use lsqminnorm, which will be expensive. But lsqminnorm actually does implicitly choose some constant for you, It just does so in a way that won't be terribly useful to you. Instead, the correct way to solve the problem is to indeed enforce that constant to take on some specific value. Just pick one corner of the domain, and set u(0,0)=0, for example as another boundary condition. That will immediately resolve the problem.
If you think about it, a singular system of equations in linear algebra typically means there are infinitely many equivalent solutions. The singularity just reflects back to you what you already know. But it is also what gives the solver a case of the hiccups, since it will not resolve that singularity for you.
This is a common issue faced by students in engineering courses. They might be told to build a truss model of a bridge that minimizes the potential energy of the structure. But unless they build in a boundary condition that properly fixes the structure at some location in space, the equations become singular. Essentially the entire structure can move as a whole anywhere in the universe, as long as it has the same fundamental shape.

请先登录，再进行评论。

类别

在 Help Center 和 File Exchange 中查找有关 Structural Mechanics 的更多信息

标签

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by