Laplace expansion in fortran

1. Jun 25, 2015

kinbeat

Hi, i'm dealing with an exam problem (which i have 8 days to solve, with no outside help limitations), which is to write a fortran95 program that solves linear systems.
the first part asks to find the determinant of a NxN matrix with the laplace expansion, implementing it as a recursive function.
this is what i've done so far

module prec
integer, parameter :: rk= selected_real_kind(6)
end module prec

module lap
use prec
implicit none

contains
recursive function det(a) result (d)
real(kind=rk), intent(in), dimension(:,:) :: a
real(kind=rk), allocatable,dimension(:,:) :: b
real(kind=rk) :: d
integer ::i,n
n=0
b(1:i-1,1:n) = a(1:i-1,2:n)
b(i:n,1:n) =a(i+1:n,2:n)

if (size(a) > 2) then
do i=1,n
d= ((-1)**(i+1))*a(i,1)*det(b)
end do
else
d = a(1,1)*a(2,2)-a(1,2)*a(2,1)

end if
end function det
end module lap

program sistema
use lap
real(kind=rk), allocatable, dimension(:,:) :: c
real(kind=rk) :: z
integer ::j,m
print*,'Inserisci la dimensione della matrice'
read*,m
allocate(c(m,m))
print*, 'inserisci i coefficienti della matrice'
read*, c

do j=1,m
print*,c(j,:)
end do

z=det(c)

print*, z

end program sistema

but the compiler gives me this error:
error 112, reference to undefined variable, array element of function result (/undef)
LAP!DET - in file matrice.f95 at line 18 [+00aa]
main - in file matrice.f95 at line 51 [+0521]

what am i not seeing?

2. Jun 25, 2015

Staff: Mentor

Hi kinbeat.

Code (Text):
integer ::i,n
n=0
b(1:i-1,1:n) = a(1:i-1,2:n)
b(i:n,1:n) =a(i+1:n,2:n)
Does i have a value at the stage of being used on the 3rd line here?

3. Jun 25, 2015

kinbeat

ehi there!
yeah, i got a bit further now, i realized i didn't declare the n.
now i get the "nonconformant arrays" error, at that line, "
b(1:i-1,1:n) = a(1:i-1,2:n)

what i'm trying to tell the compiler is that the matrix b is a minus the i row and first column.

now i changed that part so it looks like this:

recursive function det(a,n) result (d)
real(kind=rk), intent(in), dimension(:,:) :: a

real(kind=rk) :: d
integer ::i
integer, intent(in)::n

if (size(a) > 2) then
do i=1,n

b(1:i-1,1:n) = a(1:i-1,2:n)
b(i:n,1:n) =a(i+1:n,2:n)
d= ((-1)**(i+1))*a(i,1)*det(b,n)
end do
else
d = a(1,1)*a(2,2)-a(1,2)*a(2,1)

end if

Last edited: Jun 25, 2015
Know someone interested in this topic? Share this thread via Reddit, Google+, Twitter, or Facebook

Have something to add?
Draft saved Draft deleted