Questions about AdBF-QM/MM in cp2k

Ivan Gladich

unread,

Feb 2, 2017, 12:59:23 AM2/2/17

to cp2k

Dear Cp2K users,

I am a newbie of CP2K and I am trying to use Adaptive Buffered force QM-MM (AdBF-QM/MM) in cp2k

https://www.ncbi.nlm.nih.gov/pmc/articles/PMC4351341/

to study a radical specie in water

I have a couple of questions.

1) Using

...

&QMMM

&QM_KIND Br

MM_INDEX 1

&END QM_KIND

&QM_KIND O1

MM_INDEX 2

&END QM_KIND

&FORCE_MIXING

MAX_N_QM 500

R_CORE 0.0 0.1 #dimension of the core region for the reduced QM/MM calculation. Could be 0 if you have the force field

R_QM 3.0 3.5 #Dimension of the dynamical core region

R_BUF 3.0 3.5 #dimensione of the buffere region

QM_KIND_ELEMENT_MAPPING Br Br #atoms that can be adaptively selectefor the QM region

QM_KIND_ELEMENT_MAPPING O O #atoms that can be adaptively selected for the QM region

QM_KIND_ELEMENT_MAPPING H H #atoms that can be adaptively selectedfor the QM region

EXTENDED_DELTA_CHARGE 0 # this is the charge of the QM_dynamical + buffer region.

MOMENTUM_CONSERVATION_REGION QM #only to dynamical core, you can also choose BUFFER for QM+BUFFER

MOMENTUM_CONSERVATION_TYPE EQUAL_A #equal acceleration on each atom, or EQUAL_F, equal force

&END FORCE_MIXING

....

And I got in my output

QMMM FORCE MIXING final count (not including links): N_QM core_list 2 N_QM core 0 N_QM extended 12 N_QM buffered 63

...

Is N_QM core = 0 the number of atoms in my dynamical QM region?

Moreover, reading the paper N_QM extended "is constructed adding a buffer region around the QM dynamical core"

But N_QM extended =12 and N_QM buffered is 63, so it cannot be (N_QM should be the same of the buffer plus the dynamical core)...

2) How does the R_BUF and R_QM works?

R_BUF is a radius around each atoms belonging to the QM dynamical core region

OR

R_BUF is centered around the QM_CORE (in my case determined by Br and O position) and, thus, I should impose always R_QM =< R_BUF ?

3) Is it possible to write, at each time step of my run, a label or the indexes of the atoms in the dynamical and extended QM/MM region?

I tried using

&FORCE_MIXING_LABELS

FORMAT XYZ

&EACH

MD 1

&END EACH

COMMON_ITERATION_LEVELS 0

&END FORCE_MIXING_LABELS

but it seems that it labels only the atoms in the QM CORE, i.e., *fmlabels-1_23.xyz

6482

i = 23, time = 11.500, E = -21.9207476721

BR 10.0000000000 10.0000000000 10.0000000000

O 10.0000000000 10.0000000000 10.0000000000

O 0.0000000000 0.0000000000 0.0000000000

H 0.0000000000 0.0000000000 0.0000000000

...

Thank you for any possible help or comment

Best Regards

Ivan

Noam Bernstein

unread,

Feb 2, 2017, 3:11:57 PM2/2/17

to cp2k

No, I apologize for the somewhat confusing output. "N_QM core" is the
number added to the core list by the action of R_CORE. I.e. The core
region consists of a fixed list (the atoms specified by the
conventional QM/MM directives) + an additional region, which is every
atom within R_CORE of the fixed list. The core region (combining
those two, the list and the R_CORE based additions) is QM in both
reduced and extended calculations whose forces will be mixed. The QM
dynamics region includes all those core atoms, plus whatever is
included by R_QM (based on distance from core atoms).

>
> Moreover, reading the paper N_QM extended "is constructed adding a buffer
> region around the QM dynamical core"
>
> But N_QM extended =12 and N_QM buffered is 63, so it cannot be (N_QM should
> be the same of the buffer plus the dynamical core)...

I think here the notation of the code and the paper are unfortunately
not quite consistent. N_QM extended refers to what the paper calls
the dynamical QM region (I'm pretty sure it's actually just the part
that's in addition to the core).

It doesn't help that there are dynamically selected atoms in the core
region (as well as in the QM dynamics region and buffer region), but
that shouldn't be confused with the region that follows QM force
dynamics.

>
>
> 2) How does the R_BUF and R_QM works?
>
>
> R_BUF is a radius around each atoms belonging to the QM dynamical core
> region
>
>
> OR
>
>
> R_BUF is centered around the QM_CORE (in my case determined by Br and O
> position) and, thus, I should impose always R_QM =< R_BUF ?

More like the former, but not exactly, at least not if you're using
the terminology of either the paper or the code. As the code
documentation says, it's applied around every atom that should follow
QM forces, i.e. the core list + the core distance (R_CORE) based + the
QM dynamical region distance (R_QM) based.

>
>
> 3) Is it possible to write, at each time step of my run, a label or the
> indexes of the atoms in the dynamical and extended QM/MM region?
>
>
> I tried using
>
>
> &FORCE_MIXING_LABELS
>
> FORMAT XYZ
>
> &EACH
>
> MD 1
>
> &END EACH
>
> COMMON_ITERATION_LEVELS 0
>
> &END FORCE_MIXING_LABELS
>
>
> but it seems that it labels only the atoms in the QM CORE, i.e.,
> *fmlabels-1_23.xyz
>
>
> 6482
>
> i = 23, time = 11.500, E = -21.9207476721
>
> BR 10.0000000000 10.0000000000 10.0000000000
>
> O 10.0000000000 10.0000000000 10.0000000000
>
> O 0.0000000000 0.0000000000 0.0000000000
>
> H 0.0000000000 0.0000000000 0.0000000000
>
> H 0.0000000000 0.0000000000 0.0000000000
>

Good question. I'll take a look at the source and remind myself what
kind of functionality is available.

Noam

Noam Bernstein

unread,

Feb 2, 2017, 3:40:31 PM2/2/17

to cp2k

I looked at the source, and I see no reason why this should be the
case. As far as I can tell it should write all the labels, values
4-10. How sure are you that there aren't some atoms with other label
settings (> 0)? I think I'd need a full example (all the input files)
that shows this problem to usefully investigate further. Ideally one
that runs quickly.

Noam

Rolf David

unread,

Feb 3, 2017, 5:41:34 AM2/3/17

to cp2k

For the force-mixing labels, for me it works

I have 10.0 for list/core, 7.0 for extended and 5.0 for buffered. (for CP2K version 3.0)

I have a parell question if I may:

How can we restart a ad-qmmm calculations ?:

I run a first dynamics, at the last step at the end i got:

QMMM FORCE MIXING final count (not including links): N_QM core_list 11 N_QM core 0 N_QM extended 42 N_QM buffered 111

And restarting it, the very first step :

QMMM FORCE MIXING final count (not including links): N_QM core_list 11 N_QM core 0 N_QM extended 33 N_QM buffered 87

+ the WFN: READ RESTART : WARNING : DIFFERENT natom, returning 131 164

So for the minimal QM region, I got only core_list, so it doesen't change and it can read the WFN restart files

But for the extended, since it is buffer in buffer out related, it varies.

Is there a way to let CP2K now which atoms where in the extended/buffered ? (Except putting the min/max the same for R_*** distances)

For info, these test were run on the 2.x (rev 16205) and i used:

R_CORE [angstrom] 0.0 0.0

R_QM [angstrom] 2.5 3.0

R_BUF [angstrom] 2.5 3.0

A part for the restart issues, which can be penalising when using national/eu clusters where there is a wall time, it runs smoothly

Thanks

Ivan Gladich

unread,

Feb 5, 2017, 1:00:13 AM2/5/17

to cp...@googlegroups.com

Dear Noam,

Dear Rolf,

thank you for your reply, they helped me to better understand the meaning of R_QM, R_CORE and EXTENDED and the differences with to the paper.

Concerning the labels I still have some unclear point. For example, I got for the last step

> QMMM FORCE MIXING final count (not including links): N_QM core_list 2 N_QM core 0 N_QM extended 12 N_QM buffered 66

but, if I look at my non-zero labels I see

6482

i = 50, time = 25.000, E = -108.8420993135

BR 10.0000000000 10.0000000000 10.0000000000

O 10.0000000000 10.0000000000 10.0000000000

O 5.0000000000 5.0000000000 5.0000000000

H 5.0000000000 5.0000000000 5.0000000000

O 7.0000000000 7.0000000000 7.0000000000

H 7.0000000000 7.0000000000 7.0000000000

O 5.0000000000 5.0000000000 5.0000000000

H 5.0000000000 5.0000000000 5.0000000000

O 5.0000000000 5.0000000000 5.0000000000

H 5.0000000000 5.0000000000 5.0000000000

O 5.0000000000 5.0000000000 5.0000000000

H 5.0000000000 5.0000000000 5.0000000000

So, two questions

a) 10 should be the core list atoms but which is the meaning of labels 7 and 5?

b) I cannot match the numbers because the labels tells me that I have 17 QM atoms in the extended or extended plus buffer.

The Output points toward much higher number (66)

You can find my inputs at the following link

goo.gl/C15cBY

Thank you again for your help, it was extremely useful.

Ivan

--
You received this message because you are subscribed to a topic in the Google Groups "cp2k" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/cp2k/u4xcztECsE0/unsubscribe.
To unsubscribe from this group and all its topics, send an email to cp2k+unsubscribe@googlegroups.com.
To post to this group, send email to cp...@googlegroups.com.
Visit this group at https://groups.google.com/group/cp2k.
For more options, visit https://groups.google.com/d/optout.

Noam Bernstein

unread,

Feb 6, 2017, 4:17:55 PM2/6/17

to cp2k

On Sun, Feb 5, 2017 at 1:00 AM, Ivan Gladich <igla...@gmail.com> wrote:
> Dear Noam,

>
>
> So, two questions
>
> a) 10 should be the core list atoms but which is the meaning of labels 7 and
> 5?

10.0 is core list
7.0 is QM dynamics (what's called, misleadingly, "N_QM extended" above),
5.0 is buffer

> b) I cannot match the numbers because the labels tells me that I have 17 QM
> atoms in the extended or extended plus buffer.
> The Output points toward much higher number (66)

If those lines are the only ones that are non-zero, then I agree
there's a mismatch, but it's hard for me to imagine how that happens,
since the files are written from the same internal data structures
that are used to generate the counts that are printed above. I'll try
to run with your input files.

Noam

Noam Bernstein

unread,

Feb 7, 2017, 11:06:14 AM2/7/17

to cp2k

On Mon, Feb 6, 2017 at 4:17 PM, Noam Bernstein
<no...@bollweevil.gdbg.org> wrote:
>
> If those lines are the only ones that are non-zero, then I agree
> there's a mismatch, but it's hard for me to imagine how that happens,
> since the files are written from the same internal data structures
> that are used to generate the counts that are printed above. I'll try
> to run with your input files.

I ran with the input files on google drive (except only 5 time steps),
and do not see any problems. My output is entirely consistent:

tin 1203 : fgrep N_QM stdout