# Is this the correct procedure for calculating matrix spectrum?

I am not sure if my question is on topic but I have a piece of Fortran code that is used to perform successive over relaxation. Prior to performing successive over relaxation the author is calculating the spectrum of the matrix.There is very little documentation so I am just trying to figure out all the details for documentation. It appears from the code comments that the following code calculates matrix spectrum. Would it be possible to get a mathematical background for this code ? Here is the equation for the successive over relaxation as requested.

The linearized second order partial differential equation(Laplacian) that I am trying to solve is this

$$\Delta_x(A.\Delta_x\psi) + \Delta_y(B.\Delta_y\psi) + \Delta_z(C.\Delta_z\psi) + S = 0$$

where $$\Delta_x(A.\Delta_x\psi) = A_{i+1/2,j,k} .(\psi_{i+1,j,k}-\psi_{i,j,k}) - A_{i-1/2,j,k} .(\psi_{i,j,k}-\psi_{i-1,j,k})$$

$$\Delta_y(B.\Delta_y\psi) = B_{i,j+1/2,k} .(\psi_{i,j+1,k}-\psi_{i,j,k}) - B_{i,j-1/2,k} .(\psi_{i,j,k}-\psi_{i,j-1,k})$$

$$\Delta_z(C.\Delta_z\psi) = C_{i,j,k+1/2} .(\psi_{i,j,k+1}-\psi_{i,j,k}) - A_{i,j,k-1/2} .(\psi_{i,j,k}-\psi_{i,j,k-1})$$

$$S_{i,j,k} = \Delta_x .\Delta_y . \Delta_z * \sigma_{i,j,k}$$

$$\psi_i^{(n+1)} = \omega\left( b_i - \sum_{j=1}^{i-1} A_{ij} \psi_j^{(n+1)} - \sum_{j=i}^{m}A_{ij}\psi_j^{(n)}\right) + (1-\omega)\psi_i^{(n)}.$$

  i=nx/2
specx=4.*a(i,0,0)/
>      (2.*a(i,0,0)+b(i,0,0)+b(i,1,0)+c(i,0,0)+c(i,0,1))
specy=2.*(b(i,0,0)+b(i,1,0))/
>      (2.*a(i,0,0)+b(i,0,0)+b(i,1,0)+c(i,0,0)+c(i,0,1))
specz=2.*(c(i,0,0)+c(i,0,1))/
>      (2.*a(i,0,0)+b(i,0,0)+b(i,1,0)+c(i,0,0)+c(i,0,1))
do k=1,nz-2
do j=1,ny-2

helpx=4.*a(i,j,k)/
>            (2.*a(i,j,k)+b(i,j,k)+b(i,j-1,k)+c(i,j,k)+c(i,j,k-1))
if (helpx.gt.specx) specx=helpx

helpy=2.*(b(i,j,k)+b(i,j+1,k))/
>            (2.*a(i,j,k)+b(i,j,k)+b(i,j-1,k)+c(i,j,k)+c(i,j,k-1))
if (helpy.gt.specy) specy=helpy

helpz=2.*(c(i,j,k)+c(i,j,k+1))/
>            (2.*a(i,j,k)+b(i,j,k)+b(i,j-1,k)+c(i,j,k)+c(i,j,k-1))
if (helpz.gt.specz) specz=helpz

enddo
enddo


The following code fragments calculate the a,b,c matrices

  do i=0,nx
do j=0,ny
do k=0,nz
a(i,j,k)=dy*dz*rhoref(2*k)/(dx*coriol(i,j))
if (j.lt.ny) then
b(i,j,k)=dx*dz*rhoref(2*k)/
>                     (dy*0.5*(coriol(i,j)+coriol(i,j+1)))
else
b(i,j,k)=0.
endif
if (k.lt.nz) then
c(i,j,k)=dx*dy*rhoref(2*k+1)*coriol(i,j)/
>                     (dz*nsq(2*k+1))
else
c(i,j,k)=0.
endif
enddo
enddo
enddo

• What are $a$, $b$, and $c$, and what linear system is it solving? – Kirill Jun 29 '15 at 5:24
• No, I meant what are the matrices precisely? Give details, expressions, anything! I don't know that equation. How do the matrices enter into it? What is the equation exactly, and how does the code try to solve it? – Kirill Jun 29 '15 at 5:41
• I converted your formula from JPEG to LaTeX (BTW, compare it with the similar one at en.wikipedia.org/wiki/Successive_over-relaxation), but what has it got to do with the code at hand? See my first comment. – Kirill Jun 29 '15 at 11:26
• @Kirill - Better late than never. I answered your comment. I hope :-) – gansub Sep 7 '16 at 15:18

The spectrum of a matrix is a set of numbers, so since there isn't a set of numbers being computed here, it's almost certainly not the spectrum. Also relevant is that you have three tensors ($a$, $b$, and $c$), so which spectrum do you think it would calculate? Calculating the spectrum of a matrix numerically also cannot be done in a finite number of rational operations and requires an iterative algorithm, unless something very special is known about the matrix.