Exceeding of array bounds

[ferry2]

Hi, fellows.

I'm writing a scientific program for my thesis. But I have a little problem. I'm getting error message: severe <161>: Program Exception - array bounds exceeded for this line:

r(1) = alpha

Here is the whole subprogram:

subroutine Tridiagonal_Solv(n, a, b, alpha, beta, x, y, y_deriv, error)
	implicit none 
	
	integer n, i 
 	double precision h 
	double precision s(n), t(n) 
	double precision down(n), main(n), up(n) 
	double precision exact_sol(n + 2) 
	double precision r(n) 
	double precision x(n + 2) 
	double precision y(n), y_deriv(n) 
	double precision error

	common /DIM/ D(1)
	integer D
	common /PAR/P(4)
	double precision P

	n = D(1)
	a = P(1)
	b = P(2)
	alpha = P(3)
	beta = P(4)


	h = Step(n, a, b) 
	x(1) = a
	r(1) = alpha
	exact_sol(1) = ExactSolution(x(1))
	error = r(1)**2
	do i = 2, n - 1
		x(i) = a + (i - 1)*h
		r(i) = 2.0d0*dsin(x(i))
		exact_sol(i) = ExactSolution(x(i))
		error = error + r(i)**2
		!print *, "r(", i, ") = ", r(i)
	end do
	r(n) = beta
	error = error + r(n)**2
	x(n) = b
	exact_sol(n) = ExactSolution(x(n))
	x(n + 1) = a + n*h
	exact_sol(n + 1) = ExactSolution(x(n + 1))
	x(n + 2) = a + (n + 1)*h
	exact_sol(n + 2) = ExactSolution(x(n + 2))

	do i = 1, n
		y_deriv(i) = (-3.0d0*exact_sol(i) + 4.0d0*exact_sol(i + 1) - exact_sol(i + 2))/2.0d0*h
	end do

	do i = 2, n - 1
		down(i) = -1.0d0/(h*h)
	end do
	down(n) = 2.0d0/h

	main(1) = -3.0d0/(2.0d0*h)
	do i = 2, n - 1
		main(i) = (2.0d0/(h*h)) + 1.0d0
	end do
	main(n) = 3.0d0/(2.0d0*h)

	up(1) = 2.0d0/h
	do i = 2, n - 1
		up(i) = -1.0d0/(h*h)
	end do
	
	s(1) = -(main(1)*up(2) - main(2)*up(1))/(down(1)*up(2) - down(2)*up(1))
	t(1) = (up(2)*r(1) - up(1)*r(2))/(down(1)*up(2) - down(2)*up(1))
	s(2) = -(up(1)*down(2) - up(2)*down(1))/(main(1)*down(2) - main(2)*down(1))
	t(2) = (down(2)*r(1) - down(1)*r(2))/(main(1)*down(2) - main(2)*down(1))

	do i = 3, n - 2
		s(i) = -up(i)/(main(i) + down(i)*s(i - 1))
		t(i) = (r(i) - down(i)*t(i - 1))/(main(i) + down(i)*s(i - 1))
	end do

	! Çàïèñâàìå ðåçóëòàòèòå îò ðåøåíèåòî íà ñèñòåìàòà âúâ âåêòîðà y
	
	y(n) = ((r(n - 1) - t(n - 2)*down(n - 1))*(s(n - 2)*down(n) + main(n)) - (r(n) - t(n - 2)*down(n))*(s(n - 2)*down(n - 1) + main(n - 1)))/(up(n - 1)*(down(n)*s(n - 2) + main(n)) - up(n)*(down(n - 1)*s(n - 2) + main(n - 1)))
	y(n - 1) = (r(n)*up(n - 1) - t(n - 2)*down(n)*up(n - 1) - r(n - 1)*up(n) + t(n - 2)*down(n - 1)*up(n))/(s(n - 2)*down(n)*up(n - 1) + main(n)*up(n - 1) - s(n - 2)*down(n - 1)*up(n) - main(n - 1)*up(n))

	do i = n - 2, 3, -1
		y(i) = s(i)*y(i + 1) + t(i)
	end do

	y(2) = s(2)*y(3) + t(2)
	y(1) = s(1)*y(2) + t(1)
	
	!print *
	!do i = 1, n
	!	print *, "y(", i, ") = ", y(i) 
	!end do

end subroutine Tridiagonal_Solv

I already check that the variable n is properly readed form the data file.

 

[gummibaer]

What is the actual value of n at the time you subroutine is called ?
If it is zero, then r would be declared as an array of zero elements.

In fact, you set n = D(1) after your arrays are declared in this sub and then, later on you reference r and other array elements. This is bound to create trouble if D(1) is bigger than n at the time your sub is called.

BTW: seems you pass n to your sub twice, one as n in the argument list and then the other as D(1) in the common. If both are equal, you do not need to do this, if they are not, you are going to have trouble.


Norbert


The optimist believes we live in the best of all possible worlds - the pessimist fears this might be true.

 

[salgerman]

Yeah, ferry, you may be stepping on your toes by overwriting 'n'...which is it?

By the way, I use g95 for compiling my fortran programs and
I typically use the "-fbounds-check" flag so that if something
happens at runtime with array indeces, etc...it lets me know
with more detail which array, which dimension and index value.

 

[gummibaer]

Just looking a little more deeply into your code, especially your argument's list and data declaration. I guess you have to pay a lot of attention to keep your data declarations and dimensions coherent.

You could simplify this by using modules. this would look something like this:

module mydata
    integer, parameter :: n = 100   ! .. or whatever you need as dimension
    double precision a                       ! I would add an explanation what a is here
    double precision b
    double precision r(n)
    .
    .
    .
end module

You reference this module in all your program-units that need these data as

subroutine Tridiagonal_Solv (...)    ! only those data in the argument list that are not in the module
use mydata
.
.
.

and you do not have to worry about dimensions again.

Norbert


The optimist believes we live in the best of all possible worlds - the pessimist fears this might be true.

 

[melmacianalf]

BTW, why do you have too lengthy lines?

y(n) = ((r(n - 1) - t(n - 2)*down(n - 1))*(s(n - 2)*down(n) + main(n)) - (r(n) - t(n - 2)*down(n))*(s(n - 2)*down(n - 1) + main(n - 1)))/(up(n - 1)*(down(n)*s(n - 2) + main(n)) - up(n)*(down(n - 1)*s(n - 2) + main(n - 1)))
y(n - 1) = (r(n)*up(n - 1) - t(n - 2)*down(n)*up(n - 1) - r(n - 1)*up(n) + t(n - 2)*down(n - 1)*up(n))/(s(n - 2)*down(n)*up(n - 1) + main(n)*up(n - 1) - s(n - 2)*down(n - 1)*up(n) - main(n - 1)*up(n))

 

[ferry2]

@gummibaer the actual value of n is 5. @melmacianalf because the formulas are too long.

 

[gummibaer]

Ferry,

the vital point of my question is '...at the time your subroutine is called'. If n really is 5 at the time your type declaration takes place, this error would never occur. Change your prog either to

program main
.
.
.
write (*,*) 'n = ', n
call Tridiagonal_Solv(n, a, b, alpha, beta, x, y, y_deriv, error)
.
.

or

subroutine Tridiagonal_Solv(n, a, b, alpha, beta, x, y, y_deriv, error)
	implicit none 
	
	integer n, i 
 	double precision h 
	double precision s(n), t(n) 
	double precision down(n), main(n), up(n) 
	double precision exact_sol(n + 2) 
	double precision r(n) 
	double precision x(n + 2) 
	double precision y(n), y_deriv(n) 
	double precision error

	common /DIM/ D(1)
	integer D
	common /PAR/P(4)
	double precision P

         write (*,*) 'n = ', n   ! important: before you set n to D(1)

	n = D(1)
	a = P(1)
	b = P(2)

I guess you should reconsider how you pass your data to your sub. Having a argument list and then setting the values of the arguments by what is passed via the common-statements is not just a little strange.

Norbert




The optimist believes we live in the best of all possible worlds - the pessimist fears this might be true.