【发布时间】：2014-09-29 16:09:44
【问题描述】：
我一直在 FORTRAN 中编写代码，但在使用 lapack dsyevr 时遇到问题：
http://netlib.sandia.gov/lapack/double/dsyevr.f
我遇到的问题似乎与内存分配问题有关，特别是我认为与 dsyevr 生成的输出数组有关（包括作为输入和输出的 A）。
我尝试编写一个简化的代码来演示我遇到的问题。请让我知道是否需要澄清。代码名称为prof1.f90，调用dsyevr函数：
PROGRAM prog1
implicit none
real(kind=8), allocatable :: W(:)
real(kind=8), allocatable :: Z(:,:)
real(kind=8), allocatable :: A(:,:)
integer(kind=8) :: n, info, il, iu, m, lwork, liwork
integer(kind=8) :: i, k, p, q, nu
real(kind=8) :: abstol, vl, vu
real(kind=8), allocatable :: work(:)
integer, allocatable :: isuppz(:), iwork(:)
n = 3
allocate(W(3),Z(3,3),A(n,n),stat=info)
if (info .ne. 0) stop "error allocating arrays"
A(1,1)=3.78136524999999994E-003
A(1,2)=0.0000000000000000
A(1,3)=-7.92918150000000038E-004
A(2,1)=0.0000000000000000
A(2,2)=5.20293929999999984E-003
A(2,3)=0.0000000000000000
A(3,1)=-7.92918150000000038E-004
A(3,2)=0.0000000000000000
A(3,3)=3.78136524999999994E-003
vl = 1.06451084056294826E-313
vu = 0.0
il = 4294967297
iu = 8839891
m = 140733655445712
W(1) = 2.98844710000000001E-003
W(2) = 4.57428340000000030E-003
W(3) = 5.20293929999999984E-003
Z(1,1) = 8.65587596665713699E-317
Z(1,2) = 8.65587596665713699E-317
Z(1,3) = 1.58101006669198894E-322
Z(2,1) = 1.58101006669198894E-322
Z(2,2) = 0.0000000000000000
Z(2,3) = 8.65569415049946741E-317
Z(3,1) = 4.24400777097956191E-314
Z(3,2) = 4.79243676466009148E-322
Z(3,3) = 3.51391740150311405E-316
lwork = -1
liwork = -1
abstol = 1d-5
allocate(work(1),iwork(1),isuppz(6))
call dsyevr('V','A','U',n,A,n,vl,vu,il,iu,abstol,m,W,Z,n,isuppz,work,lwork,iwork,liwork,info)
if (info .ne. 0) stop "error obtaining work array dimensions"
lwork = work(1)
liwork = iwork(1)
deallocate(work,iwork)
allocate(work(lwork),iwork(liwork),stat=info)
if (info .ne. 0) stop "error allocating work arrays"
call dsyevr('V','A','U',n,A,n,vl,vu,il,iu,abstol,m,W,Z,n,isuppz,work,lwork,iwork,liwork,info)
if (info .ne. 0) stop "error diagonalizing the hamiltonian"
deallocate(A,work,iwork,isuppz)
END PROGRAM prog1
在上面的代码中，dsyevr 函数被调用了两次，第一次只是为了获取工作的尺寸等...矩阵正确运行但是第二次调用它时返回以下错误
*** glibc detected *** ./PROGRAM: munmap_chunk(): invalid pointer: 0x000000000134cc20 ***
我还可以提供一个 Backtrace 和 MemoryMap。
如果有用的话，我一直在使用的 makefile 如下所示。该程序是使用以下行创建的：
make PROGRAM
生成文件：
FC = gfortran
FCFLAGS = -g -fbounds-check
FCFLAGS = -O2
FCFLAGS += -I/usr/include
%: %.o
$(FC) $(FCFLAGS) -o $@ $^ $(LDFLAGS)
%.o: %.f90
$(FC) $(FCFLAGS) -c $< -fno-range-check
%.o: %.F90
$(FC) $(FCFLAGS) -c $<
.PHONY: clean veryclean
PROGRAM: prog1.f90 prog1.o
$(FC) $(FCFLAGS) -o $@ prog1.o $(LIBS) -Wl,--start-group -L$(MKLROOT)/lib/intel64 -lmkl_gf_ilp64 -lmkl_core -lmkl_sequential -Wl,--end-group -lpthread
clean:
rm -f *.o *.mod *.MOD
veryclean: clean
rm -f *~ $(PROGRAM)
其中 $MKLROOT 是 /opt/intel/composer_xe_2011_sp1.8.273/mkl
如果有用的话我用过valgrind：
valgrind --tool=memcheck --db-attach=yes ./PROGRAM
我发现了以下错误：
==31069==
==31069== Invalid write of size 8
==31069== at 0x57B9F92: mkl_lapack_dsyevr (in /opt/intel/composer_xe_2011_sp1.8.273 /mkl/lib/intel64/libmkl_core.so)
==31069== by 0x4D9A580: DSYEVR (in /opt/intel/composer_xe_2011_sp1.8.273/mkl/lib /intel64/libmkl_gf_ilp64.so)
==31069== by 0x401163: MAIN__ (in /home/j/workbook/Test4/PROGRAM)
==31069== by 0x401F09: main (in /home/j/workbook/Test4/PROGRAM)
==31069== Address 0x6afe0b0 is 0 bytes inside a block of size 4 alloc'd
==31069== at 0x4A069EE: malloc (vg_replace_malloc.c:270)
==31069== by 0x400FE4: MAIN__ (in /home/j/workbook/Test4/PROGRAM)
==31069== by 0x401F09: main (in /home/j/workbook/Test4/PROGRAM)
==31069==
我不确定“Invalid write of size 8”是指某些整数或实数的大小（如 Jonathan Dursi 所述），还是指传递给 dsyevr 的数组的大小。
【问题讨论】：
-
Invalid write of size 8表示写入 8 个字节的单个元素。它可以是整数或实数。 -
我尝试将所有参数设置为 kind=4，但是当我这样做时，我收到一条错误消息：“MKL 错误：进入 DSYEVR 时参数 20 不正确”
-
难怪你用lp64 MKL，整数必须是8字节。
-
所以 MKL 要求所有 (?) 整数为 8 字节。但是 gfortran 编译器要求整数很长，即 4 字节（如 Jonathan Dursi 所述）。这是一个正确的。有没有办法解决这个问题？
-
当然，使用 8 字节类型的整数。