eFF - ProkopHapala/FireCore GitHub Wiki
Quick start
Python scrips
examples of running eFF library over python API are located in tests/tEFF/
- First it is necessary to compile the C++ library
libeFF_lib.so
. This can be done withcmake
inFireCore/cpp/Build/
following the Readme.md.- Note that
libeFF_lib.so
does not depend on OpenCL,OpenGL,nor SDL so both option can be set off (cmake .. -DWITH_SDL=OFF -DWITH_OPENCL=OFF
) and no additional libraries have to be installed.
- Note that
- running
tests/tEFF/run.sh
will automatically recompile theFireCore/cpp/Build/libs/Molecular/libeFF_lib.so
and than runstests/tEFF/run_dynamics.py
simulation of impact of free-flying electron of water molecule. - call function
relax_mol("H2O")
from withing the scripttests/tEFF/run.py
will relax H2O molecule starting from geometry in stored inFireCore/cpp/sketches_SDL/Molecular/data/H2O.fgo
- it produces final relaxed structure in
tests/tEFF/H2O_relaxed.fgo
and relaxation trajectorytests/tEFF/H2O_relax.xyz
.H2O_relax.xyz
can be visualized with jmol- electrons are represented by heavy atoms (92 U for spin +1, and 109 Mt for spin -1)
- Jmol visualizes 92U=spin_up as blue and 109Mt=spin_down as red.
- Electron radius saved as 6th columns in .xyz is also properly visualized in Jmol
- number of atoms and electrons as well as components of total energy are stored in the comment (2nd) line of .xyz
- see function EFF::to_xyz() in cpp/common/molecular/eFF.h for more information about the output .xyz trajectroy
- electrons are represented by heavy atoms (92 U for spin +1, and 109 Mt for spin -1)
- it produces final relaxed structure in
Python API
from pyBall import eFF as eff
loads eFF library. You may need to setup pythonpath first. For example assuming you rate in tests/tEFF/
do sys.path.append("../../"); from pyBall import eFF as eff
To initialize eFF (allocate arrays) you may run eff.init(natoms,nelectrons)
but the arrays will be uninitalized. Typically you rather run
eff.load_fgo("molecule.fgo")
instead to initialize the problem from .fgo file.
Python API allows to share internal C data-arrays by wrapping them as numpy arrays. This allows both read and modify internal state of eFF library on the fly without boiler-plate functions. This can be done by calling eff.eff.getBuffs( )
. Then you can read internal C buffers from python e.g. by calling print(eff.eforce)
. Or you can modify them. E.g. eff.eforce[-1,1]+=50.0
increase velocity of last electron along y-axis by 50.0 [A/fs].
And example how to use python API see /tests/tEFF/run.py.
For example to relax molecular geometry (both atoms and electrons) call:
eff.load_fgo("initial_geometry.fgo" ) # load molecule in .fgo format (i.e. floating-gaussian-orbital)
eff.initOpt(dt=0.03,damping=0.1 ) # initialize optimizer/propagator
eff.run( 10000, Fconv=1e-3, ialg=2 ) # run simuation for maximum 1000 time steps intil it converge to |F|<1e-3, ialg=2 is FIRE http://users.jyu.fi/~pekkosk/resources/pdf/FIRE.pdf https://www.sciencedirect.com/science/article/pii/S0927025620300756
eff.save_fgo(name+"_relaxed.fgo" ) # save final relaxed geometry to .fgo format (i.e. floating-gaussian-orbital).
To run electron dynamics simulation do
eff.load_fgo("data/H2O_shoot.fgo", bVel_=True ) # load H2O moleule with electron-bullet in .fgo format (i.e. floating-gaussian-orbital) including initial velocities
eff.initOpt(dt=0.001,damping=0.0, bMass=True) # initialize optimizer/propagator
eff.invSmass[:] = 0 # fix size of electrons by setting invMass=0 (mass=innfinity)
eff.setTrjName("trj_bullet_vs_H2O.xyz") # setup output .xyz file to save trajectory of all atoms and electrons at each timespep (comment-out to ommit .xyz and improve performance )
eff.run( 1000, ialg=-1 ) # run eFF for 10 iterations, using leap-frog progrgator (ialg=-1)
.fgo file format (floating gaussian orbitals)
See function EFF::loadFromFile_fgo()
in cpp/common/molecular/eFF.h for reference about reading .fgo format.
Typical .fgo format is as follows:
natoms nelectrons 1 0
x y z Q sQ sP cP [vx vy vz] # atom Q=core charge, sQ=core charge smearing radius, sP=core electron radius, cP=Pauli repulsion stregth of core electrons, velocities [vx vy vz] are optional
x y z s 1.0 ispin [vx vy vz vs] # electron with position (x,y,z) and size s and spin ispin={+1 or -1}. velocities [vx vy vz vs] are optional
example H atom with 1 electron spin up (see cpp/sketches_SDL/Molecular/data/H_eFF.fgo)
1 1 1 0
+0.0 0.0 0.0 -1. 0.0 0.0 0.0
+0.0 0.0 0.0 0.5 1.0 +1
Frozen core approximation
Original formulation of eFF is able to simulate all electron (at least for light elements with Z=1..9). But due to strong interaction of inner shell (s1) with nucleus these electrons oscillate rapidly, and require order of magnitude shorter time-step. For example, relaxation of free electron pairs in H2O molecule together with core electron is very painful and requires >10000 iterations.
For this reason and for improved numerical efficiency I implemented frozen core version of eFF which fix both position and radius of the electron in the inner shells. Instead of treating core electrons as independent particles, I replace them by corresponding in (pseudo-)potential of atomic nucleus comprising from both electrostatic and Pauli terms. All core electrons are currently represented by single Gaussian with single width and corresponding charge. Currently I tested only 2nd period elements (Li,Be,B,C,N,O) which have 1 spin pair (2 electrons) in core. Thanks to this modification relaxation of molecules such as H2O,NH3,CH4 is order of magnitude faster and significantly more numerically stable.
Form comparison this is H2O molecule with frozen core (core radius is 0.2A and contains 1.0 of electron pairs, leading to effective charge -6.0):
FireCore/cpp/sketches_SDL/Molecular/data/H2O_fixcore.fgo
3 8 1 0
0.0 0.0 0.0 -8. 0.0 0.2 1.0
+1.0 0.0 1.0 -1. 0.0 0.0 0.0
-1.0 0.0 1.0 -1. 0.0 0.0 0.0
+1.0 0.0 0.5 1.0 1.0 +1
+1.0 0.0 0.5 1.0 1.0 -1
-1.0 0.0 0.5 1.0 1.0 +1
-1.0 0.0 0.5 1.0 1.0 -1
0.0 +1.0 -0.5 1.0 1.0 +1
0.0 +1.0 -0.5 1.0 1.0 -1
0.0 -1.0 -0.5 1.0 1.0 +1
0.0 -1.0 -0.5 1.0 1.0 -1
In contrast this is equivalent H2O molecule without frozen core, which has additional 2 electron with size 0.2A. Both systems should evaluate same total energy, although this energy is partitioned differently between componens Eee, EeePauli, Eae, EaePauli
.
FireCore/cpp/sketches_SDL/Molecular/data/H2O.fgo
3 10 1 0
0.0 0.0 0.0 -8. 0.0 0.0 0.0
+1.0 0.0 1.0 -1. 0.0 0.0 0.0
-1.0 0.0 1.0 -1. 0.0 0.0 0.0
+1.0 0.0 0.5 0.5 1.0 +1
+1.0 0.0 0.5 0.5 1.0 -1
-1.0 0.0 0.5 0.5 1.0 +1
-1.0 0.0 0.5 0.5 1.0 -1
0.0 +1.0 -0.5 0.5 1.0 +1
0.0 +1.0 -0.5 0.5 1.0 -1
0.0 -1.0 -0.5 0.5 1.0 +1
0.0 -1.0 -0.5 0.5 1.0 -1
0.0 0.0 0.0 0.2 1.0 +1
0.0 0.0 0.0 0.2 1.0 -1
Effect of core radius on molecular geometry
With frozen core electrons, the core electron radius is a tuneable parameter which allows to fine-tune equlibrium distance of valence electrons from nuclei which has significant effect on resulting chemical structures of molecules - e.g. on bond lengths. Non-trivial dependence of bond length in CH4,NH3 and H2O si presented in following chart.
As you can see the lenght of C-H bond in CH4 monotoneously and roughly linearly increase with core radius, while in case of H2O molecule there is roughly parabolic minimum for core radius 0.2 [A]. How can we understand this?
- In case of CH4 larger core of C atom repels more the bonding electrons therefore increasing bond lenghs.
- In case of H2O we have 3 different kinds of electrons (1. core, 2. free electron pairs, 3. sigma bonds). The effect of core radius affect considerably more the (2.) free electron pairs than the (3.) bond electrons. This is because bonding electrons are affected by both nuclei, while free electron pairs depend on single nucleus. Small core radius allows free electron pairs to concentrate much closer to nucleus, which in effect repels the bonding electrons further from nucleus. This is the reason for contra-intuitive increase of bond-lengths for smaller core radius (<0.2A). For larger core radius (>0.2A) the situation is similar to to CH4 - bonding electrons are repelled by core electrons causing increase of bondlenghs with core size.
Standalone graphical program
- First you have to compile the C++ subfolder (
FireCore/cpp/
) according to Readme.md. Do not forget activate-DWITH_SDL=ON
since this is graphical application. - Then run
/FireCore/cpp/Build/sketches_SDL/Molecular/test_eFF_MD
. - Hit
space
key to unleash the simulations. - It presents H2O molecule hit by fast electron with speed ~50 [A/fs].
- In real time are plotted components of energy and force acting of electrons and atoms.
Where to find relevant source code files ?
- Actual electron forcefield implementation /cpp/common/molecular/eFF.h
- Integrals between gaussian functions and its derivatives cpp/common/molecular/InteractionsGauss.h
- Dynamic liberty
libeFF_lib.so
which expose functionality of eFF.h as pure C functions cpp/libs/Molecular/eFF_lib.cpp - Python API wraper for
libeFF_lib.so
/pyBall/eFF.py - graphical example programs using eFF in
/FireCore/cpp/Build/sketches_SDL/Molecular/
- relaxation test_eFF.cpp
- semiclasical electron dynamics test_eFF_MD.cpp