Ross,
I can only confirm that I have never experienced a (unexplainable) sander
bomb while using an orthorhombic cell.
Brian
On Thu, Oct 6, 2011 at 1:11 PM, Ross Walker <ross.rosswalker.co.uk> wrote:
> > 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
>
--
================================ Current Address =======================
Brian Radak : BioMaPS
Institute for Quantitative Biology
PhD candidate - York Research Group : Rutgers, The State
University of New Jersey
University of Minnesota - Twin Cities : Wright-Rieman Hall 101
Graduate Program in Chemical Physics : 610 Taylor Road,
Department of Chemistry : Piscataway, NJ
08854-8066
radak004.umn.edu :
radakb.biomaps.rutgers.edu
====================================================================
Sorry for the multiple e-mail addresses, just use the institute appropriate
address.
_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber
Received on Thu Oct 06 2011 - 11:00:04 PDT
