FlucDens

class AmberFD.FlucDens(*args)

This class impliments the fluctuating monopole density interaction. Similar to a fluctuating charge model, this force uses atom centered monopole densities described by slater functions for the electrons and point-particles for the nuclei. As such, FlucDens includes effects of electron overlap between atomic sites.

To use this force, create a FlucDens force object by specifying the number of charge sites and the parameters for all sites. These parameters can be altered after the creation of the force object. Next, each site must be assigned a molecular “fragment” Each fragment, defined by it’s site indicies, are groups of charge sites whose dynamic desity coefficients will sum to zero when the electrostatic energy is minized. As a result, charge transfer is not allowed between fragments, but each fragment will however reduce to it’s frozen electron populations when complety separated from all other framgents. Each fragment is also used to define which frozen density - dynamic density interactions are included. Frozen densities on each site do not interact with the dynamic densities on each fragment, allowing for this limmiting behavior to occur.

Next, set the model’s pther settings by calling their respective setter function. This includes the polarization dampening strengths, charge-transfer coefficients, and additional hardness values to be applied to the Coulomb matrix. Polarization dampening is applied between a polarized density site and it’s interacing frozen density site and decays exponentially as a function of their distance. The strength and rate of this polarization can be adjusted, as well as whether this is applied to the linear (default) or the quadratic matrix terms (old method). Charge transfer is not explicitly allowed in this force, and instead the energy associated with this interation is estimated by a scalar multiple of the total polarization energy.

Additional hardness values can be applied to the Coulomb matrix. These values are not overrides, but instead add a constant value to the self-interaction between two densities and hinder (or enhance) the ability of that site to gain or loose electrons. These values can be either positive or negative, with positive values increasing the hardness. Care must be taken that the system remains positive-definite if negative values are used.

Aconstant electric field can also be applied to each system. Simply call set_external_field() with the desired field components to apply this field. Periodic boundary conditions are also possible, and can be set by calling set_use_PBC(). If PBC is used, then a long-range cutoff distance must be used (the default is to use half the minimum periodic box lengths). Turning this off and using PBC will have undefined behavior.

This module also computes all interaction energies, including minimizing the total energy w.r.t. the dynamic density coefficients. Do run this computation, first call calc_energy() to compute all frozen-frozen energy components and to create the polarization matrix components. Then, call solve_minimization() to perform the actual matrix algebra associated with the minimization procedures. Energy components can then be returned calling get_energies() or any of the particular energy getter functions. Dipoles (including the polarized density) can also then be called with the updated electron populations.

All values in this module use atomic units unless otherwise specified.

To construct a new FlucDens force object, define the object with

Parameters

n_sites – number of particles (including cirtual sites)
frozen_charges_in – the frozen charges of all sites
nuclei_in – the nuclei of each site. Core electrons will be remove automatically
frozen_exp – the Slater exponents of all frozen densities
dynamic_exp – the Slater exponents of all dynamic densities

__init__(*args)

Methods

`__init__`(*args)
`add_del_frz_exclusion`(delta_i, frz_j)	Exclude the frozen electrostatic interactions between a dynamic and frozen site
`add_fragment`(site_idx_list)	Create a new fragment whoes dynamic density whould sum to zero.
`add_frz_frz_exclusion`(frz_i, frz_j)	Exclude the interaction of two frozen densities
`apply_field_to_system`(coords)	Apply the field to the system.
`calc_density`(points, pos, density_type)	Compute the electron density at the specified coordinates
`calc_energy`(coords[, calc_frz, calc_pol])	Calculate the energy and forces of the system.
`calc_frz_ext_field_energy`(coords, forces)	Compute the energy and forces of the frozen charges (nuclei and frozen electrons) with the external electric field
`calc_one_electro`(deltaR, i, j, calc_pol, ...)	Compute one pair of interactions between two atomic sites
`calc_one_frozen`(coords, i, j)	Compute one pair of frozen interactions between two atomic sites (does not include polarization)
`calc_overlap`(coords)	Calculate the total electron overlap
`create_frz_exclusions_from_bonds`(bonds, ...)	Create frozen-frozen exclusions based on each site's bonds.
`dot3Vec`(coords, i, j)	Perform the dot product between two coordinates
`elec_elec_energy`(inv_r, a, b, exp_ar, ...)	Compute the density-density interaction, either frozen or dynamic
`elec_nuclei_energy`(inv_r, a, exp_ar, dEdR)	Compute the energy between a nuclei and an electron density
`frz_frz_overlap`(inv_r, a, b, exp_ar, exp_br)	Overlap between frozen-frozen dsnsity
`get_constraints`()	Get constraint matrix applied to the dyanmic densities
`get_ct_coeff`()	Get the charge-transfer coefficient
`get_ct_energy`()	Get the estimated CT energy as a multiple of the polarization energy
`get_cutoff_distance`()	Get the cutoff distance (in bohr)
`get_dampening`()	Return the polarization dampening parameters
`get_del_frz_exclusions`(particle1)	Get the indicies of all frozen sites that do not interact with a single dynamic density
`get_delta_rho`()	Get the delta rhos of each site that are computed after minimizing the energy
`get_dipole`(coords, density_dype)	Compute the dipole for a particule type of charge
`get_dipoles`(coords)	Return the electric dipoles of the system.
`get_energies`()	Return all energy components computed
`get_forces`()	Get the forces exerted on all particles after an energy computation is performed
`get_fragments`()	Return the partitioned fragments and their site indicies for each
`get_frozen_energy`()	Get the energy associated with frozen-frozen dsitity interactions
`get_frz_frz_exclusions`(particle1)	Get the all sites that are excluded from interacting from a particular site
`get_num_constraints`()	Get the number of total constrains applied to the energy minimization
`get_num_fragments`()	Get the num fragments
`get_num_frz_frz_exclusions`()	Get the number of frz-frz exclusions
`get_param_names`()	Get the name of the parameters that can be return
`get_params_by_name`(param_name)	Get an array of parameters by it's name
`get_polarization_energy`()	Get the polarization energy
`get_rho_coulomb_mat`()	Get the delta-rho Coumomb J-matrix
`get_rho_pot_vec`()	Get the electroc potentia vector (linear term)
`get_site_params`(index)	Get the parameters of an atomic site
`get_total_time`()	Not implimented yet! Return the total wall time used for computations
`get_use_PBC`()	Inquire if PBC is used or not
`get_use_SR_cutoff`()	Inquire of SR cutoffs are being used
`get_use_cutoff`()	Inquire of nonbonded cutoffs are being used
`initialize_calculation`()	initizlize arrays before an energy computation is performed.
`print_params`(message, param_name)	Print site parameters to the console.
`set_additional_hardness`(*args)	Overload 1:
`set_calc_forces`(calculate_forces)	Set wether or not forces should also be computed.
`set_ct_coeff`(coeff)	Set the charge-transfer coefficient.
`set_cutoff_distance`(distance_in_nm)	Set the cutoff distance (in bohr)
`set_dampening`(*args)	Set the polarization dampening parameters
`set_dyn_exp`(*args)	Overload 1:
`set_external_field`(field_x, field_y, field_z)	Apply a constant electric field to the entire system
`set_frag_constraints`(constr_frags)	Choose whether constraints are used on all fragments or not.
`set_frz_exp`(index, value)	Set the Slater exponent of a particular frozen density site
`set_site_params`(index, frz_chg, frz_exp, dyn_exp)	Change the parameters for one charge site
`set_use_PBC`(*args)	Overload 1:
`set_use_SR_cutoff`(use_cutoff)	Set the use SR cutoff.
`set_use_cutoff`(useCutoff)	Set whether or not to use a nonbonded cutoff.
`solve_minimization`(forces)	Minimize the total energy of the system w.r.t.

Attributes

`A_mat_save`	The quadratic A-matrix that was used in the last energy minimization.
`All`	The total density of the system (nuclei, frozen, and dynamic)
`B_vec_save`	The linear B-vector used in the last energy minimization.
`Delta`	Dynamic electron density only
`Frozen`	Frozen electron density only
`Linear`	Dampening will be applied to the linear terms
`Nuclei`	Nuclei charge
`Quadratic`	Dampening will be applied to the quadratic terms (old method)
`thisown`	The membership flag

property thisown: The membership flag

Linear = 1: Dampening will be applied to the linear terms

Quadratic = 2: Dampening will be applied to the quadratic terms (old method)

All = 0: The total density of the system (nuclei, frozen, and dynamic)

Frozen = 1: Frozen electron density only

Delta = 2: Dynamic electron density only

Nuclei = 3: Nuclei charge

calc_density(points, pos, density_type)

Compute the electron density at the specified coordinates

Parameters

points (vec_d) – A 3K array of K number of points to evaluate the density at
pos (vec_d) – The 3N array of positions for the charge sites of the system
density_type (int) – The type of density to compute

Return type

vec_d

Returns

vec_d A K-length array of density evaluations

set_site_params(index, frz_chg, frz_exp, dyn_exp)

Change the parameters for one charge site

Parameters

index (int) – The index of the site to adjust the parameters
frz_chg (float) – The frozen charge to set
frz_exp (float) – The Slater exponent of the frozen density
dyn_exp (float) – The Slater exponent of the dynamic density

get_site_params(index)

Get the parameters of an atomic site

Parameters

index (int) – The index of the site to adjust the parameters
frz_chg (float) – The frozen charge of the site
frz_exp (float) – The Slater exponent of the frozen density
dyn_exp (float) – The Slater exponent of the dynamic density

set_dyn_exp(*args)

Overload 1:

Set the Slater exponent of a particular dynamic density site

Parameters

index (int) – The index of the site
value (float) – The value of the Slater exponent

Overload 2:

Set the Slater exponents of all dynamic density sites at once

Parameters: exponents (vec_d) – An array of all exponents in order of their indicies

set_frz_exp(index, value)

Set the Slater exponent of a particular frozen density site

Parameters

index (int) – The index of the site
value (float) – The value of the Slater exponent

set_ct_coeff(coeff)

Set the charge-transfer coefficient. Charge transfer is extimated as a multiple of the polarization energy between two inducable sites. This does not apply to the polarization due to external fields

Parameters: coeff (float) – The scalar multiple of the polarization energy used for the charge transfer energy

get_ct_coeff()

Get the charge-transfer coefficient

Return type: float
Returns: double The coefficient

set_dampening(*args)

Set the polarization dampening parameters

Parameters

coeff (float) – The strength of the dampening
exponent (float) – The exponential decay of the dampening
damp (int) – Tye type of dampening, enumerated by DampType

get_dampening()

Return the polarization dampening parameters

Parameters

coeff (float) – The strength of the dampening
exponent (float) – The exponential decay of the dampening

set_additional_hardness(*args)

Overload 1:

Add additional atomic hardness to the diagonal J-matrix elements. This is added to the value already computed by the self interaction in the matrix, not an override.

Parameters

index (int) – The particle to add the hardness to.
value (float) – The value of the additional hardness

Overload 2:

Return the additional hardness values of all sites

Parameters: values (vec_d) – An array of the hardness values

add_fragment(site_idx_list)

Create a new fragment whoes dynamic density whould sum to zero. Every site must be assigned to a fragment, else an exception will be thrown.

Parameters: site_idx_list (std::vector< int,std::allocator< int > >) – A list of site indicies that make up this fragment

get_fragments()

Return the partitioned fragments and their site indicies for each

Return type: std::vector< vec_i,std::allocator< vec_i > >
Returns: An array of fragment indicies. Each element correpsonds to a fragment and contains a list of indicies that make up that fragment

get_num_fragments()

Get the num fragments

Return type: int
Returns: The number of fragments

add_del_frz_exclusion(delta_i, frz_j)

Exclude the frozen electrostatic interactions between a dynamic and frozen site

Parameters

delta_i (int) – The index of the dynamic density site
frz_j (int) – The index of the frozen density site

get_del_frz_exclusions(particle1)

Get the indicies of all frozen sites that do not interact with a single dynamic density

Parameters: particle1 (int) – The index of the dynamic density
Return type: std::set< int,std::less< int >,std::allocator< int > >
Returns: std::set<int> The indicies of the frozen sites that do not interact with this dynamic density

set_frag_constraints(constr_frags)

Choose whether constraints are used on all fragments or not. NOT FULLY IMPLIMENTED YET!

Parameters: constr_frags (boolean) – If true, fragment constraints will be used. If false, then only the total dynamic density of the entire system will be summed to zero, not per fragment.

get_num_constraints()

Get the number of total constrains applied to the energy minimization

Return type: int
Returns: The total number of constraints

get_constraints()

Get constraint matrix applied to the dyanmic densities

Return type: std::vector< std::vector< int,std::allocator< int > >,std::allocator< std::vector< int,std::allocator< int > > > >
Returns: std::vector<std::vector<int>>

add_frz_frz_exclusion(frz_i, frz_j)

Exclude the interaction of two frozen densities

Parameters

frz_i (int) – The index of the first frozen site
frz_j (int) – The index of the second frozen site

get_frz_frz_exclusions(particle1)

Get the all sites that are excluded from interacting from a particular site

Parameters: particle1 (int) – The frozen site whoes exclusions are being inquired
Return type: std::set< int,std::less< int >,std::allocator< int > >
Returns: std::set<int> All indicies that are excluded

get_num_frz_frz_exclusions()

Get the number of frz-frz exclusions

Return type: int
Returns: int The total number of exclusions in the system

create_frz_exclusions_from_bonds(bonds, bond_cutoff)

Create frozen-frozen exclusions based on each site’s bonds.

Parameters

bonds (std::vector< std::pair< int,int >,std::allocator< std::pair< int,int > > >) – An vector of pair of indicies for each bond
bond_cutoff (int) – All indicies that are this number of bonds away will be excluded

set_external_field(field_x, field_y, field_z)

Apply a constant electric field to the entire system

Parameters

field_x (float) – The x-direction Cartesian component of the field
field_y (float) – The y-direction Cartesian component of the field
field_z (float) – The z-direction Cartesian component of the field

apply_field_to_system(coords)

Apply the field to the system. This must be called before any minimization of the total energy is performed.

Parameters: coords (vec_d) – 1D array of positions of the all sites

get_dipoles(coords)

Return the electric dipoles of the system. This includes nuclear, frozen electrons, dynamic electrons, and total dipoles

Parameters: coords (vec_d) – 1D array of positions of the all sites
Return type: std::vector< Vec3,std::allocator< Vec3 > >
Returns: std::vector<Vec3> The dipoles that are computed. Indicies are determind by DensityType enumeration

get_dipole(coords, density_dype)

Compute the dipole for a particule type of charge

Parameters

coords (vec_d) – 1D array of positions of the all sites
density_dype (int) – DensityType: The type of charges to be used in the computation

Return type

vec_d

Returns

The dipole components

calc_frz_ext_field_energy(coords, forces)

Compute the energy and forces of the frozen charges (nuclei and frozen electrons) with the external electric field

Parameters

coords (vec_d) – 1D array of positions of the all sites
forces (std::vector< Vec3,std::allocator< Vec3 > >) – Array of forces to be updated

Return type

float

Returns

The value of the interaction energy

set_use_PBC(*args)

Overload 1:

Turn on or off periodic boundary conditions

Parameters: is_periodic (boolean) – If true, use PBC, if false, do not use them

Overload 2:

Set the periodic boundaries

Parameters

is_periodic (boolean) – If true, use PBC, if false, do not use them
x (float) – The x-direction of the periodic box
y (float) – The y-direction of the periodic box
z (float) – The z-direction of the periodic box

get_use_PBC()

Inquire if PBC is used or not

Return type: boolean
Returns: true PBC is being used
Return type: boolean
Returns: false PBC is not being used

set_use_cutoff(useCutoff)

Set whether or not to use a nonbonded cutoff. Care must be taked to ensure that the polarization J-matrix remains is positive definite.

Parameters: useCutoff (boolean) – If true, enable the use of a cutoff distance between all interactions. If false, do not use any cutoff distance.

get_use_cutoff()

Inquire of nonbonded cutoffs are being used

Return type: boolean
Returns: true cutoffs are being used
Return type: boolean
Returns: false cutoffs are not being used

set_cutoff_distance(distance_in_nm)

Set the cutoff distance (in bohr)

Parameters: distance_in_nm (float) – The cutoff distance (in bohr)

get_cutoff_distance()

Get the cutoff distance (in bohr)

Return type: float
Returns: distance_in_nm The cutoff distance (in bohr)

set_use_SR_cutoff(use_cutoff)

Set the use SR cutoff. If enabled, electron overlap will be ignored beyond a particular distance.

Parameters: use_cutoff (boolean) – If true, enable the use this apprximation

get_use_SR_cutoff()

Inquire of SR cutoffs are being used

Return type: boolean
Returns: true SR cutoffs are being used
Return type: boolean
Returns: false SR cutoffs are not being used

calc_overlap(coords)

Calculate the total electron overlap

Parameters: coords (vec_d) – the positions of all sites
Return type: float
Returns: the total overlap in e^2 /bohr^3

calc_energy(coords, calc_frz=True, calc_pol=True)

Calculate the energy and forces of the system.

Parameters

coords (vec_d) – A 1D array of all the positions (in bohr)
calc_frz (boolean) – Whether or not to compute frozen electrostatic energies
calc_pol (boolean) – Whether or not to compute polarization

Return type

float

Returns

The total frozen energy computed. Must cal solve_minimization() in order to get polarization energy.

calc_one_electro(deltaR, i, j, calc_pol, calc_frz, energies, forces, thread_num=0)

Compute one pair of interactions between two atomic sites

Parameters

deltaR (DeltaR) – The distance object between the two sites
i (int) – The index of the first site
j (int) – The index of the second site
calc_pol (boolean) – Whether or not to include polarization
calc_frz (boolean) – Whether or not to include frozen electrostatic energies
energies (Energies) – Returns the energies computed
forces (std::vector< Vec3,std::allocator< Vec3 > >) – Returns the forces between both sites
thread_num (int) – the index of the OpenMP thread to compute the energies and forces on

calc_one_frozen(coords, i, j)

Compute one pair of frozen interactions between two atomic sites (does not include polarization)

Parameters

coords (vec_d) – 1D array of all the positions (in bohr)
i (int) – The index of the first site
j (int) – The index of the second site

Return type

Energies

Returns

Energies energies of the computed pair

initialize_calculation(): initizlize arrays before an energy computation is performed. This is called automatically by calc_energy(), but must be called manually to clear all energy and force arrays if calling calc_one_electro() or calc_one_frozen().

solve_minimization(forces)

Minimize the total energy of the system w.r.t. the fluctuating densities. This function solve the matrix equations that minimizes the quadratic polarization energy, under the constraint that all fragment dynamic densities sum to zero.

Parameters: forces (std::vector< Vec3,std::allocator< Vec3 > >) – array of forces to be updated based on the polarization interactions

get_frozen_energy()

Get the energy associated with frozen-frozen dsitity interactions

Return type: float
Returns: double frozen-frozen energy

get_polarization_energy()

Get the polarization energy

Return type: float
Returns: double polarization energy

get_ct_energy()

Get the estimated CT energy as a multiple of the polarization energy

Return type: float
Returns: double

get_energies()

Return all energy components computed

Return type: FlucDensEnergies
Returns: FlucDensEnergies Object that constains the energy components

get_forces()

Get the forces exerted on all particles after an energy computation is performed

Return type: std::vector< vec_d,std::allocator< vec_d > >
Returns: std::vector<vec_d> Nx3 Array of forces on all sites

frz_frz_overlap(inv_r, a, b, exp_ar, exp_br)

Overlap between frozen-frozen dsnsity

Parameters

inv_r (float) – inverse distance between the two sites
a (float) – frozen exponent of the first site
b (float) – frozen exponent of the second site
exp_ar (float) – pre-computed exp(-a*r) using the first exponent
exp_br (float) – pre-computed exp(-a*r) using the secpnd exponent

Return type

float

Returns

double the overlap in e^2/bohr^3

elec_elec_energy(inv_r, a, b, exp_ar, exp_br, dEdR)

Compute the density-density interaction, either frozen or dynamic

Parameters

inv_r (float) – inverse distance between the two sites
a (float) – frozen exponent of the first site
b (float) – frozen exponent of the second site
exp_ar (float) – pre-computed exp(-a*r) using the first exponent
exp_br (float) – pre-computed exp(-a*r) using the secpnd exponent
dEdR (float) – The derivative of the energy w.r.t. the distance

Return type

float

Returns

The elec-elec energy in a.u.

elec_nuclei_energy(inv_r, a, exp_ar, dEdR)

Compute the energy between a nuclei and an electron density

Parameters

inv_r (float) – inverse distance between the two sites
a (float) – frozen exponent of the density site
exp_ar (float) – pre-computed exp(-a*r) using the density exponent
dEdR (float) – The derivative of the energy w.r.t. the distance

Return type

float

Returns

The nuclei-density energy in q.u.

get_rho_coulomb_mat()

Get the delta-rho Coumomb J-matrix

Return type: vec_d
Returns: 1D array of the symmetric matrix elements. Contains NxN elements

get_rho_pot_vec()

Get the electroc potentia vector (linear term)

Return type: vec_d
Returns: vec_d Nx1 array

get_delta_rho()

Get the delta rhos of each site that are computed after minimizing the energy

Return type: vec_d
Returns: vec_d Nx1 array

get_total_time()

Not implimented yet! Return the total wall time used for computations

Parameters: calculate_forces –

set_calc_forces(calculate_forces)

Set wether or not forces should also be computed. If set to false, some computation time can be saved

Parameters: calculate_forces (boolean) – If true, then forces will be computed

print_params(message, param_name)

Print site parameters to the console.

Parameters

message (string) – string: A comment message to be printed
param_name (string) – string: the name of the parameter to be printed

get_param_names()

Get the name of the parameters that can be return

Return type: std::vector< std::string,std::allocator< std::string > >
Returns: std::vector<std::string> A list of possible parameters to enquire about

get_params_by_name(param_name)

Get an array of parameters by it’s name

Parameters: param_name (string) – string: the name of the parameter to return
Return type: vec_d
Returns: vec_d The value of the parameter for each site

static dot3Vec(coords, i, j)

Perform the dot product between two coordinates

Parameters

coords (vec_d) – A 3N array of coordinates for ALL sites
i (int) – index of the first first site
j (int) – index of the second site

Return type

float

Returns

the dot product

property A_mat_save: The quadratic A-matrix that was used in the last energy minimization. This includes The constraint matrix in it’s off-diagonal elements

property B_vec_save: The linear B-vector used in the last energy minimization.