Hi everyone,

Below is a subroutine I've created (the variable fit_coeffs is defined externally and is allocatable, it's not the problem). When I compile I get memory allocation errors, that appear when I assign fit_coeffs, due to the matmul(ATA,AT) line. I know this from inserting a bunch of print statements. Also, both error checking statements are invoked from the dgetrf and dgetri routines. WHy might this be happening? I'm compiling using the command:

gfortran -Wall -cpp -std=f2003 -ffree-form -L/home/***USER***/lapack-3.4.0/ file.f -llapack -lrefblas.

Thanks in advance!

subroutine polynomial_fit(x_array, y_array, D)

integer, intent(in) :: D

real, intent(in), dimension(:) :: x_array, y_array

real, allocatable, dimension(:,:) :: A, AT, ATA

real, allocatable, dimension(:) :: work

integer, dimension(:), allocatable :: pivot

integer :: l, m, n, lda, lwork, ok

l = D + 1

lda = l

lwork = l

allocate(fit_coeffs(l))

allocate(pivot(l))

allocate(work(l))

allocate(A(size(x_array),l))

allocate(AT(l,size(x_array)))

allocate(ATA(l,l))

do m = 1,size(x_array),1

do n = 1,l,1

A(m,n) = x_array(m)**(n-1)

end do

end do

AT = transpose(A)

ATA = matmul(AT,A)

call dgetrf(l, l, ATA, lda, pivot, ok)

! ATA is now represented as PLU (permutation, lower, upper)

if (ok /= 0) then

write(6,*) "HERE"

end if

call dgetri(l, ATA, lda, pivot, work, lwork, ok)

! ATA now contains the inverse of the matrix ATA

if (ok /= 0) then

write(6,*) "HERE"

end if

fit_coeffs = matmul(matmul(ATA,AT),y_array)

deallocate(pivot)

deallocate(fit_coeffs)

deallocate(work)

deallocate(A)

deallocate(AT)

deallocate(ATA)

end subroutine polynomial_fit