[CP2K:9624] The charge set for QM/MM MD simulation

[Teodoro Laino]

Dear Junbo,

first of all many thanks for having posted your request here: more people will be able to benefit of the on-going discussion.

Let me explain my point of view regarding your problem, with the understanding that they are very much in-line with Juerg’s comment.

First of all as you already realised, the ADBF QM/MM needs to have a constant charge and the charge is not allowed to change during the run. There is a very good reason for that: if you change the charge of the QM region you are changing the zero reference for your calculations during MD and this makes the analysis of your simulations less trivial.

Therefore, this is not about an "unimplemented feature” (*ignore the recommandation received in a separate list, to implement this functionality*), rather a protection to prevent users from running non-sense simulations.

Once we agree that we need to use the ADBF system with this “limitation", let’s move on to your questions.

(1) QM/MM runs do not have to be necessarily neutral. You can run QM/MM of QM charged systems (both fully periodic and non-periodic) being extremely careful to how you are decoupling the periodic images (in case of fully non-periodic calculations). You may have the corresponding counterion in the classical (MM) portion and the interaction would still be computed correctly (again be careful to decouple the QM part properly), in case QM+MM is neutral. You could also have the QM+MM charged and this would still work ok.

Option 1: If you believe that the Cl- is not important for your simulation you can freely decide to remove it from the entire system and run instead a simulation that is overall (QM+MM) charged (2+).

Option 2: If you decide to treat the Cl- classically, be sure that the Cl- is always far from your Zn2+ center enough (at least 2 solvation shells) to avoid any interaction between the two ions that could not be treated classically (electrostatic+non-bonded). In this case your overall system (QM+MM) is neutral, with your QM region being 2+ and the MM one 2-. Please, use a decoupling to be sure that the strong electrostatic coupling between the QM periodic replicas is dumped. 

(2) Using ADBF QM/MM for both option 1 and 2 is trivial as your QM charge is not changing (and again, for fundamental reasons it should not!).
Option 3: If you start with Option 2, you may notice instead that Cl- may play a fundamental role in your simulations or you may already know before starting your simulations that you are interested in mapping the interaction between the Cl- and  Zn2+.
In this case, there is only one thing that you can do. If you believe that the Cl- is an important part of the QM center, then you’ll have to include the Cl- into the QM region since the very beginning of your simulations. This will have the major drawback that Cl- may leave the Zn2+ center with the drawback of being forced to use very large QM boxes. 
Still, there is nothing lost or wasted in life (and in simulations too): you can apply a physical restraint (i.e. spherical wall) to prevent the Cl- from flying far and subsequently  reweight the calculated free energy (assuming you are interested in free energies) for the additional constraint. This will help you in containing the size of  your QM box even for extended simulation times (several hundreds of ps).

Personally, I like option 3.. but depending on your understanding of the system you may decide (equally good) to go for option 2 or option 1. … Possibly even for another option not necessarily considered here….

Hope this helps and good luck!
Best,
Teo



> On 4 Nov 2017, at 17:07, luju... at gmail.com wrote:
> 
> Hi, I want to do a simulation of hydration of Zn2+ by QM/MM MD of CP2K. But there are several problems for my simulation:
> (1). I do not know whether we should neutralize the system by minus ion in CP2K, such as Cl-. We can also not neutralize the system, is it right?
> (2). If we add some minus ions into system, how can we define the charge of dynamic QM region in adbf QM/MM MD calculation? When a charged molecule enter the QM region, such as Cl- ion, the total charge of dynamic QM region will change. I have read the CP2K manual, I find it can not change the charge of QM region when doing MD. So this is a problem.
> Can you help me for these two problems? I will be very appreciated!
> With my best regards,
> Junbo
>  <>
> 
> -- 
> You received this message because you are subscribed to the Google Groups "cp2k" group.
> To unsubscribe from this group and stop receiving emails from it, send an email to cp2k+uns... at googlegroups.com <mailto:cp2k+uns... at googlegroups.com>.
> To post to this group, send email to cp... at googlegroups.com <mailto:cp... at googlegroups.com>.
> Visit this group at https://groups.google.com/group/cp2k <https://groups.google.com/group/cp2k>.
> For more options, visit https://groups.google.com/d/optout <https://groups.google.com/d/optout>.

-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://lists.cp2k.org/archives/cp2k-user/attachments/20171108/450364a2/attachment.htm>