[CP2K-user] [CP2K:21170] Re: CP2K: How to Define a Molecule and Restrict Lennard-Jones Interactions Within It?

Marcella Iannuzzi marci.akira at gmail.com
Tue Feb 18 14:34:47 UTC 2025


Hi 

If that is the problem , a straightforward way would be to pass the 
connectivity information through a connectivity file. 
The connectivity generated by CP2K can also be tuned from 
- FORCE_EVAL <https://manual.cp2k.org/trunk/CP2K_INPUT/FORCE_EVAL.html> 
- SUBSYS <https://manual.cp2k.org/trunk/CP2K_INPUT/FORCE_EVAL/SUBSYS.html> 
- TOPOLOGY 
<https://manual.cp2k.org/trunk/CP2K_INPUT/FORCE_EVAL/SUBSYS/TOPOLOGY.html> 
- GENERATE
 
Coordinates passed in PDB format, where molecule name and molecule number 
are specified, might also work. 

Regards
Marcella

On Tuesday, February 18, 2025 at 3:06:02 PM UTC+1 Planck wrote:

> Hi, 
>                I am trying to simulate molecular dynamics for a 
> hypothetical chain of carbon atoms (this is a minimum working example 
> designed to illustrate a larger issue I am facing). Each molecule consists 
> of 5 carbon atoms, and my simulation contains 4 such molecules.
>
> ```
>   20
>   #Built with Packmol
>   C          264.316076      114.220486       50.721755
>   C          264.392382      114.277533       50.950866
>   C          264.525927      114.377373       51.351839
>   C          264.659472      114.477213       51.752813
>   C          264.735777      114.534260       51.981923
>   C          327.956067      237.050267       49.303207
>   C          327.861562      237.033691       49.074379
>   C          327.696165      237.004681       48.673899
>   C          327.530769      236.975671       48.273418
>   C          327.436264      236.959095       48.044590
>   C          433.879103      142.141879      115.467902
>   C          433.814625      141.918476      115.381286
>   C          433.701780      141.527492      115.229696
>   C          433.588936      141.136507      115.078106
>   C          433.524458      140.913105      114.991490
>   C          234.869472      121.984410      116.856352
>   C          234.879704      121.977472      116.608530
>   C          234.897610      121.965330      116.174809
>   C          234.915516      121.953188      115.741088
>   C          234.925748      121.946251      115.493266
> ```
> I have defined a Lennard-Jones potential between the carbon atoms. My 
> input file looks like this:
>
> ```
> &GLOBAL
>   PRINT_LEVEL HIGH ! Low verbosity
>   PROJECT MD_simulation
>   RUN_TYPE MD
> &END GLOBAL
>
> &MOTION
>   &MD
>     ENSEMBLE NVT ! Constant volume
>     TEMPERATURE 420 
>     TIMESTEP 4.0 ! 4 fs
>     STEPS 100000 ! 
>
>     &THERMOSTAT
>       TYPE CSVR
>     &END THERMOSTAT
>     
>     &PRINT
>   &ENERGY
>           &EACH
>           MD 25 ! 25*TIMESTEP
>     &END EACH
>   &END ENERGY
>           &PROGRAM_RUN_INFO
>           &EACH
>             MD 25 ! 25*TIMESTEP
>     &END EACH
>   &END PROGRAM_RUN_INFO
>     &END PRINT
>   &END MD
>   &PRINT
>     &TRAJECTORY
>     &EACH
>       MD 25 ! 25*TIMESTEP
>       &END EACH
>     &END TRAJECTORY
>   &END PRINT
> &END MOTION
>
> &FORCE_EVAL
>   METHOD FIST
>   &MM
>     &FORCEFIELD
>       DO_ELECTROSTATICS F ! Ignore electrostatics (potential from QC)
>       IGNORE_MISSING_CRITICAL_PARAMS T     ! Ignore bonded parameters 
> (potential from QC)
>
>       &SPLINE
>         EMAX_SPLINE 5000
>         EPS_SPLINE 1.00000000E-8 ! Increase spline accuracy
>       &END SPLINE
>       
>       &NONBONDED
>
>         &LENNARD-JONES
>           atoms C C
>           EPSILON [kcalmol]  0.08
>           SIGMA   [angstrom] 3.1507
>           RCUT    [angstrom] 11.4
>           RMIN 1.0
>         &END LENNARD-JONES
>
>       &END NONBONDED
>
>
>     &END FORCEFIELD
>   
>     &POISSON
>       PERIODIC NONE
>       &EWALD
>         EWALD_TYPE NONE
>       &END EWALD
>     &END POISSON
>   
>   &END MM
>
>   &SUBSYS
>   
>     &CELL
>       PERIODIC XYZ
>       ABC 800 800 800
>     &END CELL 
>
>     &TOPOLOGY
>       COORD_FILE_FORMAT XYZ
>       COORD_FILE_NAME  C5_20.xyz
>     &END TOPOLOGY
>   
>   &END SUBSYS
>   
>   STRESS_TENSOR NUMERICAL
> &END FORCE_EVAL
>
> ```
>
>
>  but when I run the simulation, I encounter the following error:
> ```
>  SPLINE_INFO| Generating 1 splines for NONBONDED interactions 
>               Due to 1 different atomic kinds
>
>
>  *******************************************************************************
>  *   ___                                                                   
>     *
>  *  /   \                                                                 
>      *
>  * [ABORT]  SPLINE_INFO| Number of points:  2346797 obtained accuracy 
>  43.1289 *
>  *  \___/         . MM SPLINE: no convergence on required accuracy (adjust 
>     *
>  *    |                            EPS_SPLINE and rerun)                   
>     *
>  *  O/|                                                                   
>      *
>  * /| |                                                                   
>      *
>  * / \                                                   
>  pair_potential.F:426 *
>
>  *******************************************************************************
>
>
>  ===== Routine Calling Stack ===== 
>
>             7 spline_nonbond_control
>             6 force_field_pack_splines
>             5 force_field_pack
>             4 force_field_control
>             3 fist_init
>             2 fist_create_force_env
>             1 CP2K
> ```
> I suspect this is happening because I did not explicitly specify that each 
> 5-atom unit is a single molecule, so CP2K is applying the Lennard-Jones 
> potential between atoms that should be bonded. Since the atoms within each 
> molecule are very close together, CP2K treats them as interacting via 
> nonbonded forces, which leads to numerical instability in the Lennard-Jones 
> calculations.
>
> **If my assumption is correc**, How can I properly define a molecule in 
> CP2K using the 5 carbon atoms?
> How can I exclude Lennard-Jones interactions **within the same molecule**, 
> while keeping them for interactions between different molecules?
>
> N.B:
> I understand that the carbon atoms are very close together, but this is a 
> minimum working example meant to illustrate the issue.
> I cannot change the atomic positions, as the actual problem requires 
> keeping them fixed. Also I want to ignore the bond vibrations, rotations, 
> torsions etc. 
>
> Any help or guidance would be greatly appreciated!
>
>

-- 
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+unsubscribe at googlegroups.com.
To view this discussion visit https://groups.google.com/d/msgid/cp2k/c8965abb-ad8a-41c7-98b1-ca8674207b75n%40googlegroups.com.
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://lists.cp2k.org/archives/cp2k-user/attachments/20250218/88ec2c4b/attachment-0001.htm>


More information about the CP2K-user mailing list