> That is good to know. Hopefully none of my suggestions (which are born
> mostly out of necessary habits from my CHARMM days) would make things
> worse. ^^^^^^
|
I rest my case. ---------------------------- ;-)
> Last I looked at the code (or maybe this is just our local development
> version), I think the coordinates are actually shifted by the QM region
> center of geometry (not that I can imagine that making much of a
> difference
> except for maybe really big QM regions). In any event, I am still
Here's the relevant code from qm_mm.f
subroutine qm_fill_qm_xcrd_periodic(x,natom, &
iqmatoms,scaled_mm_chrgs,real_scratch)
...
...
...
!Move the first QM atom to be at the origin.
! iqm_one = iqmatoms(1)
! frac(1:3) = x(1,iqm_one)*recip(1,1:3)+x(2,iqm_one)*recip(2,1:3)+ &
! x(3,iqm_one)*recip(3,1:3)
!Moving the first QM atom can be an issue if the box size is small
!since it could end up with atoms that can stick outside the box.
!Thus we can change this to calculate the center of the coordinates of
!the QM region and move this instead. This could cause problems if the QM
!region diffuses significantly but if that happens we are in a whole
!world of problems anyway and it should trigger the cutoff too large error.
... NOTE WE JUST USE THE FIRST QM atom position and move that to the origin
and assume given that none of the QM atoms stick out of the box. This should
be true unless the box is crazy small.
The QM atoms are then shiften by the 3D translation calculated for moving
the first atom to the origin.
...
...
do j=1,qmmm_struct%nquant
xbnd0=min(xbnd0,real_scratch(1,j))
xbnd1=max(xbnd1,real_scratch(1,j))
ybnd0=min(ybnd0,real_scratch(2,j))
ybnd1=max(ybnd1,real_scratch(2,j))
zbnd0=min(zbnd0,real_scratch(3,j))
zbnd1=max(zbnd1,real_scratch(3,j))
enddo
xbnd0=xbnd0-qmmm_nml%qmcut
ybnd0=ybnd0-qmmm_nml%qmcut
zbnd0=zbnd0-qmmm_nml%qmcut
xbnd1=xbnd1+qmmm_nml%qmcut
ybnd1=ybnd1+qmmm_nml%qmcut
zbnd1=zbnd1+qmmm_nml%qmcut
...
...
The code then finds x,y,zbnd0,1 and which give the maximum and minimum
coordinates of the QM region relative to having the first QM atom at the
origin. The qmcut is subtracted from the minimums and added to the maximums
to give you the extent of the QM+MM region in the X, Y and Z directions.
...
!------ Check if QM region plus cutoff around it is too large for this box
!--- The sphere method is the simplest check, but we can get more
!--- sophisticated using the distance between parallel faces later
if we need to...mfc
!--- sphere is calculated in ew_box.f and used here.
if( (xbnd1-xbnd0 > two*sphere) .or. (ybnd1-ybnd0 > two*sphere) .or.
(zbnd1-zbnd0 > two*sphere) )then
write(6,*) " ****************************************************"
write(6,*) " ERROR: QM region + cutoff larger than box dimension:"
write(6,'(2X,"QM-MM Cutoff = ",f8.4)') qmmm_nml%qmcut
write(6,*) " Coord Lower Upper Size Radius of largest
sphere inside unit cell"
write(6,'(5X,"X",4(2X,f8.3))') xbnd0, xbnd1, xbnd1-xbnd0, sphere
write(6,'(5X,"Y",4(2X,f8.3))') ybnd0, ybnd1, ybnd1-ybnd0, sphere
write(6,'(5X,"Z",4(2X,f8.3))') zbnd0, zbnd1, zbnd1-zbnd0, sphere
write(6,*) " ****************************************************"
call sander_bomb("QM_CHECK_PERIODIC<qm_mm.f>", &
"QM region + cutoff larger than box", &
"cannot continue, need larger box.")
endif
...
The code then basically calculates the largest sphere that will fit inside
the box and then checks if the QM region + MM atoms within qmcut would fit
within the diameter of this sphere ((2 x sphere) which is the radius). In
the case of Kaustubh calculation when it fails he has:
ERROR: QM region + cutoff larger than box dimension:
QM-MM Cutoff = 4.0000
Coord Lower Upper Size Radius of largest sphere inside unit
cell
X -9.745 9.611 19.357 11.528
Y -5.919 11.077 16.996 11.528
Z -12.188 11.060 23.248 11.528
****************************************************
So here his QMatoms + MM atoms extend from -9.745 to +9.611 in the X
direction, a distance of 19.4 angstroms. The largest sphere to fit in his
box is 11.528 x 2 = 23 angstroms diameter so he is fine in the X direction,
same in the Y direction, but the Z direction he is larger than the sphere,
so this violates the minimum image convention and the code quits. Simple
answer is that his box is too small.
> confused
> by how the upper/lower values in the sander bomb message are able to
> differ
> by their sign when visualization of the trajectory appears to show the
> QM
> atoms quite close together (maybe I used some really misleading
> re-wrapping?).
Visualization is in the original coordinate frame. The QM calculation is
done with the origin shifted to be that of the first QM atom. If you want to
visualize with respect to the coordinates used internally in the QM code
then run the trajectory through ptraj and have reimage for you with origin
center being defined as the first atom of your QM region.
> I'm not experiencing my old problems any more, but I
> still
> don't know what exactly caused them other than, apparently, triclinic
> boundary conditions.
Truncated octahedron boxes could cause issues. The original code that I
wrote only supported alpha=beta=gamma = 90 degrees. Non-orthogonal box
support was added with this commit:
commit 38e98745a98d3cc5fdb73d9d30177b7987464e1d
Author: rcw <rcw>
Date: Sat Aug 18 01:58:37 2007 +0000
RCW: Add support for non orthogonal boxes in QMMM.
By me of all people but that was a pretty hectic time for me so it is
possible there are issues with it, especially as a system approaches the
edge of the box. If you have a good test case to show the issue then I can
try to look into what is going wrong. It would good to know if you only see
problems with non-orthogonal boxes though.
All the best
Ross
/\
\/
|\oss Walker
---------------------------------------------------------
| Assistant Research Professor |
| San Diego Supercomputer Center |
| Adjunct Assistant Professor |
| Dept. of Chemistry and Biochemistry |
| University of California San Diego |
| NVIDIA Fellow |
|
http://www.rosswalker.co.uk |
http://www.wmd-lab.org/ |
| Tel: +1 858 822 0854 | EMail:- ross.rosswalker.co.uk |
---------------------------------------------------------
Note: Electronic Mail is not secure, has no guarantee of delivery, may not
be read every day, and should not be used for urgent or sensitive issues.
_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber
Received on Thu Oct 06 2011 - 10:30:05 PDT
