QMCFCMD is a program for combined Quantum Mechanics/Molecular Mechanics calculations using Molecular Dynamics by sampling in an NVT or in a NVE ensemble including a charge field embedding. The default mode is to run only Molecular Mechanics Molecular Dynamics simulation. Upon user request, however, additional ab initio forces within a user defined sphere –the hotspot cell– may be calculated utilizing an external program such as. The total force is then calculated by the equation
VMD screenshot of a QMCFC simulation box. 499 water molecules surrounding a Al(III) ion. The Hotspot area is highlighted as a sphere.
where F_{MM, system} are the MD forces of the simulation system. F_{core} and F_{layer}$ are the quantum mechanical forces of the hotspot cell, including for the layer part the noncoulombic interaction with the water sourrounding . Finaly are the MD forces concerning the solvent sourrounding the spherical QMmechanicaly treated hotspot zone. The smoothing function is mandatory to ensure a continous change of forces between the QM area (the hotspot zone) and the mm solvent part oft the simulation box.
Especially here it smoothes the Forces between the layer and the MM part of the simulation box.
The input file is read either from standard input or by giving the input file name as the first argument to the program. By default the output, energy, trajectory, and velocity files (present from a previous run) are not overwritten. If called with h the program shows a short help message and terminates immediately.
The input file is in strict format keyword = value for the keyword; Note that comment lines starting with a # are allowed, empty lines too. A complete list of available keywords is given in the keyword list
Parameters describing the potential function for the MM part are loaded from the data file guff.dat The potential function for each atom pair follows the general form.
Potential for coulombic forces
Potential for noncoulombic forces
r denotes the distance between pairs of atoms; c,0 − 22 are the special set of parameters for the potential of this two atoms involved. Parameters are usually taken from already published potential fit of older QM/MM simulation. The data for the parametrisation and the literature citations are given in this coefficients tabular. Missing pair interactions were paremtrisated by gaussion or turbomol at HF level. For every molecule, or atom, used in the simulation there has to be a set of parameters for interaction, if not the simulation will not start. Terms which are not used in the potential function are “switched off by setting the appropriate parameters to zero. These sets of parameters for pair interaction are stored in the file guff.dat.
Additionally a threebody potential might be specified to augment the twobody potential. This is necessary especially for doubly and higher charged ions because the pure twobody potentials are too less repulsive. Since a threebody potential is fully repulsive a function like the following proofs to be very useful
where r_{12} is the distances between ion and atom 1, r_{13} is the distances between ion and atom 2, r_{23} is the distance between atom 1 and atomoxygen ; A,B, and C are potential parameters and D,E, and k are optional potential parameters. The last two terms make sure that not only the function but also it’s first derivative are zero at the cutoff limit r_{limit}. Currently the forces are evaluated only for the ion and the two oxygens, hydrogens stay unaffected. Parameters for the 3bd correction terms are listed in the file threebody.dat.
Potential function
From TCWiki
In order to perform a qmcfc simulation potential funktions for each pair interaction in the mm part is needed. To generate those potential functions of the general form:
Potential for coulombic forces
Potential for noncoulombic forces
a reasonable straight forward methodology is supported on our tci desktops.
Contents[hide] 
if (window.showTocToggle) { var tocShowText = “show”; var tocHideText = “hide”; showTocToggle(); }
[edit] Background
In performing QMCF MD or MMMD simulations, the simulation box has areas where interactions are approximated with ab initio constructed twobody and threebody potentials. These potentials reflect the noncoloumbic and coulombic interactions present in the bulk of the solvation box. In many particle systems that have an inherent complexity (e.g. interacting molecules of water exhibit hydrogen bonding), the nature of these potentials would be difficult to reproduce, and it is thus necessary to provide an approximation through numerical means. For solvated metal ion systems, this may be accomplished by considering the variation of interaction energy between a metal ion and a water molecule (for constructing twobody potentials, and a metal ion and two water molecules (for threebody interactions) at varying spatial orientations in cartesian space.
In more practical terms, the process is done in the following manner:
1) Compute for the SCF (self consistent field) single point energies of free ion, E_{ion}, and free water molecule, at a suitable computational level (normally, HartreeFock).
2) Scan for the energy of the ion with the water molecule(s), E_{sys}, at different configurations in cartesian space
This energy is computed at each point in a parameter space defining the positions of each chemical species (ion and water molecule) with each other. Thus for considering twobody interactions depictable from the figure below:
Thus, the hypersurface is specified through the set {r, θ, φ}, for twobody scans, and the set {r1, r2, τ, φ}, for threebody scans.
3) Once E_{sys} is computed for each point in the hypersurface, the interaction energy, E_{inter}, is obtained accordingly as:
n = 1 for twobody interactions is then applied to obtain an analytical approximation to the behavior of the potential observed in the hypersurface. The analytical fit represents the noncoloumbic interaction augmented to the coloumbic attraction present in the system.
4) A fitting algorithm The procedure can be easily extended for considering threebody interactions:
[edit] Computation of energies: E_{ion}, and E_{sys}
SCF (self consistent field) calculations on free ion, free water, and interaction energies of the ion with the water molecule(s) can be done with any standard abinitio software package like Turbomole using appropriate basis sets for the ion and for water. Energy values are often expressed in kcal/mole as . Because basis sets are used in SCF methods to approximate the BornOppenheimer limit, the choice of basis sets requires careful consideration so as to arrive at reliable results. Normally, basis sets with appropriate ECPs (effective core potentials) are used for metal ions, and full valence doublezeta basis set with polarization for the water molecule. In performing a scan, not all points are significant and computational time may be minimized by giving fine sampling only to a region where a global minimum is expected (the interaction for, say, a metal ion and water is at an optimum). The scan through the energy hypersurface may therefore be constructed as a finite nonuniform grid of {r, θ, φ} or of {r1, r2, τ, φ} values with most of the sampling points lying in the region of the expected minimum (see Figure below). This is analogous to riding a car over a terrain and speeding up on flat ground while slowing down over dips and bumps.
[edit] Approximation of noncoloumbic interaction through numerical fitting
Once the energy hypersurface is generated, application of a nonlinear least squares fitting algorithm (e.g. LevenbergMarquardt) over the hypersurface is done to give an analytical form to the behavior of noncoloumbic interaction. Nonlinear least squares fitting algorithms minimize the residual between empirical dataset and a model equation having the set of parameters that are iteratively modified. Normally, a polynomial approximation of the form:
is sufficient in describing the noncoloumbic twobody potentials. The exercise is then finding the most suitable combination of coefficients (A, B, C, …) and exponents (a, b, c,…) that would give the lowest enumeration statistic, χ2.
A similar principle applies for the consideration of threebody interaction except that the model equation does not directly represent a potential function, but rather, correction terms to the twobody potentials. The three body correction terms take the form:
The exact forms of the model equations for twobody and three body potentials may be seen in the documentation of the QMCF program.
The success of iterative methods is intimately tied to convergence; least squares fitting algorithms are no exception. Several factors must be considered in ensuring convergence: the starting values for the parameters to the model equation, the termination criteria (how much fluctuation in the residuals can be tolerated), and a cutoff (a value above which fitting is not applied). The parameters to the model equation may all be minimized simultaneously, or if knowledge on the values of all or any given set is known, it is possible to simply proceed in fixing some parameters while leaving the rest to get modified in the iterations. A subsequent inspection of enumeration statistics would then provide a cursory means for choosing a workable analytical approximation from the fit.
[edit] Implementation
The manner in which the classical potentials are constructed leads to a natural development of tools to implement each of the steps in the foregoing discussion. Following is a listing of these tools and how they work.
[edit] Twobody energy scan
The script tm_2bd_scan generates a dataset of energy values for each of the points in the hypersurface defined by {r, θ, φ} via successive single point energy calculations using Turbomole. For fitting these calculated values of the 2body hypersurface the use of the fit2bd programm is recommended
[edit] Threebody energy scan
The script tm_3bd_scan extends the principle involved two body interaction by another molecule. These additional interaction provides a correction value for the interactions based only on two body potentials. Again the values of the scan can be fitted by the fit3bd programm installed on the TCIDesktop PCs.
[edit] Successes and Pitfalls
TrajToMov is a tool to convert old QMMMMD trajectories to the standardized xyz format (mov).
TrajToMov trajectory_file.traj
Tm 2bd scan
From TCWiki
The tm_2bd_scan is a program for performing a ab initio scan for the creation of a two body pair potential used in many mm simulation approaches (e.g. potential function) The pearl script tm_2bd_scan generates a dataset of energy values for each of the points in the hypersurface defined by {r, θ, φ} via successive single point energy calculations using Turbomole.
syntax ./tm_2bd_scan tm_2bd_scan.in
A file named tm_2bd_scan.crd, is also read in for the specification of the internal geometry of water necessary for the single point energy calculation for the system. The user must also provide a tm_define.in file for the ab initio calculation done from the pearl script.
The program generates two output files: tm_2bd_scan.log, containing a log of the scan, and fit2bd.in, containing the list of obtained energy values versus each configuration in parameter space {r, θ, φ}. It is the latter file that becomes input to the numerical fitting procedure implemented in the program, fit2bd.
Fit2bd
From TCWiki
The programm fit2bd provides to fit an analytical expression on the datapoints calculated by the tm_2bd_scan script.
[edit] Fitting procedure for the two body scan
The program, fit2bd, applied to the dataset in fit2bd.in (a data file from a run of the tm_2bd_scan, generates an output file containing values for the coefficients A, B, C, and D of the analytical approximation. It is important to note that the program does not simultaneously minimize coefficients and exponents of the model approximation, but rather proceeds by minimizing the set of coefficients given a fixed set of exponents. This is accomplished via compile time modification to the source code. To be able to perform minimization over a defined set of exponents, one can either manually modify the source code or employ a script to perform the modification, with a subsequent rebuild and execution of the program.
Because iterative methods of approximation are, at best, a directed search, prior knowledge or experience on classical potentials that can be assigned to the system under study is necessary to avoid naively testing all combinations of fixed parameters to the model equation. This, and the inherent effect of the energy cutoff (upperbound of energy values to consider for the fit) defined in fit2bd.in to the convergence of the fitting procedure, must be kept well in mind for any useful result to be had.
Currently, the program is only applicable to twobody potentials involving ionwater interaction The output (depending on your changes to the fitted formula) contains the coefficients, for the pair interactions.
Syntax: ./fit2bd <fit2bd.in > output_file.dat
Tm 3bd scan
From TCWiki
The script, tm_3bd_scan, extends the principle involved in the twobody scan to the parameter space {r1, r2, τ, φ}. However, because of the increase in the number of input parameters, the format of the input file and the behavior of the script is less direct than that in the case of twobody scans. Several input files must be prepared, each corresponding to a value of τ and a range of values of φ. The program places its output to a fit3bd.in file corresponding to the parameters specified in the input file. This must however be stored and renamed accordingly for later consolidation into one fit3bd.in file that would contain the full hypersurface as input to the fitting procedure. As with the twobody scan, a log file can also be generated for inspection on the success/failure of the scan.
Syntax: "./tm_3bd_scan < input_file.in > log_file.log"
Fit3bd
From TCWiki
The program fit3bd behaves differently from fit2bd, as no compile time modifications must be made prior to execution. The important distinction is that because the threebody fitting procedure is for construction of correction terms to the two body potentials, the energy cutoff is not treated as an upperbound but a lowerbound of energies. This also means that the generated coefficients, A1A3, implicitly describe threebody interaction involving water.
As with fit2bd, the fit3bd program generates a list of values for the parameters to the model equation that each have a corresponding enumeration statistics.
Syntax: "./fit3bd < input_file.in > output_file.log"
Coulomb’s law
From Wikipedia, the free encyclopedia
Coulomb’s law, sometimes called the Coulomb law, is an equation describing the electrostatic force between electric charges. It was studied and first published in the 1780s by French physicist Charles Augustin de Coulomb and was essential to the development of the theory of electromagnetism. Nevertheless, the dependence of the electric force with distance (inverse square law) had been proposed previously by Joseph Priestley^{[1]} and the dependence with both distance and charge had been discovered, but not published, by Henry Cavendish, prior to Coulomb’s works.
Coulomb’s law may be stated in scalar form as follows:
 The magnitude of the electrostatic force between two point electric charges is directly proportional to the product of the magnitudes of each of the charges and inversely proportional to the square of the total distance between the two charges.
Contents[hide] 
//<![CDATA[
if (window.showTocToggle) { var tocShowText = “show”; var tocHideText = “hide”; showTocToggle(); }
//]]>
[edit] Scalar form
The scalar form of Coulomb’s law will only describe the magnitude of the electrostatic force between two electric charges. If direction is required, then the vector form is required as well. The magnitude of the electrostatic force (F) on a charge (q_{1}) due to the presence of a second charge (q_{2}), is given by
where r is the distance between the two charges and k_{e} a proportionality constant. A positive force implies a repulsive interaction, while a negative force implies an attractive interaction.^{[2]}
The proportionality constant k_{e}, called Coulomb’s constant (sometimes called Coulomb’s force constant) is related to the properties of space and can be calculated exactly:^{[3]}
In SI units the speed of light in vacuum, denoted c_{0}^{[4]} is defined as 299,792,458 m·s^{−1},^{[5]} and the magnetic constant (μ_{0}), is defined as 4π × 10^{−7} H·m^{−1},^{[6]} leading to the definition for the electric constant (ε_{0}) as ε_{0} = 1/(μ_{0}c20) ≈ 8.854187817×10^{−12} F·m^{−1}.^{[7]} In cgs units, the unit charge, esu of charge or statcoulomb, is defined so that this Coulomb force constant is 1.
This formula says that the magnitude of the force is directly proportional to the magnitude of the charges of each object and inversely proportional to the square of the distance between them. The exponent in Coulomb’s Law has been found to differ from −2 by less than one in a billion.^{[8]}
When measured in units that people commonly use (such as SI—see International System of Units), the electrostatic force constant (k_{e}) is numerically much much larger than the universal gravitational constant (G).^{[4]} This means that for objects with charge that is of the order of a unit charge (C) and mass of the order of a unit mass (kg), the electrostatic forces will be so much larger than the gravitational forces that the latter force can be ignored. This is not the case when Planck units are used and both charge and mass are of the order of the unit charge and unit mass. However, charged elementary particles have mass that is far less than the Planck mass while their charge is about the Planck charge so that, again, gravitational forces can be ignored. For example, the electrostatic force between an electron and a proton, which constitute a hydrogen atom, is almost 40 orders of magnitude greater than the gravitational force between them.^{[9]}
Coulomb’s law can also be interpreted in terms of atomic units with the force expressed in Hartrees per Bohr radius, the charge in terms of the elementary charge, and the distances in terms of the Bohr radius.
[edit] Electric field
It follows from the Lorentz Force Law that the magnitude of the electric field (E) created by a single point charge (q) at a certain distance (r) is given by:
For a positive charge, the direction of the electric field points along lines directed radially away from the location of the point charge, while the direction is the opposite for a negative charge. The SI units of electric field are volts per meter or newtons per coulomb.
[edit] Vector form
In order to obtain both the magnitude and direction of the force on a charge, q_{1} at position , experiencing a field due to the presence of another charge, q_{2} at position , the full vector form of Coulomb’s law is required.
where r is the separation of the two charges. Note that this is simply the scalar definition of Coulomb’s law with the direction given by the unit vector, , parallel with the line from charge q_{2} to charge q_{1}.^{[9]}
If both charges have the same sign (like charges) then the product q_{1}q_{2} is positive and the direction of the force on q_{1} is given by ; the charges repel each other. If the charges have opposite signs then the product q_{1}q_{2} is negative and the direction of the force on q_{1} is given by ; the charges attract each other.
[edit] System of discrete charges
The principle of linear superposition may be used to calculate the force on a small test charge, q, due to a system of N discrete charges:
where q_{i} and are the magnitude and position respectively of the i^{th} charge, is a unit vector in the direction of (a vector pointing from charge q_{i} to charge q), and R_{i} is the magnitude of (the separation between charges q_{i} and q).^{[9]}
[edit] Continuous charge distribution
For a charge distribution an integral over the region containing the charge is equivalent to an infinite summation, treating each infinitesimal element of space as a point charge dq.
For a linear charge distribution (a good approximation for charge in a wire) where gives the charge per unit length at position , and is an infinitesimal element of length,
 .^{[10]}
For a surface charge distribution (a good approximation for charge on a plate in a parallel plate capacitor) where gives the charge per unit area at position , and is an infinitesimal element of area,
For a volume charge distribution (such as charge within a bulk metal) where gives the charge per unit volume at position , and is an infinitesimal element of volume,
 ^{[9]}
The force on a small test charge at position is given by
[edit] Graphical representation
Below is a graphical representation of Coulomb’s law, when q_{1}q_{2} > 0. The vector is the force experienced by q_{1}. The vector is the force experienced by q_{2}. Their magnitudes will always be equal. The vector is the displacement vector between two charges (q_{1} and q_{2}).
[edit] Electrostatic approximation
In either formulation, Coulomb’s law is fully accurate only when the objects are stationary, and remains approximately correct only for slow movement. These conditions are collectively known as the electrostatic approximation. When movement takes place, magnetic fields are produced which alter the force on the two objects. The magnetic interaction between moving charges may be thought of as a manifestation of the force from the electrostatic field but with Einstein’s theory of relativity taken into consideration.
[edit] Table of derived quantities
Particle property  Relationship  Field property  
Vector quantity 



Relationship  
Scalar quantity 


[edit] See also
 BiotSavart Law
 Method of image charges
 Electric field
 Electric constant
 Coulomb, the SI unit of electric charge named after Charles Augustin de Coulomb
 Electromagnetic force
 Molecular modelling
[edit] Notes
 ^ Robert S. Elliott (1999), Electromagnetics: History, Theory, and Applications, ISBN 9780780353848, http://eu.wiley.com/WileyCDA/WileyTitle/productCd0780353846.html
 ^ Coulomb’s law, Hyperphysics
 ^ Coulomb’s constant, Hyperphysics
 ^ ^{a} ^{b} Current practice is to use c_{0} to denote the speed of light in vacuum according to ISO 31. In the original Recommendation of 1983, the symbol c was used for this purpose. See NIST Special Publication 330, Appendix 2, p. 45
 ^ http://physics.nist.gov/cuu/Units/meter.html]
 ^ http://physics.nist.gov/cuu/Units/ampere.html]
 ^ http://physics.nist.gov/cgibin/cuu/Value?ep0
 ^ Williams, Faller, Hill (1971), “New Experimental Test of Coulomb’s Law: A Laboratory Upper Limit on the Photon Rest Mass“, Physical Review Letters 26: 721–724, doi:, http://prola.aps.org/abstract/PRL/v26/i12/p721_1
 ^ ^{a} ^{b} ^{c} ^{d} Coulomb’s law, University of Texas
 ^ Charged rods, PhysicsLab.org
[edit] References
 Griffiths, David J. (1998). Introduction to Electrodynamics (3rd ed.). Prentice Hall. ISBN 013805326X.
 Tipler, Paul (2004). Physics for Scientists and Engineers: Electricity, Magnetism, Light, and Elementary Modern Physics (5th ed.). W. H. Freeman. ISBN 0716708108.
[edit] External links
 Coulomb’s Law on Project PHYSNET.
 Electricity and the Atom — a chapter from an online textbook
 A maze game for teaching Coulomb’s Law—a game created by the Molecular Workbench software
 The inverse cube law The inverse cube law for dipoles (PDF file) by Eng. Xavier Borg
 Electric Charges, Polarization, Electric Force, Coulomb’s Law Walter Lewin, 8.02 Electricity and Magnetism, Spring 2002: Lecture 1 (video). MIT OpenCourseWare. License: Creative Commons AttributionNoncommercialShare Alike.