Welcome to our community

Be a part of something great, join today!

Output-Nan..is it wrong??

evinda

Well-known member
MHB Site Helper
Apr 13, 2013
3,720
Hi :) I have written a code in matlab to implement the Jacobi method.I calculated the error [tex] \left \| x_{k}-D \right \|[/tex] (where D is the real solution) of the last iteration for different values of [tex] n [/tex] .For [tex] n\geq 250 [/tex] ,I get this as result: NaN.
Does this mean that I have done something wrong??
 

Ackbach

Indicium Physicus
Staff member
Jan 26, 2012
4,193
Yes, there's probably something wrong if that's supposed to be an error. NaN means "Not a Number". MATLAB responds with that if the result is infinite, or if there's a category error like you ask it to do something with strings that you normally only do with numbers. Can you post your code?
 

evinda

Well-known member
MHB Site Helper
Apr 13, 2013
3,720
Yes, there's probably something wrong if that's supposed to be an error. NaN means "Not a Number". MATLAB responds with that if the result is infinite, or if there's a category error like you ask it to do something with strings that you normally only do with numbers. Can you post your code?
That is my code:
Code:
it=1;
x2=-inv(F)*W*x1+inv(F)*b
while ((it<=maxit)&&(norm(x2-x1)<e)&&(norm(b-A*x2))<e))
        error=norm(x-D);
        x1=x2;
        x2=-inv(F)*W*x1+inv(F)*b
        it=it+1;
end
I get Nan as result when I apply the Jacobi method at the Hilbert matrix,for [tex] n\geq 250 [/tex].Why does this happen???
 

Ackbach

Indicium Physicus
Staff member
Jan 26, 2012
4,193
There are a lot of variables here that I don't know what they are. Could you please provide a complete cheat sheet for all your variables? Ideally, you could have all their definitions. I especially want to know the dimensions of all variables (whether matrix or vector). Thanks!
 

evinda

Well-known member
MHB Site Helper
Apr 13, 2013
3,720
There are a lot of variables here that I don't know what they are. Could you please provide a complete cheat sheet for all your variables? Ideally, you could have all their definitions. I especially want to know the dimensions of all variables (whether matrix or vector). Thanks!
Code:
it=1;
x2=-inv(F)*W*x1+inv(F)*b
while ((it<=maxit)&&(norm(x2-x1)<e)&&(norm(b-A*x2))<e))
        error=norm(x-D);
        x1=x2;
        x2=-inv(F)*W*x1+inv(F)*b
        it=it+1;
end
F:the diagonal matrix of the Hilbert matrix
W:the Hilbert matrix - F
b:Ax=b
it:number of iteration
x2:approximation of solution
x1: previous approximation of solution
maxit:maximum number of iterations
e:a small number given from the user and
error:the error we are looking for
 

Ackbach

Indicium Physicus
Staff member
Jan 26, 2012
4,193
Ok, I have a few comments/thoughts:

1. You could probably speed this code up just a bit by setting Finv = inv(F), and referencing Finv instead of inv(F). As you are not changing F in the loop, it would be faster simply to recall the contents of memory somewhere than to evaluate the inverse of F multiple times in the loop.

2. Refering to your OP, is n the number of iterations? If so, could you please give me a plot of several errors on the y axis versus the corresponding number of iterations on the x axis? I'm wondering if you might not be getting an underflow. That is, 250 iterations of the Jacobi method might make your error so small that you run up against the machine epsilon, and the computer cannot represent an error that small.

3. There might be a logical error in your program: you have the line
error=norm(x-D);
Should you be saying, instead, that
error=norm(x1-D)? Or error=norm(x2-D)? (Probably the latter, due to first-time-through issues.)
I mean, does the computer know what x is? It looks to me like you only know what x1 and x2 are (x isn't in your list). That indicates to me that that it's probably a typo in your code.

4. The while statement has unbalanced parentheses. Try enclosing the last test in its own parentheses.

See if one of these things, especially 2 or 3, doesn't help.
 

evinda

Well-known member
MHB Site Helper
Apr 13, 2013
3,720
Ok, I have a few comments/thoughts:

1. You could probably speed this code up just a bit by setting Finv = inv(F), and referencing Finv instead of inv(F). As you are not changing F in the loop, it would be faster simply to recall the contents of memory somewhere than to evaluate the inverse of F multiple times in the loop.

2. Refering to your OP, is n the number of iterations? If so, could you please give me a plot of several errors on the y axis versus the corresponding number of iterations on the x axis? I'm wondering if you might not be getting an underflow. That is, 250 iterations of the Jacobi method might make your error so small that you run up against the machine epsilon, and the computer cannot represent an error that small.

3. There might be a logical error in your program: you have the line
error=norm(x-D);
Should you be saying, instead, that
error=norm(x1-D)? Or error=norm(x2-D)? (Probably the latter, due to first-time-through issues.)
I mean, does the computer know what x is? It looks to me like you only know what x1 and x2 are (x isn't in your list). That indicates to me that that it's probably a typo in your code.

4. The while statement has unbalanced parentheses. Try enclosing the last test in its own parentheses.

See if one of these things, especially 2 or 3, doesn't help.

1)Ok..And...could I also write it like that: x2=-F\W*x1+F\b ??

2)I wrote this code in a function and with n I mean the maxit I give each time when I call the function...The underflow,has it to do with the Hilbert matrix??Because,when I ran my code,using other matrices,I didn't get Nan as result... :confused:

3)Oh sorry,I accidentally wrote it like that. :eek: At the code I ran ,I wrote it like that: error=norm(x2-D) ..
 

Ackbach

Indicium Physicus
Staff member
Jan 26, 2012
4,193
I mean something more like this:

Code:
it=1;
Finv=inv(F);
x2=-Finv*W*x1+Finv*b
while ((it<=maxit)&&(norm(x2-x1)<e)&&(norm(b-A*x2)<e))
        error=norm(x2-D);
        x1=x2;
        x2=-Finv*W*x1+Finv*b
        it=it+1;
end
 

evinda

Well-known member
MHB Site Helper
Apr 13, 2013
3,720
I mean something more like this:

Code:
it=1;
Finv=inv(F);
x2=-Finv*W*x1+Finv*b
while ((it<=maxit)&&(norm(x2-x1)<e)&&(norm(b-A*x2)<e))
        error=norm(x2-D);
        x1=x2;
        x2=-Finv*W*x1+Finv*b
        it=it+1;
end
I wrote it like that right now..but I get the same result again(NaN)..Have I done maybe an other error?? :confused:
 

Ackbach

Indicium Physicus
Staff member
Jan 26, 2012
4,193
There's something rather puzzling about your code: you check whether x1 is close to x2, which you wouldn't necessarily expect early in your iterations. So I can easily imagine a scenario where you wouldn't even enter the loop. However, presumably the loop is where you want the iterations to get closer and closer to the final answer. I don't have the time right this second to think about the proper loop condition, but think carefully about what comes right after the 'while' keyword.
 

evinda

Well-known member
MHB Site Helper
Apr 13, 2013
3,720
There's something rather puzzling about your code: you check whether x1 is close to x2, which you wouldn't necessarily expect early in your iterations. So I can easily imagine a scenario where you wouldn't even enter the loop. However, presumably the loop is where you want the iterations to get closer and closer to the final answer. I don't have the time right this second to think about the proper loop condition, but think carefully about what comes right after the 'while' keyword.
I checked if,giving the Hilbert matrix with dimension 250,maxit=250,my code enters the loop,and it does...it makes 250 iterations...Have you an idea what I could have done wrong??
 

evinda

Well-known member
MHB Site Helper
Apr 13, 2013
3,720
Also,the method does not converge...Has it maybe to do something with that?? :eek: :confused:
 

evinda

Well-known member
MHB Site Helper
Apr 13, 2013
3,720
But if I want to find the error of the last iteration of the Gauss-Seidel method,using the Hilbert Matrix,I get a number and not Nan as result...And,in this case the method does also not converge...

Why does this happen???
 

evinda

Well-known member
MHB Site Helper
Apr 13, 2013
3,720
I know where the error is now...Using matlab,the determinant of the matrix equals to zero..But,what could I do,so that I have more precision of digits,to have a determinant [tex] \neq 0[/tex] ???
 

Ackbach

Indicium Physicus
Staff member
Jan 26, 2012
4,193
I know where the error is now...Using matlab,the determinant of the matrix equals to zero..But,what could I do,so that I have more precision of digits,to have a determinant [tex] \neq 0[/tex] ???
Having more precision of digits isn't necessarily going to help with the determinant being zero.

So... if the determinant is zero, then you most likely have infinite solutions with that particular matrix. Are you sure you are using the correct matrix? Can you double-check all the entries (I know there are a lot of them!)? If you're going to solve $Ax=b$, then having the right $A$ is definitely of the essence.
 

evinda

Well-known member
MHB Site Helper
Apr 13, 2013
3,720
Having more precision of digits isn't necessarily going to help with the determinant being zero.

So... if the determinant is zero, then you most likely have infinite solutions with that particular matrix. Are you sure you are using the correct matrix? Can you double-check all the entries (I know there are a lot of them!)? If you're going to solve $Ax=b$, then having the right $A$ is definitely of the essence.
How could I check this??? :confused: I use the command hilb(n)..
 

Ackbach

Indicium Physicus
Staff member
Jan 26, 2012
4,193
Hmm. Hilbert matrices are technically invertible, and hence have non-zero determinant. I wonder if there is a MATLAB setting you can change to use more precision. We are fast approaching the limits of my knowledge.
 

evinda

Well-known member
MHB Site Helper
Apr 13, 2013
3,720
Hmm. Hilbert matrices are technically invertible, and hence have non-zero determinant. I wonder if there is a MATLAB setting you can change to use more precision. We are fast approaching the limits of my knowledge.
Where could I find such a command?? :confused:
 

evinda

Well-known member
MHB Site Helper
Apr 13, 2013
3,720
And...something else...is this right that both of the methods,the Jacobi and Gauss-Seidel method,do not converge for the Hilbert matrix???Because,no matter with dimension and which small number ε I give,the methods do not converge...
 

Opalg

MHB Oldtimer
Staff member
Feb 7, 2012
2,708
I know where the error is now...Using matlab,the determinant of the matrix equals to zero..But,what could I do,so that I have more precision of digits,to have a determinant [tex] \neq 0[/tex] ???
The determinant of the Hilbert matrix is not zero, but it tends to zero exponentially fast as the size of the matrix increases. I know nothing about matlab, but I doubt whether it could possibly have enough precision to distinguish $\det H$ from zero for a $250\times 250$ Hilbert matrix. As is mentioned here, the Hilbert matrix is notoriously ill-conditioned for numerical computation.
 

evinda

Well-known member
MHB Site Helper
Apr 13, 2013
3,720
The determinant of the Hilbert matrix is not zero, but it tends to zero exponentially fast as the size of the matrix increases. I know nothing about matlab, but I doubt whether it could possibly have enough precision to distinguish $\det H$ from zero for a $250\times 250$ Hilbert matrix. As is mentioned here, the Hilbert matrix is notoriously ill-conditioned for numerical computation.
I wrote the command det(H) and the result is that it equals to zero :confused:
 

evinda

Well-known member
MHB Site Helper
Apr 13, 2013
3,720
I have searched in Google and found this:
Hilbert matrices are notoriously ill-conditioned: the computation of the determinant is subject to severe cancellation effects. The following results, both with HardwareFloats as well as with SoftwareFloats, are marred by numerical roundoff:

A := linalg::hilbert(15):
float(numeric::det(A, Symbolic)),
numeric::det(A, HardwareFloats),
numeric::det(A, SoftwareFloats)

How could I apply this,to check the result??Because I get an error message,when I write it...
 

evinda

Well-known member
MHB Site Helper
Apr 13, 2013
3,720
The determinant of the Hilbert matrix is not zero, but it tends to zero exponentially fast as the size of the matrix increases. I know nothing about matlab, but I doubt whether it could possibly have enough precision to distinguish $\det H$ from zero for a $250\times 250$ Hilbert matrix. As is mentioned here, the Hilbert matrix is notoriously ill-conditioned for numerical computation.
Hasn't it maybe to do with the determinant,but with the condition number of the matrix?:confused:
 

evinda

Well-known member
MHB Site Helper
Apr 13, 2013
3,720
Anyway... :) I also calculated the spectral radius of the hilbert matrix with dimension 250.Is this right that it is 217.3320,using the Jacobi method?
 

Opalg

MHB Oldtimer
Staff member
Feb 7, 2012
2,708
Anyway... :) I also calculated the spectral radius of the hilbert matrix with dimension 250. Is this right that it is 217.3320,using the Jacobi method?
The Hilbert matrices are positive, so the spectral radius is equal to the norm, which is $\leqslant \pi$. Asymptotically, as the size of the matrix increases, the norm converges to $\pi.$ See Man-Duen Choi's paper "Tricks or Treats with the Hilbert Matrix".