Actions to define on the Input parameters.
This module defines some of the most common actions that a BigDFT user might like to
perform on the input file. Such module therefore set some of the keys of the input
dictionary to the values needed to perform the operations.
Users might also inspire to the actions performed in order to customize the runs in a different way.
All the functions of this module have as first argument inp
, the dictionary of the input parameters.
Many other actions are available in BigDFT code. This module only regroups the most common.
Any of these functionalities might be removed from the input file by the remove()
function.
Note
Any of the action of this module, including the remove()
function, can be also applied
to an instance of the BigDFT.Inputfiles.Inputfile
class, by removing the first argument (inp
).
This adds extra flexibility as the same method may be used to a dictionary instance or to a BigDFT input files.
See the example Example.
Note
Each of the actions here must have default value for the arguments (except the input dictionary inp
).
This is needed for a good behaviour of the function remove.
remove (inp, action) |
Remove action from the input dictionary. |
set_xc (inp[, xc]) |
Set the exchange and correlation approximation |
set_hgrid (inp[, hgrids]) |
Set the wavelet grid spacing. |
set_rmult (inp[, rmult, coarse, fine]) |
Set the wavelet grid extension by modifying the multiplicative radii. |
set_atomic_positions (inp[, posinp]) |
Insert the atomic positions as a part of the input dictionary |
set_mesh_sizes (inp[, ngrids]) |
Constrain the number of grid points in each direction. |
optimize_geometry (inp[, method, nsteps]) |
Optimize the geometry of the system |
spin_polarize (inp[, mpol]) |
Add a collinear spin polarization to the system. |
charge (inp[, charge]) |
Charge the system |
charge_and_polarize (inp) |
Charge the system by removing one electron. |
set_symmetry (inp[, yes]) |
Set the symmetry detection for the charge density and the ionic forces and stressdef set_symmetry(inp,yes=True): |
apply_electric_field (inp[, elecfield]) |
Apply an external electric field on the system |
set_random_inputguess (inp) |
Input orbitals are initialized as random coefficients |
write_orbitals_on_disk (inp[, format]) |
Set the code to write the orbitals on disk in the provided format |
read_orbitals_from_disk (inp) |
Read the orbitals from data directory, if available |
write_density_on_disk (inp) |
Write the charge density on the disk after the last SCF convergence |
calculate_dipole (inp) |
Extract the dipole momenet from the total charge density. |
use_gpu_acceleration (inp) |
Employ gpu acceleration when available, for convolutions and Fock operator |
change_data_directory (inp[, name]) |
Modify the name of the data- directory. |
connect_run_data (inp[, log]) |
Associate the data of the run of a given logfile to the input by retrieving the data directory name of the logfile. |
add_empty_SCF_orbitals (inp[, norbs]) |
Insert norbs empty orbitals in the SCF procedure |
extract_virtual_states (inp, nvirt[, davidson]) |
Extract a given number of empty states after the scf cycle. |
set_electronic_temperature (inp[, kT, T]) |
Define the electronic temperature, in AU (kT ) or K (T ) |
calculate_tddft_coupling_matrix (inp[, tda, …]) |
Perform a casida TDDFT coupling matrix extraction. |
BigDFT.InputActions.
add_empty_SCF_orbitals
(inp, norbs=10)[source]¶Insert norbs
empty orbitals in the SCF procedure
Parameters: | norbs (int) – Number of empty orbitals |
---|
Warning
In linear scaling case, this is only meaningful for the direct minimization approach.
BigDFT.InputActions.
apply_electric_field
(inp, elecfield=[0, 0, 0.001])[source]¶Apply an external electric field on the system
Parameters: | electric (list, float) – Values of the Electric Field in the three directions. Might also be a scalar. |
---|
BigDFT.InputActions.
calculate_dipole
(inp)[source]¶Extract the dipole momenet from the total charge density.
Note
This function is useful for the linear scaling setup as the cubic scaling approach always calculates the charge density multipoles.
BigDFT.InputActions.
calculate_tddft_coupling_matrix
(inp, tda=False, rpa=True, fxc=True)[source]¶Perform a casida TDDFT coupling matrix extraction.
Parameters: |
|
---|
Note
The arguments fxc
and rpa
should not be simultaneously False
.
Warning
Presently the LR-TDDFT casida fxc is only available for LDA functionals in ABINIT flavour.
BigDFT.InputActions.
change_data_directory
(inp, name='')[source]¶Modify the name of the data-
directory.
Useful to grab the orbitals from another directory than the run name
BigDFT.InputActions.
charge
(inp, charge=-1)[source]¶Charge the system
Parameters: | charge (int,float) – value of the charge in units of e (the electron has charge -1). Also accept floating point numbers. |
---|
BigDFT.InputActions.
charge_and_polarize
(inp)[source]¶Charge the system by removing one electron. Assume that the original system is closed shell, thus polarize.
BigDFT.InputActions.
connect_run_data
(inp, log=None)[source]¶Associate the data of the run of a given logfile to the input by retrieving the data directory name of the logfile.
Parameters: | log (Logfile) – instance of a Logfile class |
---|
BigDFT.InputActions.
extract_virtual_states
(inp, nvirt, davidson=False)[source]¶Extract a given number of empty states after the scf cycle.
Parameters: | davidson (bool) – If set to True activates davidson calculation, otherwise Trace Minimization of the Hamiltonian is employed. |
---|
BigDFT.InputActions.
optimize_geometry
(inp, method='FIRE', nsteps=50)[source]¶Optimize the geometry of the system
Parameters: |
|
---|
BigDFT.InputActions.
read_orbitals_from_disk
(inp)[source]¶Read the orbitals from data directory, if available
BigDFT.InputActions.
remove
(inp, action)[source]¶Remove action from the input dictionary.
Remove an action from the input file, thereby restoring the default value, as if the action were not specified.
Parameters: |
|
---|
Example
>>> from Calculators import SystemCalculator as C
>>> code=C()
>>> inp={}
>>> set_xc(inp,'PBE')
>>> write_orbitals_on_disk(inp)
>>> log=code.run(input=inp) # perform calculations
>>> remove(write_orbitals_on_disk) #remove the action
>>> read_orbitals_from_disk(inp)
>>> log2=code.run(input=inp) #this will restart the SCF from the previous orbitals
BigDFT.InputActions.
set_SCF_method
(inp, method='dirmin', mixing_on='density', mixing_scheme='Pulay')[source]¶Set the algorithm for SCF.
Parameters: |
|
---|
Warning
Only the FOE method exhibit asymptotic linear scaling regime.
Todo
Check if the linear scaling case needs another input variable for the mixing of the potential (density)
BigDFT.InputActions.
set_atomic_positions
(inp, posinp=None)[source]¶Insert the atomic positions as a part of the input dictionary
BigDFT.InputActions.
set_electronic_temperature
(inp, kT=0.001, T=0)[source]¶Define the electronic temperature, in AU (kT
) or K (T
)
BigDFT.InputActions.
set_hgrid
(inp, hgrids=0.4)[source]¶Set the wavelet grid spacing.
Parameters: | hgrid (float,list) – list of the grid spacings in the three directions. It might also be a scalar, which implies the same spacing |
---|
BigDFT.InputActions.
set_mesh_sizes
(inp, ngrids=64)[source]¶Constrain the number of grid points in each direction. This is useful when performing periodic system calculations with variable cells which need to be compared each other. In this way the number of degrees of freedom is kept constant throughout the various simuilations.
Parameters: | ngrids (int,list) – list of the number of mesh points in each direction. Might be a scalar. |
---|
BigDFT.InputActions.
set_random_inputguess
(inp)[source]¶Input orbitals are initialized as random coefficients
BigDFT.InputActions.
set_rmult
(inp, rmult=None, coarse=5.0, fine=8.0)[source]¶Set the wavelet grid extension by modifying the multiplicative radii.
Parameters: |
|
---|
BigDFT.InputActions.
set_symmetry
(inp, yes=True)[source]¶Set the symmetry detection for the charge density and the ionic forces and stressdef set_symmetry(inp,yes=True):
Parameters: | yes (bool) – If False the symmetry detection is disabled |
---|
BigDFT.InputActions.
set_xc
(inp, xc='PBE')[source]¶Set the exchange and correlation approximation
Parameters: | xc (str) – the Acronym of the XC approximation |
---|
Todo
Insert the XC codes corresponding to libXC
conventions
BigDFT.InputActions.
spin_polarize
(inp, mpol=1)[source]¶Add a collinear spin polarization to the system.
Parameters: | mpol (int) – spin polarization in Bohr magneton units. |
---|
BigDFT.InputActions.
use_gpu_acceleration
(inp)[source]¶Employ gpu acceleration when available, for convolutions and Fock operator
Todo
Verify what happens when only one of the functionality is enabled at compile-time
BigDFT.InputActions.
write_cubefiles_around_fermi_level
(inp, nplot=1)[source]¶Writes the nplot
orbitals around the fermi level in cube format
Parameters: | nplot (int) – the number of orbitals to print around the fermi level. |
---|
Warning
This is presently meaningful only for a empty states calculation.
Warning
This would work only for the cubic scaling code at present.
BigDFT.InputActions.
write_density_on_disk
(inp)[source]¶Write the charge density on the disk after the last SCF convergence
BigDFT.InputActions.
write_orbitals_on_disk
(inp, format='binary')[source]¶Set the code to write the orbitals on disk in the provided format
Parameters: | format (str) – The format to write the orbitals with. Accepts the strings: * ‘binary’ * ‘text’ * ‘etsf’ (requires etsf-io enabled) |
---|
Todo
Verify if this option works for a linear scaling calulation.