Why are the solutions zeros for a linear equations ?

18 views (last 30 days)
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 Comment
Pierre Benoit
Pierre Benoit on 17 Sep 2014
Edited: Pierre Benoit on 17 Sep 2014
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.

Sign in to comment.

Accepted Answer

Pierre Benoit
Pierre Benoit on 17 Sep 2014
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 Comment
rong
rong on 17 Sep 2014
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. "

Sign in to comment.

More Answers (2)

Matt J
Matt J on 17 Sep 2014
Edited: Matt J on 17 Sep 2014
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 on 17 Sep 2014
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 Comment
Matt J
Matt J on 17 Sep 2014
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

Sign in to comment.

Categories

Find more on Mathematics and Optimization in Help Center and File Exchange

Community Treasure Hunt

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

Start Hunting!