[gmx-users] Reading XTC files from fortran90
Jones de Andrade
johannesrs at gmail.com
Sat Sep 29 06:00:12 CEST 2007
Hi all.
I'm writing this beucase I'm having terrible problems in trying to deal with
gromacs trajectory files with my own programs, originally written on
fortran90.
The codes themselves are enoughly big to make it out of question to try to
"translate" them. What leads to some sort of mixed language compilation.
My first guess was to try to use the guidelines from
http://xray.bmc.uu.se/~spoel/md/online/xtc.html, but it absolutelly didn't
worked. Could not find a compatible library. Tried to recompile gromacs in
order to do so, but instead of a libsdrf.a library file, it created just
common .lo and .o object files.
Looked inside concoord program, as there are some mentions in the list to
look at the code. Was a really good aprentize, could understand better the
logics involved when dealing with .xtc files, BUT I could not use the
library files that came with that. Somehow, it wasn't recognized as a valid
file.
My last try up to now was to try to link my test-read-program with the
libxdrf object file. It yelded the following error:
CruNumMac src/own/B/9/xtc 238% ifort libxdrf.lo xtciof.f90 test_xtc.f90
IPO: WARNING: no IR in object file libxdrf.lo; was the source file compiled
with -ipo?
IPO Error: unresolved : xdrfint_
Referenced in /tmp/ipo_ifortZNvZHh.o
IPO Error: unresolved : xdrfclose_
Referenced in /tmp/ipo_ifortZNvZHh.o
IPO Error: unresolved : xdrfopen_
Referenced in /tmp/ipo_ifortZNvZHh.o
IPO Error: unresolved : xdrffloat_
Referenced in /tmp/ipo_ifortZNvZHh.o
IPO Error: unresolved : xdrf3dfcoord_
Referenced in /tmp/ipo_ifortZNvZHh.o
IPO: performing multi-file optimizations
IPO: generating object file /tmp/ipo_ifortZNvZHh.o
test_xtc.f90(17) : (col. 3) remark: LOOP WAS VECTORIZED.
test_xtc.f90(25) : (col. 8) remark: LOOP WAS VECTORIZED.
xtciof.f90(112) : (col. 3) remark: LOOP WAS VECTORIZED.
libxdrf.lo: file not recognized: File format not recognized
If I try to link to the .o objetc:
CruNumMac src/own/B/9/xtc 240% ifort libxdrf.o xtciof.f90 test_xtc.f90
IPO: WARNING: no IR in object file libxdrf.o; was the source file compiled
with -ipo?
IPO Error: unresolved : xdrfint_
Referenced in /tmp/ipo_ifort2f5KMw.o
IPO Error: unresolved : xdrfclose_
Referenced in /tmp/ipo_ifort2f5KMw.o
IPO Error: unresolved : xdrfopen_
Referenced in /tmp/ipo_ifort2f5KMw.o
IPO Error: unresolved : xdrffloat_
Referenced in /tmp/ipo_ifort2f5KMw.o
IPO Error: unresolved : xdrf3dfcoord_
Referenced in /tmp/ipo_ifort2f5KMw.o
IPO Error: unresolved : ffclose
Referenced in libxdrf.o
IPO: performing multi-file optimizations
IPO: generating object file /tmp/ipo_ifort2f5KMw.o
test_xtc.f90(17) : (col. 3) remark: LOOP WAS VECTORIZED.
test_xtc.f90(25) : (col. 8) remark: LOOP WAS VECTORIZED.
xtciof.f90(112) : (col. 3) remark: LOOP WAS VECTORIZED.
/tmp/ipo_ifort2f5KMw.o(.text+0x1e1): In function `MAIN__':
/tmp/ipo_ifort2f5KMw.f: undefined reference to `xdrfopen_'
/tmp/ipo_ifort2f5KMw.o(.text+0x4ea):/tmp/ipo_ifort2f5KMw.f: undefined
reference to `xdrfclose_'
/tmp/ipo_ifort2f5KMw.o(.text+0x517):/tmp/ipo_ifort2f5KMw.f: undefined
reference to `xdrfopen_'
/tmp/ipo_ifort2f5KMw.o(.text+0x92b): In function `xtciof_mp_xtcheader_':
/tmp/ipo_ifort2f5KMw.f: undefined reference to `xdrfint_'
/tmp/ipo_ifort2f5KMw.o(.text+0x941):/tmp/ipo_ifort2f5KMw.f: undefined
reference to `xdrfint_'
/tmp/ipo_ifort2f5KMw.o(.text+0x95e):/tmp/ipo_ifort2f5KMw.f: undefined
reference to `xdrfint_'
/tmp/ipo_ifort2f5KMw.o(.text+0x97b):/tmp/ipo_ifort2f5KMw.f: undefined
reference to `xdrffloat_'
/tmp/ipo_ifort2f5KMw.o(.text+0xb8e): In function `xtciof_mp_xtcio_':
/tmp/ipo_ifort2f5KMw.f: undefined reference to `xdrffloat_'
/tmp/ipo_ifort2f5KMw.o(.text+0xbde):/tmp/ipo_ifort2f5KMw.f: undefined
reference to `xdrffloat_'
/tmp/ipo_ifort2f5KMw.o(.text+0xc2e):/tmp/ipo_ifort2f5KMw.f: undefined
reference to `xdrffloat_'
/tmp/ipo_ifort2f5KMw.o(.text+0xc7e):/tmp/ipo_ifort2f5KMw.f: undefined
reference to `xdrffloat_'
/tmp/ipo_ifort2f5KMw.o(.text+0xcce):/tmp/ipo_ifort2f5KMw.f: more undefined
references to `xdrffloat_' follow
/tmp/ipo_ifort2f5KMw.o(.text+0xe6b): In function `xtciof_mp_xtcio_':
/tmp/ipo_ifort2f5KMw.f: undefined reference to `xdrf3dfcoord_'
/tmp/ipo_ifort2f5KMw.o(.text+0xfc6): In function `xtciof_mp_xtcopen_':
/tmp/ipo_ifort2f5KMw.f: undefined reference to `xdrfopen_'
/tmp/ipo_ifort2f5KMw.o(.text+0x1292):/tmp/ipo_ifort2f5KMw.f: undefined
reference to `xdrfclose_'
/tmp/ipo_ifort2f5KMw.o(.text+0x12b3):/tmp/ipo_ifort2f5KMw.f: undefined
reference to `xdrfopen_'
/tmp/ipo_ifort2f5KMw.o(.text+0x134e): In function `xtciof_mp_xtccoord_':
/tmp/ipo_ifort2f5KMw.f: undefined reference to `xdrffloat_'
/tmp/ipo_ifort2f5KMw.o(.text+0x139f):/tmp/ipo_ifort2f5KMw.f: undefined
reference to `xdrffloat_'
/tmp/ipo_ifort2f5KMw.o(.text+0x13f0):/tmp/ipo_ifort2f5KMw.f: undefined
reference to `xdrffloat_'
/tmp/ipo_ifort2f5KMw.o(.text+0x1441):/tmp/ipo_ifort2f5KMw.f: undefined
reference to `xdrffloat_'
/tmp/ipo_ifort2f5KMw.o(.text+0x1492):/tmp/ipo_ifort2f5KMw.f: undefined
reference to `xdrffloat_'
/tmp/ipo_ifort2f5KMw.o(.text+0x14e3):/tmp/ipo_ifort2f5KMw.f: more undefined
references to `xdrffloat_' follow
/tmp/ipo_ifort2f5KMw.o(.text+0x162d): In function `xtciof_mp_xtccoord_':
/tmp/ipo_ifort2f5KMw.f: undefined reference to `xdrf3dfcoord_'
libxdrf.o(.text+0x3d86): In function `xdrclose':
libxdrf.c: undefined reference to `ffclose'
libxdrf.o(.text+0x3da8):libxdrf.c: undefined reference to `ffclose'
I'm really not getting what can be the source of error. I placed the two
test wources from above trials below, in case it can be of any help. But I
really not getting what am I doing wrong, and also could not find any other
clue in gromacs mail lists, neither google.
Sorry it this is a much simpler question than it seems to me, but I'm really
astonished not finding the source of error neither a workaround or an
altertive way which I could follow and reach some results.
Thanks a lot in advance,
Sincerally yours,
Jones
********Attachments:
*****test_xtc.f90:
program test_xtc
use xtciof
implicit none
! Declaracao de Variaveis:
integer :: xd, natoms, step, ret
integer :: nat, maxat
! real*4 :: time, box(9), x(*), prec
real :: time, box(9), prec
real, dimension(:), allocatable :: x
character*1 :: stat
character(len=1) :: file_name
maxat = 10000
stat = 'a'
file_name = 'teste.xtc'
allocate(x(30000))
x(:) = 0.
! Aqui comeca o programa de teste:
call xtcopen(xd, file_name, stat, maxat, nat)
write(6,*) 'xtc file openned'
call readxtc(xd, natoms, step, time, box, x, prec, ret)
write(6,*) 'xtc file readed'
call writextc(xd, natoms, step, time, box, x, prec, ret)
write(6,*) 'xtc file written'
deallocate(x)
end program test_xtc
******xtciof.f90:
module xtciof
contains
subroutine xtccheck(ret,ivar)
implicit none
!
! Passed variables
!
integer, intent(in) :: ret, ivar
if ( ret .eq. 0 ) then
if (ivar == 1) write(6,*) '> XTC-Error reading/writing natoms'
if (ivar == 2) write(6,*) '> XTC-Error reading/writing step'
if (ivar == 3) write(6,*) '> XTC-Error reading/writing time'
if (ivar == 4) write(6,*) '> XTC-Error reading/writing box'
if (ivar == 5) write(6,*) '> XTC-Error reading/writing x'
! stop
endif
end subroutine xtccheck
subroutine xtcheader(xd,magic,natoms,step,time,ret)
implicit none
!
! Passed variables
!
integer :: xd, natoms, step, ret, magic
! real*4 :: time
real :: time
call xdrfint(xd, magic, ret)
if (ret == 0) return
call xdrfint(xd, natoms, ret)
call xtccheck(ret, 1)
call xdrfint(xd, step, ret)
call xtccheck(ret, 2)
call xdrffloat(xd, time, ret)
call xtccheck(ret, 3)
end subroutine xtcheader
subroutine xtccoord(xd,natoms,box,x,prec,ret)
implicit none
!
! Passed variables
!
integer :: xd, natoms, ret
! real*4 :: box(9), x(*), prec
real :: box(9), x(*), prec
!
! local
!
integer :: i
do i=1, 9
call xdrffloat(xd, box(i), ret)
call xtccheck(ret, 4)
enddo !i
call xdrf3dfcoord(xd, x, natoms, prec, ret)
call xtccheck(ret, 5)
end subroutine xtccoord
subroutine xtcio(xd, natoms, step, time, box, x, prec, mode, ret)
implicit none
!
! Passed variables
!
integer :: xd, natoms, step, mode, ret
! real*4 :: box(9)
! real*4 :: box(9)
! real*4 :: x(*)
real :: time, prec
real :: x(*)
real :: box(9)
!
! local variables
!
integer :: xtcmagic
integer :: magic
xtcmagic=1995
if (mode .eq. 0) magic = xtcmagic
call xtcheader(xd, magic, natoms, step, time, ret)
if (ret == 0) return
if (mode == 1) then
if (magic == xtcmagic) write(6,*) '> Fatal error, magic number read as
',magic, ' should be ',xtcmagic
endif !(mode == 1)
call xtccoord(xd, natoms, box, x, prec, ret)
end subroutine xtcio
subroutine readxtc(xd, natoms, step, time, box, x, prec, ret)
implicit none
integer :: xd, natoms, step, ret
integer :: nat
! real*4 :: time, box(9), x(*), prec
real :: time, box(9), prec
real, dimension(:), allocatable :: x
allocate (x(natoms))
x(:) = 0.
call xtcio(xd,natoms,step,time,box,x,prec,1,ret)
end subroutine readxtc
subroutine writextc(xd, natoms, step, time, box, x, prec, ret)
implicit none
integer :: xd, natoms, step, ret
real :: time, box(9), x(*), prec
call xtcio(xd, natoms, step, time, box, x, prec, 0, ret)
end subroutine writextc
subroutine xtcopen (xd, file_name, stat, maxat, nat)
implicit none
integer :: xd, ret, nat, maxat, magic, step
! real*4 :: time
real :: time
character*1 :: stat
character*(*) :: file_name
integer, parameter :: xtcmagic = 1995
call xdrfopen(xd, file_name, stat, ret)
if (ret == 0) write(6,*) '> Fatal error, could not open ',file_name(1:20)
call xtcheader(xd, magic, nat, step, time, ret)
if (magic /= xtcmagic) then
write(6,*) '> Fatal error, magic number read as ',magic,' should be
',xtcmagic
stop
endif
if (nat < 0) then
write(6,*) '> Fatal error, negative number of atoms'
stop
endif
if (nat > maxat) then
write(6,*) '> Fatal error, number of atoms (',nat,') ','larger than
allowed (',maxat,')'
stop
endif
call xdrfclose(xd, ret)
call xdrfopen(xd, file_name, stat, ret)
end subroutine xtcopen
!end contains
end module xtciof
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://maillist.sys.kth.se/pipermail/gromacs.org_gmx-users/attachments/20070929/83ecf9c6/attachment.html>
More information about the gromacs.org_gmx-users
mailing list