
.. _dacapo netcdf: http://oldwww.fysik.dtu.dk/CAMPOS/Documentation/Dacapo/dacapoNETCDF/dacapoNETCDF.html

.. _dacapo pseudopotential page: http://www.fysik.dtu.dk/CAMPOS/Documentation/Dacapo/PseudoPotentialOverView/index.html

.. _dacapo Ascii file: http://www.fysik.dtu.dk/CAMPOS/Documentation/Dacapo/dacapoNETCDF/node61.html

.. _DFT Calculator: http://www.fysik.dtu.dk/campos/ASE/manual/manual.html#dft-calculator-interface

.. _ListOfAtoms: http://www.fysik.dtu.dk/campos/ASE/manual/manual.html#listofatoms

.. _ASE manual: http://www.fysik.dtu.dk/campos/ASE/manual/manual.html

.. _filter: http://www.fysik3.dtu.dk/campos/ASE/manual/manual.html#filter

.. _dynamics: http://www.fysik3.dtu.dk/campos/ASE/manual/manual.html#dynamics

.. _trajectory: http://www.fysik3.dtu.dk/campos/ASE/manual/manual.html#trajectory

.. _Nudged Elastic Band method: http://www.fysik.dtu.dk/campos/ASE/manual/manual.html#nudged-elastic-band-method

.. _h2o example: http://www.fysik.dtu.dk/campos/Dacapo/examples/examples.html#combining-several-constrain-filters

.. _Wannier example: http://www.fysik.dtu.dk/campos/Dacapo/examples/examples.html#c2h4

.. _plot wavefunction: http://www.fysik.dtu.dk/campos/Dacapo/examples/examples.html#plotwavefunction

.. _tutorials: http://www.fysik.dtu.dk/campos/Dacapo/tut/tutorials.html

.. _Nudged Elastic Band example: http://www.fysik.dtu.dk/campos/Dacapo/examples/examples.htm#neb

.. _Harris calculation example: http://www.fysik.dtu.dk/campos/Dacapo/examples/examples.htm#harris


------------
Introduction
------------

Dacapo is a total energy program based on density functional theory. It uses a
plane wave basis for the valence electronic states and describes the
core-electron interactions with Vanderbilt ultrasoft pseudo-potentials. The
program performs self-consistent calculations for both Local Density
Approximation (LDA) and various Generalized Gradient Approximation (GGA)
exchange-correlations potentials, using state-of-art iterative algorithms. The
code may perform molecular dynamics / structural relaxation simultaneous with
solving the Schrodinger equations within density functional theory. The
program may be compiled for seriel as well as parallel execution and the code
has been ported to many hardware platforms.

The manual describes how to use the dacapo program from the Campos 
Atomistic Simulation enviroment, or for short ASE.  

.. _SetListOfAtoms:

-----------------------------------
Setting up the atomic configuration 
-----------------------------------

The atoms for Dacapo are defined via the ASE ``ListOfAtoms``
object. 

A ``ListOfAtoms`` object is a collection of atoms.  The ASE
implementation of a ListOfAtoms can be used like this:

>>> from ASE import Atom, ListOfAtoms
>>> atoms = ListOfAtoms([Atom('C', (0, 0, 0)),
...                  Atom('O', (0, 0, 1.0))])
>>> atoms.SetUnitCell((10,10,10))

This will define a CO molecule, with a CO distance of 1A, 
in a cubic 10Ax10Ax10A unitcell.  

For more details of the ``ListOfAtoms`` class see
`listofatoms`_  in the ASE manual.

The ``atoms`` defined here can be attached to a ``Dacapo`` calculator 
by using: 

   atoms.SetCalculator(calculator)

alternatively, one can attach the ``atoms`` to a calculator by using: 
 
   calculator.SetListOfAtoms(atoms) 
  
The next chapter describes how to define the calculator for dacapo. 

-------------------------------
Dacapo ASE calculator interface
-------------------------------

As a minimum for defining a calculation in python you need to specify the atoms 
as described above, and you need to set the planewave cutoff and the number of
bands,  this will look like this using the python interface: 

>>> from Dacapo import Dacapo 
>>> calculator = Dacapo()
>>> calculator.SetNumberOfBands(6)
>>> calculator.SetPlaneWaveCutoff(340)
>>> atoms.SetCalculator(calculator)
>>> calculator.GetPotentialEnergy() 

The rest of this chapter describes the different methods you 
can use for the ``calculator`` object.
More detailed explanation of the python construction found above can be 
found in the `tutorials`_.

For the most important of the method described one can use ``keyword`` 
arguments in the construction of the ``Dacapo`` calculator, 
the table below list the ``keywords``:

============================   =================== 
 Method                        keyword
============================   ===================
`SetPlaneWaveCutoff`_           ``planewavecutoff``
`SetDensityCutoff`_             ``densitycutoff``
`SetNumberOfBands`_            ``nbands``
`SetXCFunctional`_             ``xc``
`SetBZKPoints`_                ``kpts``
`SetSpinPolarized`_            ``spin``
`SetNetCDFFile`_               ``out``
`SetTxtFile`_                  ``txtout``
`DipoleCorrection`_            ``dipole``
`SetScratch`_		       ``scratchdir``
============================   ===================


So the example above can be written more compact: 

>>> calculator = Dacapo(
...                     nbands=6,
...                     planewavecutoff=340) 


.. _SetNumberOfBands: 

Number of bands and Electronic temperature 
------------------------------------------

The method '``SetNumberOfBands(nbands)``' can be used to 
define the number of bands for the calculatior. ``nbands`` should be a least 
equal to halv the number of valence electrons.

The occupation statitstics can be set using the method 
``SetOccupationStatistics(method)``.

dacapo support two methods for occupation statistics: 

``FermiDirac``: 
  Fermi-Dirac occupation statistics. The temperature is controlled 
  by the method ``SetElectronicTemperature(temp)``. Default is 
  0.1 eV. 

``MethfesselPaxton``: 
  Currently, the Methfessel-Paxton scheme (PRB 40,3616, 1989) 
  is implemented to 1st order. 


.. _SetPlaneWaveCutoff: 

.. _SetDensityCutoff:

Planewave and Density-cutoff
----------------------------

The method ``SetPlaneWaveCutoff(planewavecutoff)`` defines the plavewave cutoff
for the calculation. The unit for ``planewavecutoff`` is eV.
 
The method ``SetDensityCutoff(densitycutoff)`` sets the density cutoff. 
By default the two cutoffs are equal, corresponding to using just one grid. 
For elements like ``Cu`` using ultra-soft pseudo potentials the ultra-soft 
pseudo-wavefunction are very soft, while the density still requires a 
fine grid, in this case it can be advantagous to set a larger value for 
``densitycutoff`` that for ``planewavecutoff``.

.. _SetNetCDFFile: 

.. _SetTxtFile: 

.. _SetScratch:

Changing file names
-------------------

The following two methods can be use to change the file names
used by dacapo: 

``SetNetCDFFile(filename)``:
  This will set the NetCDF filename, default is 'out.nc'.

``SetTxtFile(filename)`` 
  This will set the Ascii text file name, default is 'out.txt'. 
  For information on the most inportant data in the Ascii file 
  see `dacapo Ascii file`_.

The method:

``SetScratch(scratchdir)``:

will set the scratch directory used for swapfiles to ``scratchdir``. 
Default for this directory is ``/scratch/<username>`` or if this 
is not found just ``/scratch``. 


    

.. _SetXCFunctional: 

Exchange-correlation functional
-------------------------------

You can use the method 
``SetXCFunctional(exc)``
to select the exchange-correlation used for the calculation. 
``exc`` is one of: 

``VWN``:
  Vosko Wilk Nusair LDA-parametrization.

``PW91``: 
  Perdew Wang 91 GGA-parametrization.

``PBE``:
  Perdew Burke Ernzerhof GGA-parametrization.


``RPBE``:
  Revised PBE GGA-parametrization, Hammer, Hansen and Nrskov.


The method 
``GetXCFunctional``
tells which of the above defined exchange-correlation functionals 
are currently  used. 

.. _SetSpinPolarized: 

Defining a spin-polarized calculation 
-------------------------------------

One can use the method 
``SetSpinPolarized(True)``
to tell dacapo that the calculation should be spin-polarized. 

One also need to specify an initial magnetic moment for the 
atoms. You can do this in the construction of the atoms, 
by using: 

>>> from ASE import Atom, ListOfAtoms
>>> atoms = ListOfAtoms([Atom('O', (0, 0, 0),magmom=1.0),
...                      Atom('O', (0, 0, 1.2),magmom=1.0)])

This can also be set using the ListOfAtoms method 
``SetMagneticMoments((1,1))`` or setting the ``SetMagneticMoment`` method on each atom. 

This will define a initial total magnetic moment of 2 Bohr magneton, 
for the O2 molecule, this being also the self-consistent result. 
You could also use, SetMagneticMoments((1,-1)), this will define a
'anti-ferro' magnetic state with total magnetic moment of zero. 
If this case probably the end result is still the same (2 Bohr 
magneton) but the number of self-consistent iterations will increase. 
 

.. _SetBZKpoints: 

K-points
--------

The method ``SetBZKPoints(kpts)`` can be used to define the k-points in the BZ. 
In general ``kpts`` are a python list of k-points defined in units of the 
reciprocal cell. You can however specify a 10x10x10 Monkhorst-Pack 3dimensional set 
simple by using  ``SetBZKpoint((10,10,10))``. 

Below is a list of the special k-points sets that can be used: 


``MonkhorstPack``:
  A 3-dimensional Monkhorst-Pack special k-point set.  
  (see H.J.Monkhorst and J.D.Pack, Physical Review B, vol. 13, page 5188, 1976. 
  You can defined this k-points set simply using, SetBZKPoints((n1,n2,n3)). This will define a n1 by n2 by n3,
  Monkhorst Pack special k-point set.  
  The size along each reciprocal vector given by  ``n1``, ``n2`` and ``n3``, respectively.  


``ChadiCohen``:
  A Chadi-Cohen (see Chadi and Cohen, Phys. Review B., vol 8,page 5747) 1-dimensional k-point can be defined 
  using the ASE Utilities module ChadiCohen: 
 
  >>> from ASE.Utilities.ChadiCohen import CC6_1x1
  >>> calc.SetBZKPoints(CC6_1x1)

  This will define a 6 k-point 1x1 Chadi-Cohen special k-point set.
  In general the Chadi-Cohen sets are named CC+'<Nkpoints>'+_+'shape'. 
  For example an 18 k point sq(3)xsq(3) is
  named 'CC18_sq3xsq3'.


You can use the method ``GetBZKPoints``, to access the k-points just defined. 
The method ``GetIBZKPoints`` returns the symmetry reduced set of k-points. 


Pseudo-potentials
-----------------

You can change or view the pseudopotentials use for an element using the 
methods: 

``GetPseudoPotential(Z=<elementnumber>)``: 
  View the pseudopotential path for element ``elementnumber``.  

``SetPseudoPotential(Z=<elementnumber>,path=path)``: 
  Set the pseudopotential path for element ``elementnumber``.

Charge-density mixing
---------------------

The method ``SetChargeMixing(value)`` can be used to enable or disable 
charge density mixing. 

Using ``SetChargeMixing(True)`` the program will use Pulay mixing for the density
(default). 

Using ``SetChargeMixing(False)`` correspond to running 
a calculation using the Harris-Foulkes functional using 
the input density. See also the `Harris calculation example`_. 


The method ``SetKerkerPreconditioning(value)`` can be used to set Kerker preconditioning 
for the density mixing.  
For value=True Kerker preconditiong is used, i.e. q is different from zero, see eq. 82 in 
Kresse/Furthmller: Comp. Mat. Sci. 6 (1996). The value of q is fix to give a damping of 20 of 
the lowest q vector.
For value=False q is zero and mixing is linear (default).


Eigenvalue solver
-----------------  

The eigenvalue solver can be set using the method ``SetEigenvalueSolver(method)``.
Here ``method`` is one of: 

A. ``eigsolve``: 
   Block Davidson algorithm (Claus Bendtsen et al). This is the default method. 

B. ``rmm-diis``: 
   Residual minimization method (RMM), using DIIS (direct inversion in the iterate subspace). 
   The implementation follows closely the algorithm outlined in Kresse and 
   Furthmuller, Comp. Mat. Sci, III.G/III.H


.. _DipoleCorrection: 

Dipole-correction
-----------------

``DipoleCorrection(MixingParameter=0.2,InitialValue=0.0)``: 
  This methods helps setting up the Dipole Correction scheme, 
  For more information see the `dacapo netcdf`_ manual. 


Setting non-default executable 
------------------------------

``GetJobType``:
  SetJobType can be used to view the current executable for 
  dacapo.

``SetJobType(executable=newexecutable)``:
  SetJobType can be used to change the default dacapo executable 
  ``dacapo.run`` to ``newexecutable``. 

Symmetry
--------

``SetSymmetryOff``:
  Do not use symmtry then reducing the k-point set. 

``SetSymmetryOn``:
  Use symmtry then reducing the k-point set. 
 


Electrostatic decoupling 
------------------------

``SetElectrostaticDecoupling(numberofgaussians=3,ecutoff=100)``:
  Add electrostatic decoupling (Jan Rossmeisl). 
  Decoupling activates the three dimensional electrostatic decoupling.  
  Based on paper by Peter E. Bloechl: JCP 103 page7422 (1995).


Stress calculation
------------------

By default dacapo will not calculate the stress acting on the unitcell, 
you can use the method:
 
``CalculateStress``:
  
to tell dacapo to calculate the stress. 


Local density of states
-----------------------

Dacapo can calculate the density of states projected onto atomic
orbitals.

.. warning:: 
   The density of states are not implemented for planewave parallel 
   calculations, 
   so the number of nodes must match the number kpoints (after 
   they have been reduce by symmetry). 
   So if you for example have 4 IBZ kpoints,  you can use 
   4 or 2 parallel nodes and still get the projected density of 
   states. 



Typically one would do: 


>>> atoms = ListOfAtoms(..) 
>>>
>>> calc = Dacapo(..)
>>> calc.CalculateAtomicDOS() 
>>>
>>> atoms.SetCalculator(calc)
>>>
>>> energy = atoms.GetPotentialEnergy()
>>> dos = calc.GetLDOS(atoms=[1],angularchannels=["s"]) 
>>> dos.GetPlot()

This will return a density of states projected onto the s-orbital 
of the first atom (atom numbers starts from 1). 
  
The method:  

  `CalculateAtomicDOS(energywindow=None)` 

tells the fortran program to calculate the  local density of states 
projected onto atomic orbitals, 
using ``energywindow`` for calculating the energy resolved 
LDOS.   

The method: 

   ``GetLDOS(**keywords )``:  

returns an instance of ``AtomProjectedDOSTool`` (see below). 

Keyword for the ``GetDOS`` method:

  * atoms
           
    List of atoms numbers (default is all atoms)

  * angularchannels

    List of angular channel names.
    The full list of names are: 
    's', 'p_x', 'p_y', 'p_z','d_zz','dxx-yy', 'd_xy', 'd_xz', 'd_yz'.                
    'p' and 'd' can be used as shorthand for all p and all d channels
    respectively. (default is all d channels)

  * spin
            
    List of spins (default all spins)

  * cutoffradius
    
    'short' or 'infinite'.
    For cutoffradius = 'short' the integrals are truncated at 1 Angstrom.
    For cutoffradius = 'infinite' the integrals are not truncated.
    (default 'infinite')

The resulting DOS are added over the members of the three list (atoms,angularchannels and spin).

GetDOS returns a instance of the class AtomProjectedDOSTool (Thanks to John Kitchen for providing 
grace plot and band moments methods). 

Methods for AtomProjectedDOSTool:

   ``GetPlot``       
    returns a GnuPlotAvartar for the energy resolved DOS
    A parent can be given too combine plot.

   ``GetIntegratedDOS`` 
    returns the integral up to the fermi energy

   ``GetData`` 
    returns the energy resolved DOS

   ``GetBandMoment(moments)``
    returns the moments of the projected DOS.
    GetBandMoment(0) return the integrated DOS. 
    GetBandMoment(1,2) returns the first and second moment of the projected DOS 
    (center and width)

   ``SaveData(filename)``
     saves the data to the filename

   ``GetGracePlot``
     makes a GracePlot (if GracePlot is installed)

.. note:: 
   The eigen values are not relative to the Fermi level, 
   you can get the fermi level for the calculation using 
   calc.GetFermilevel() 

Multicenter orbitals
````````````````````

You add the multicenter orbitals using

 >>> sim.mcen = MultiCenterProjectedDOS()

now you can add a multicenter orbital like

 >>> sim.mcen.add([(0,0,0,1.0),(1,0,0,-1.0)])

this will define an atom1_s-atom2_s orbital.

The general form is [{(atomnr,l,m,weight)}]

You can retrieve information using

 >>> sim.mcen.read(spin=0/1,cutoff="Infinite"/"short")

defaults are spin=0, cutoff="Infinite".


This will return a list of instances of the class 
MultiCenterProjectedDOSTool.
This class have the methods:

 GetBandMoment (GetBandMoment(1) will return the first moment of the band)
 GetCutoff
 GetData
 GetDescription
 GetEFermi
 GetGracePlot
 GetIntegratedDOS  (integrated up til the Fermi level)
 GetPlot (you can combine plot like GetPlot(parent = otherplot), 
 otherplot.Replot())
 GetSpin
 SaveData (SaveData(filename))

The MultiCenterProjectedDOSTool is getting all methods
(except GetDescription, GetSpin) from AtomProjectedDOSTool, which are
returned if you use::

 >>> sim.ados.GetDOS(atoms=[1],angularchannels=['s'],  spin=[spin]))

Note the numbering is here from 1.

You can add the following attributes to sim.mcen:

  sim.mcen.EnergyWindow=(-10,5)
  sim.mcen.EnergyWidth=0.25
  sim.mcen.NumberEnergyPoints=100
  sim.mcen.CutoffRadius=1.5

The direction cosines for which the spherical harmonics are set up are using
the next different atom in the list (cyclic) as direction pointer, so the
z-direction is chosen along the direction to this next atom.
At the moment the rotation matrices is only given in the text file,
you can use
grep 'MUL: Rmatrix' out_o2.txt
to get this information.


Changing convergence parameters
------------------------------- 

The method 
  ``SetConvergenceParameters(**keywords)``

can be used to set convergence parameters for the calculations.
The following keywords are defined:

* energy: 
    Absolute convergence in energy (default 1e-5) 
* density: 
    Convergence in density  (default 0.0001) 
* occupation: 
    Convergence in occupation numbers (default 0.001) 


SetConvergenceParameters(energy=0.00001,density=0.0001,occupation=0.001) 
correspond to the default settings. 

The method 
  ``GetConvergenceParameters()`` 

can be used to retrieve the convergence parameters. 





NetCDF specific methods
-----------------------

Two methods are defined for general reading and writing of ``NetCDF`` 
entries. This provides a method for reading and writing any NetCDF variable 
defined in the `dacapo netcdf`_ manual, so that any settings for the 
dacapo fortran program, not defined by the standard sets of method from 
the `DFT Calculator`_ or from the extra Dacapo method described above, 
can be defined. The methods are listed below:  


``SetNetCDFEntry(variable,att=None,value=None)``: 
     Define a general netcdf entry for this calculation.
     If the netcdf entry is allready defined, the value are
     modified.

    
``GetNetCDFEntry(variable,att=None)``:
     Returns the value of NetCDF entry <variable> or NetCDF attribute 
     ``att``, if specified.

The effect of writing: 

>>> calc = Dacapo() 
>>> calc.SetNumberOfBands(10) 

and 

>>> calc = Dacapo() 
>>> calc.SetNetCDFEntry('ElectronicBands',att='NumberOfBands',value=10) 

is the same. 
Now you can get the information by either using 

>>> print calc.GetNumberOfBands() 

or

>>> print calc.GetNetCDFEntry('ElectronicBands',att='NumberOfBands')


---------------------------------------
Getting information from the calculator
---------------------------------------

Below follow a list of the methods one can use to read information from 
the calculator: 

``GetNumberOBands``:
  Return the number of electronic bands.

``GetXCFunctional``:
  Return the XC-functional identifier ('LDA','PBE' ..)

``GetBZKPoints``
  Return the K-points in the in the Brillouin zone. The K-points are
  in units of the reciprocal space basis.

``GetSpinPolarized``:
  Is it a spin-polarized calculation.

``GetIBZKPoints``:
  Return k-points in the irreducible part of the irreducible
  Brillouin  zone.
  The K-points are in units of the reciprocal space basis

``GetIBZKPointWeights``
  Weights of the k-points. The sum of all weights is one.

``GetEigenvalues(kpt=0, spin=0)``:
  Return the eigenvalues for all bands for the given K-point ``kpt`` and the
  given spin ``spin``. Eigenvalues are relative to the Fermi-level.

``GetOccupationNumbers(kpt=0, spin=0)``:
  Return the occupation numbers for all bands for the given K-point ``kpt``
  and the given spin ``spin``. These numbers fall in the range [0,1] for a
  spin-polarized calculation, and in the range [0,2] otherwise.

``GetWavefunctionArray(band, kpt=0, spin=0)``:
  Return array of wavefunction values. The wavefunction are returned for the
  band ``band``, K-point ``kpt`` and spin ``spin``.

  The wavefunction is returned as  a three-dimensional Numerical array,
  corresponding to the grid use in the calculation.


``GetWavefunctionArrays(kpt=0, spin=0)``:
  Return array of wavefunction values for all bands for the given
  K-point ``kpt``  and the given spin ``spin``.
  The wavefunctions are returned as a three-dimensional Numerical array.

``GetDensityArray(spin=0)``:
  Return array of density values for the given spin ``spin``.



----------------------------------
Working with the Dacapo calculator
----------------------------------

reading from an old calculation
-------------------------------

If you want to restart a calculation or read specific information
from an old calculation, one have to access the information via
the atoms, using:

>>> atoms = Dacapo.ReadAtoms(filename='oldfile')

Now you can get the calculator, attached to the atoms, and
get the information, using:

>>> calc = atoms.GetCalculator()
>>> print calc.GetEigenValues()

A short-hand version of this is:

>>> print
>>> Dacapo.ReadAtoms(filename='oldfile').GetCalculator.GetEigenValues()

Often ASE methods take the atoms as the argument, using the ``GetCalculator``
method
to access the relevant information, i.e the ``VTK`` methods for 3D plotting:

>>> from ASE.Visualization.VTK import VTKPlotWaveFunction
>>> atoms = Dacapo.ReadAtoms(filename='oldfile')
>>> VTKPlotWaveFunction(atoms)


-----------------
Pseudo-potentials 
-----------------

For an overview of the pseudopotentials for dacapo, 
see the `dacapo pseudopotential page`_.


--------------------------------------
Minimization and Dynamics using python
--------------------------------------

Minimization of structures is done using the mimimization 
algorithms implemented in ASE. 
Constrains on atoms are done using ASE ``filters``. 
The next section describes the use of filters.

Setting up a constrain filter
-----------------------------

In the simplest form a filter is just a wrapper around a 
``ListOfAtoms`` object, eg. for a filters constraining 
some atoms in a ``ListOfAtoms``, the filter is just 
a shortened version of the orginal ``ListOfAtoms``. 
For a full description see the `filter`_ in the `ASE manual`_.  


Constraining one atom in a CO-molecule: 

>>> from ASE import Atom, ListOfAtoms
>>> from ASE.Filters import Subset
>>> atoms = ListOfAtoms([Atom('C', (0,   0, 0), tag=117),
...                  Atom('O', (1.2, 0, 0))])
>>> mask=[atom.GetTag() == 1 for atom in atoms]
>>> subset = Subset(atoms,mask) 

Another often used constrain is to  
fix a bond-lenght. For a example of combining the 
subset filter and the fix bond length filter, see the 
`h2o example`_. 
 

Setting up a dynamics
---------------------

Once the filter is defined on can setup a minimization or 
dynamics working that takes this filter as input.

To setup a conjugate-gradient minimizer using the filter ``subset`` defined 
above, one can use:
 
>>> from ASE.Dynamics import ConjugateGradient
>>> from ASE.Dynamics.LineMinimizers.LM1 import LM1
>>> dyn = ConjugateGradient(subset, fmax=0.05,lineMin=LM1(0.05))
>>> dyn.Converge()


For a full description of the minimization and dynamics algorithms 
available see `dynamics`_ in the `ASE manual`_.  


Generating a netcdf trajectory
------------------------------ 

To save the atomic configurations, energies, forces etc, produced from 
the mimimizer defined above,  you need to add a trajectory object to the 
minimizer: 

>>> path = NetCDFTrajectory('CO-trajectory.nc', atoms)
>>> dyn.Attach(path) 

See the `h2o example`_.
For more information on trajectories see `trajectory`_ in the `ASE manual`_.  

------------------------------
Wavefunction and density plots 
------------------------------

See the example `plot wavefunction`:.   

----------------
Wannier orbitals 
----------------

Using the ASE Wannier module one can find the most localized set of Wannier orbitals 
for a dacapo calculation. 
More details will be added.
See the `Wannier example`_. 


---------------------------------------
Nudged Elastic Band Minimum energy path
--------------------------------------- 

See the `Nudged Elastic Band method`_ for a general description. 
See the `Nudged Elastic Band example`_ for a dacapo application of this 
method.
 
