# Output-Nan..is it wrong??

#### evinda

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

#### Ackbach

##### Indicium Physicus
Staff member
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
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 $$n\geq 250$$.Why does this happen???

#### Ackbach

##### Indicium Physicus
Staff member
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
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
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
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...

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

#### Ackbach

##### Indicium Physicus
Staff member
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
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??

#### Ackbach

##### Indicium Physicus
Staff member
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
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
Also,the method does not converge...Has it maybe to do something with that??

#### evinda

##### Well-known member
MHB Site Helper
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
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 $$\neq 0$$ ???

#### Ackbach

##### Indicium Physicus
Staff member
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 $$\neq 0$$ ???
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
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??? I use the command hilb(n)..

#### Ackbach

##### Indicium Physicus
Staff member
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
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??

#### evinda

##### Well-known member
MHB Site Helper
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
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 $$\neq 0$$ ???
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
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

#### evinda

##### Well-known member
MHB Site Helper
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
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?

#### evinda

##### Well-known member
MHB Site Helper
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
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".