Why are the solutions zeros for a linear equations ?

26 次查看（过去 30 天）
rong
rong 2014-9-17
评论： Matt J 2014-9-17
Hi everyone,
I have a Ax=b need to be solved. A is a 119309*119309 matrix, and b=zeros(119309,1). I used Jacobi and Gauss-seidel methods to solve the equations and the initial values for x are zeros. But all of the solutions for x are zeros, that is zeros(119309,1). I calculate the determinant of A and is 0, that means A is a singular matrix, and the sprank (A)<n(119309). I don't know why I can not get other non-zero solutions. I have been confused for this problem for a week and cannot find the reason. I will appreciate your help if you can help me. My codes in matlab are as follows:
Jacob:
function[x,k,index]=Jacobi(A,b,ep,it_max)
step=0
if nargin<4
it_max=100;
end
if nargin<3
ep=1e-10;
end
n=length(A);
k=0;
x=zeros(n,1);
y=zeros(n,1);
index=1;
while 1
for i=1:n
y(i)=b(i);
for j=1:n
if j~=i
y(i)=y(i)-A(i,j)*x(j);
end
end
if abs(A(i,i))<1e-10 | k==it_max
index=0;
return;
end
y(i)=y(i)/A(i,i);
end
if norm(y-x,inf)<ep
break;
end
x=y;
k=k+1;
step=step+1
end
[x,k,index]=Jacobi(A,b,1e-10,10)
Gauss-seidel:
function[v,sN,vChain]=GaussSeidell(A,b,x0,errorBound,maxSp)
step=0
error=inf;
vChain=zeros(15,151);
k=1;
fx0=x0;
L=-tril(A,-1);
U=-triu(A,1);
C=diag(diag(A));
b=zeros(151,1);
while error>=errorBound & step<maxSp
x0=inv(C)*(L+U)*x0+inv(C)*b;
vChain(k,:)=x0;
k=k+1;
error=norm(x0-fx0);
fx0=x0;
step=step+1;
e=1
end
v=x0;
sN=step
end
[v,sN,vChain]=GaussSeidell(A,b,x0,1e-6,100)
  1 个评论
Pierre Benoit
Pierre Benoit 2014-9-17
编辑：Pierre Benoit 2014-9-17
First, please format your code.
One way to see why a program's behaviour is wrong is to test a small example and see what it does step by step (with inserting breakpoints for example). You could start with a matrix A 2x2 like
A = [1 2; 1 2];
where one of the solution is
x = [-2 ; 1]
And from what I gathered of the methods you're using to solve this problem, you use a starting point and iterates, but here you start with an obvious solution which is 0, so maybe your code just stop there since you already found a solution. You may want to start with a different vector x.

请先登录，再进行评论。

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

采纳的回答

Pierre Benoit
Pierre Benoit 2014-9-17
The null function in Matlab achieve what you desire (for b=0) but I suppose this is homework related to better understand these algorithms rather than an efficient implementation to solve this.
As I stated in my comment, if you start with a solution, your algorithm will end in the first iteration, so you have to start by another vector. And you have to remember the fact that such an algorithm can never converge or converge quite slowly. Trials and errors may be necessary to find a good starting point.
Also, are you sure you can use these methods to solve Ax = 0 ? In dimension 2, you will always come back to the same vector every two iteration since det(A) = 0. I don't know for higher dimensions but all the information I found was about to solve Ax=b with b different from 0.
Remember that such equation (b=0) have infinite solutions for A singular, see Kernel for more informations.
  1 个评论
rong
rong 2014-9-17
Yes, I changed the initial solution x=0.0001 and tried for a smaller matrix equations. Then I got a different solution with 0. But when I applied on the original matrix A(119309*119309). Errors just displayed: Warning: "Matrix is singular to working precision." or Error"using inv Out of memory. Type HELP MEMORY for your options. "

请先登录，再进行评论。

更多回答（2 个）

Matt J
Matt J 2014-9-17
编辑：Matt J 2014-9-17
In every pass through the loop, you reset the state variables to some unchanging quantity. In the Jacobi, you have
y(i)=b(i);
while in the Gauss Seidel, you have
x0=inv(C)*(L+U)*x0+inv(C)*b;
The algorithm can't progress if you're always resetting...
rong
rong 2014-9-17
Does Jacobi or Gauss-seidel solve Ax=b if A is singular matrix? If not I tried x=pinv(A)*b, but it only works for full matrix. In my case, A is a sparse matrix, so I used QR to decompose A, but I am not sure if it is correct?
  1 个评论
Matt J
Matt J 2014-9-17
Why not just do A\b? For example,
>> A=sparse([ones(3,2),zeros(3,5)]); b=[1;1;1];
>> x=A\b
Warning: Rank deficient, rank = 1, tol = 7.691851e-14.
x =
1
0
0
0
0
0
0

请先登录，再进行评论。

类别

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

标签

Community Treasure Hunt

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

Start Hunting!

Translated by