cgcmm.angle¶
Overview
cgcmm.angle.cgcmm |
Details
CGCMM angle potentials.
-
class
hoomd.cgcmm.angle.
cgcmm
¶ CGCMM angle potential.
The command angle.cgcmm defines a regular harmonic potential energy between every defined triplet of particles in the simulation, but in addition in adds the repulsive part of a CGCMM pair potential between the first and the third particle.
B. Levine et. al. 2011 describes the CGCMM implementation details in HOOMD-blue. Cite it if you utilize the CGCMM potential in your work.
The total potential is thus:
\[V(\theta) = \frac{1}{2} k \left( \theta - \theta_0 \right)^2\]where \(\theta\) is the current angle between the three particles and either:
\[V_{\mathrm{LJ}}(r_{13}) -V_{\mathrm{LJ}}(r_c) \mathrm{~with~~~} V_{\mathrm{LJ}}(r) = 4 \varepsilon \left[ \left( \frac{\sigma}{r} \right)^{12} - \left( \frac{\sigma}{r} \right)^{6} \right] \mathrm{~~~~for~} r <= r_c \mathrm{~~~} r_c = \sigma \cdot 2^{\frac{1}{6}}\]\[V_{\mathrm{LJ}}(r_{13}) -V_{\mathrm{LJ}}(r_c) \mathrm{~with~~~} V_{\mathrm{LJ}}(r) = \frac{27}{4} \varepsilon \left[ \left( \frac{\sigma}{r} \right)^{9} - \left( \frac{\sigma}{r} \right)^{6} \right] \mathrm{~~~~for~} r <= r_c \mathrm{~~~} r_c = \sigma \cdot \left(\frac{3}{2}\right)^{\frac{1}{3}}\]\[V_{\mathrm{LJ}}(r_{13}) -V_{\mathrm{LJ}}(r_c) \mathrm{~with~~~} V_{\mathrm{LJ}}(r) = \frac{3\sqrt{3}}{2} \varepsilon \left[ \left( \frac{\sigma}{r} \right)^{12} - \left( \frac{\sigma}{r} \right)^{4} \right] \mathrm{~~~~for~} r <= r_c \mathrm{~~~} r_c = \sigma \cdot 3^{\frac{1}{8}}\]with \(r_{13}\) being the distance between the two outer particles of the angle.
Coefficients:
- \(\theta_0\) - rest angle
t0
(in radians) - \(k\) - potential constant
k
(in units of energy/radians^2) - \(\varepsilon\) - strength of potential
epsilon
(in energy units) - \(\sigma\) - distance of interaction
sigma
(in distance units)
Coefficients \(k, \theta_0, \varepsilon\), and \(\sigma\) and Lennard-Jones exponents pair must be set for each type of angle in the simulation using
set_coeff()
.-
disable
(log=False)¶ Disable the force.
Parameters: log (bool) – Set to True if you plan to continue logging the potential energy associated with this force. Examples:
force.disable() force.disable(log=True)
Executing the disable command will remove the force from the simulation. Any
hoomd.run()
command executed after disabling a force will not calculate or use the force during the simulation. A disabled force can be re-enabled withenable()
.By setting log to True, the values of the force can be logged even though the forces are not applied in the simulation. For forces that use cutoff radii, setting log=True will cause the correct r_cut values to be used throughout the simulation, and therefore possibly drive the neighbor list size larger than it otherwise would be. If log is left False, the potential energy associated with this force will not be available for logging.
-
get_energy
(group)¶ Get the energy of a particle group.
Parameters: group ( hoomd.group
) – The particle group to query the energy for.Returns: The last computed energy for the members in the group. Examples:
g = group.all() energy = force.get_energy(g)
-
get_net_force
(group)¶ Get the force of a particle group.
Parameters: group ( hoomd.group
) – The particle group to query the force for.Returns: The last computed force for the members in the group. Examples
g = group.all() force = force.get_net_force(g)
-
get_net_virial
(group)¶ Get the virial of a particle group.
Parameters: group ( hoomd.group
) – The particle group to query the virial for.Returns: The last computed virial for the members in the group. Examples
g = group.all() virial = force.get_net_virial(g)
-
set_coeff
(angle_type, k, t0, exponents, epsilon, sigma)¶ Sets the CG-CMM angle coefficients for a particular angle type.
Parameters: - angle_type (str) – Angle type to set coefficients for
- k (float) – Coefficient \(k\) (in units of energy/radians^2)
- t0 (float) – Coefficient \(\theta_0\) (in radians)
- exponents (str) – is the type of CG-angle exponents we want to use for the repulsion.
- epsilon (float) – is the 1-3 repulsion strength (in energy units)
- sigma (float) – is the CG particle radius (in distance units)
Examples:
cgcmm.set_coeff('polymer', k=3.0, t0=0.7851, exponents=126, epsilon=1.0, sigma=0.53) cgcmm.set_coeff('backbone', k=100.0, t0=1.0, exponents=96, epsilon=23.0, sigma=0.1) cgcmm.set_coeff('residue', k=100.0, t0=1.0, exponents='lj12_4', epsilon=33.0, sigma=0.02) cgcmm.set_coeff('cg96', k=100.0, t0=1.0, exponents='LJ9-6', epsilon=9.0, sigma=0.3)
- \(\theta_0\) - rest angle