I am trying to create a Navier Stokes solver with a bunch of subroutines and multiple iterations which also uses some LAPACK subroutines on MAC through the accelerate Framework.
The program compiles correctly on the first iteration but collapses somewhere in the middle saying this:
a.out(28015,0x7fff74e78300) malloc: *** error for object
0x7f81e902a008: incorrect checksum for freed object -
object was probably modified after being freed.
*** set a breakpoint in malloc_error_break to debug
I am looking for some advice on how to try and debug this. A simpler version of the program runs smoothly but when I try to add the following subroutine to modify a variable I get the error. I am reproducing the subroutine here too. The subroutine modifies ui and vi variables:
subroutine interpforce(ui,vi,ub,vb,
1 xb,yb,x,y,m,n,uz,h,dt)
implicit double precision (a-h,o-z)
integer uz,k,q,IPIV,INFO
dimension ub(0:uz-1),vb(0:uz-1),xb(0:uz),yb(0:uz)
1 ,x(-1:m+1),y(-1:n+1),a(0:uz-1,0:uz-1),
2 fkx(0:uz-1),fky(0:uz-1),dels(0:uz-1),
3 bx(0:uz-1),by(0:uz-1),IPIV(uz),
4 fex(0:m,0:n),fey(0:m,0:n),ui(0:m,0:n),vi(0:m,0:n)
do i=0,m
do j=0,m
fex(i,j)=0.0d0
fey(i,j)=0.0d0
enddo
enddo
do k=0,uz-1
fkx(k)=0.0d0
fky(k)=0.0d0
dels(k)=sqrt((xb(k+1)-xb(k))**2 + (yb(k+1)-yb(k))**2)
do q=0,uz-1
a(k,q)=0.0d0
enddo
enddo
do k=0,uz-1
do q=0,uz-1
do i=0,m
do j=0,n
a(k,q)=a(k,q)+(umdirac(x(i)-xb(q),h)*
1 umdirac(y(j)-yb(q),h)*dels(q)*h*
2 umdirac(x(i)-xb(k),h)*umdirac(y(j)-yb(k),h))
enddo
enddo
enddo
enddo
do k=0,uz-1
bx(k)=-ub(k)/dt
by(k)=-vb(k)/dt
enddo
call DGESV(uz,1,a,uz,IPIV,bx,uz,INFO)
call DGESV(uz,1,a,uz,IPIV,by,uz,INFO)
do k=0,uz-1
fkx(k)=bx(k)
fky(k)=by(k)
enddo
do i=0,m
do j=0,n
do k=0,uz-1
fex(i,j)=fex(i,j)+fkx(k)*(umdirac(x(i)-xb(k),h))
1 *(umdirac(y(j)-yb(k),h))*dels(k)
fey(i,j)=fey(i,j)+fky(k)*(umdirac(x(i)-xb(k),h))
1 *(umdirac(y(j)-yb(k),h)*dels(k))
enddo
enddo
enddo
do i=0,m
do j=0,n
ui(i,j)=ui(i,j)+dt*fex(i,j)
vi(i,j)=vi(i,j)+dt*fey(i,j)
enddo
enddo
return
end
c ************************************ UMDIRAC FUNCTION ***************************************
function umdirac(x,h)
implicit real*8 (a-h,o-z)
if (abs(x).lt.h) then
umdirac=(h-abs(x))/h**2
else
umdirac=0.0d0
endif
return
end
This kind of error is generally the result of memory corruption from a stray write operation.
Assuming you're using
gfortran
, try adding the compiler option-fcheck=bounds
and running again.If that doesn't do it, try running your program under
valgrind
to see if it can spot your error.