I am taking my first steps with MM_PBSA and I think I have encountered some
small bugs in the script and make_crd_hg.f files. Apologies in advance if
these have already been reported or are my errors!
At our site, the executable to build trajectories is called make_crd_hg -
not make_crd. Until I changed
.GC_FIL = ("MAKE_CRD");
to
.GC_FIL = ("MAKE_CRD_HG");
in mm_pbsa.pl, I was unable to generate trajectory snapshots due to the
message:
"File MAKE_CRD is not found"
I made this substitution three other places in the .pl script. Please let
me know if this is a mistake on my part, or if we have built the mm_pbsa
suite incorrectly here.
================================================================
I believe there may also be a bug in make_crd_hg.f
I'll use the 01_Example to demonstrate this - and avoid getting into
details of my own work.
If I change mm_pbsa.in from
NUMBER_REC_GROUPS 2
RSTART 1
RSTOP 2621
RSTART 3863
RSTOP 3907
to
NUMBER_REC_GROUPS 1
RSTART 1
RSTOP 2621
I get lots of happy messages followed by:
=>> Creating input
make_crd input
=>> Creating coordinates
Executing makecrd
Something wrong with NGR 2622 2621
Broken pipe
I think my mm_pbsa.in change should be "legal" If not, please let me know
why it is not. Otherwise, read on.
I have checked the make_crd_hg.f source code and I believe the error may be
in these statements:
j=1
k=1
NGR=0
NGL=0
NGC=0
do i=1,NGREC+NGLIG
if(((NGREC.gt.0).and.(RSTART(j).lt.LSTART(k))).or.
1 (k.gt.NGLIG)) then
do l=RSTART(j),RSTOP(j)
if(x(l).lt.big) then
NGR=NGR+1
xR(NGR)=x(l)
yR(NGR)=y(l)
zR(NGR)=z(l)
if(NGCOM.gt.0) then
NGC=NGC+1
xC(NGC)=x(l)
yC(NGC)=y(l)
zC(NGC)=z(l)
endif
endif
end do
j=j+1
elseif(((NGLIG.gt.0).and.(LSTART(k).lt.RSTART(j))).or.
1 (j.gt.NGREC)) then
do l=LSTART(k),LSTOP(k)
if(x(l).lt.big) then
NGL=NGL+1
xL(NGL)=x(l)
yL(NGL)=y(l)
zL(NGL)=z(l)
if(NGCOM.gt.0) then
NGC=NGC+1
xC(NGC)=x(l)
yC(NGC)=y(l)
zC(NGC)=z(l)
endif
endif
end do
k=k+1
endif
end do
C
C Check again
C
========================================================
I don't understand all the nuances of FORTRAN :) - but I think structurally
the error could be avoided with two separate
loops over NGR and NGL so that the dependence of the variables on both
NGREC and NGLIG (number of receptor and ligand groups)
could be separated out. Also, I think we can get rid of j and k:
Something like this might do the trick - though I have not compiled/tested
this:
NGR=0
NGL=0
NGC=0
if (NGREC.gt.0)
do i=1,NGREC
do l=RSTART(i),RSTOP(i)
if(x(l).lt.big) then <<<<< Let me add - if x(l).nlt.big
don't we have an error? - and shouldn't the program exit with a message?
NGR=NGR+1
xR(NGR)=x(l)
yR(NGR)=y(l)
zR(NGR)=z(l)
if(NGCOM.gt.0) then
NGC=NGC+1
xC(NGC)=x(l)
yC(NGC)=y(l)
zC(NGC)=z(l)
endif
end do
end do
if (NGLIG.gt.0)
do i=1,NGLIG
do l=LSTART(i),LSTOP(i)
if(x(l).lt.big) then <<<<< Let me add - if x(l).nlt.big
don't we have an error? - and shouldn't the program exit with a message?
NGL=NGL+1
xL(NGL)=x(l)
yL(NGL)=y(l)
zL(NGL)=z(l)
if(NGCOM.gt.0) then
NGC=NGC+1
xC(NGC)=x(l)
yC(NGC)=y(l)
zC(NGC)=z(l)
endif
endif
end do
end do
endif
C
C Check again
C
Thanks
Chris Moth
-----------------------------------------------------------------------
The AMBER Mail Reflector
To post, send mail to amber.scripps.edu
To unsubscribe, send "unsubscribe amber" to majordomo.scripps.edu
Received on Tue Sep 23 2003 - 18:53:00 PDT
