SIRIUS 7.5.0
Electronic structure library and applications
dft_ground_state.cpp
Go to the documentation of this file.
1// Copyright (c) 2013-2018 Anton Kozhevnikov, Thomas Schulthess
2// All rights reserved.
3//
4// Redistribution and use in source and binary forms, with or without modification, are permitted provided that
5// the following conditions are met:
6//
7// 1. Redistributions of source code must retain the above copyright notice, this list of conditions and the
8// following disclaimer.
9// 2. Redistributions in binary form must reproduce the above copyright notice, this list of conditions
10// and the following disclaimer in the documentation and/or other materials provided with the distribution.
11//
12// THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR IMPLIED
13// WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A
14// PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE FOR
15// ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
16// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
17// CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR
18// OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
19
20/** \file dft_ground_state.cpp
21 *
22 * \brief Contains implementation of sirius::DFT_ground_state class.
23 */
24
25#include <iomanip>
26#include "dft_ground_state.hpp"
27#include "core/profiler.hpp"
30
31namespace sirius {
32
33void
35{
36 PROFILE("sirius::DFT_ground_state::initial_state");
37
40 if (!ctx_.full_potential()) {
41 if (ctx_.cfg().parameters().precision_wf() == "fp32") {
42#if defined(SIRIUS_USE_FP32)
45#else
46 RTE_THROW("not compiled with FP32 support");
47#endif
48
49 } else {
52 }
53 }
54}
55
56void
57DFT_ground_state::create_H0()
58{
59 PROFILE("sirius::DFT_ground_state::create_H0");
60
61 H0_ = std::shared_ptr<Hamiltonian0<double>>(new Hamiltonian0<double>(potential_, true));
62}
63
64void
66{
67 PROFILE("sirius::DFT_ground_state::update");
68
69 ctx_.update();
70 kset_.update();
73
74 if (!ctx_.full_potential()) {
76 }
77}
78
79double
80DFT_ground_state::energy_kin_sum_pw() const
81{
82 double ekin{0};
83
84 for (auto it : kset_.spl_num_kpoints()) {
85 auto kp = kset_.get<double>(it.i);
86
87 #pragma omp parallel for schedule(static) reduction(+:ekin)
88 for (int igloc = 0; igloc < kp->num_gkvec_loc(); igloc++) {
89 auto Gk = kp->gkvec().gkvec_cart<index_domain_t::local>(igloc);
90
91 double d{0};
92 for (int ispin = 0; ispin < ctx_.num_spins(); ispin++) {
93 for (int i = 0; i < kp->num_occupied_bands(ispin); i++) {
94 double f = kp->band_occupancy(i, ispin);
95 auto z = kp->spinor_wave_functions().pw_coeffs(igloc, wf::spin_index(ispin), wf::band_index(i));
96 d += f * (std::pow(z.real(), 2) + std::pow(z.imag(), 2));
97 }
98 }
99 if (kp->gkvec().reduced()) {
100 d *= 2;
101 }
102 ekin += 0.5 * d * kp->weight() * Gk.length2();
103 } // igloc
104 } // ikloc
105 ctx_.comm().allreduce(&ekin, 1);
106 return ekin;
107}
108
109double
110DFT_ground_state::total_energy() const
111{
113}
114
115json
116DFT_ground_state::serialize()
117{
118 return energy_dict(ctx_, kset_, density_, potential_, ewald_energy_, this->scf_correction_energy_);
119}
120
121/// A quick check of self-constent density in case of pseudopotential.
122json
124{
125 if (ctx_.full_potential()) {
126 return json();
127 }
128
129 auto gs0 = energy_dict(ctx_, kset_, density_, potential_, ewald_energy_);
130
131 /* create new potential */
132 Potential pot(ctx_);
133 /* generate potential from existing density */
134 bool transform_to_rg{true};
135 pot.generate(density_, ctx_.use_symmetry(), transform_to_rg);
136 /* create new Hamiltonian */
137 bool precompute_lapw{true};
138 Hamiltonian0<double> H0(pot, precompute_lapw);
139 /* initialize the subspace */
141 /* find new wave-functions */
142 ::sirius::diagonalize<double, double>(H0, kset_, ctx_.cfg().settings().itsol_tol_min());
143 /* find band occupancies */
145 /* generate new density from the occupied wave-functions */
146 bool add_core{true};
147 /* create new density */
148 Density rho(ctx_);
149 rho.generate<double>(kset_, ctx_.use_symmetry(), add_core, transform_to_rg);
150
151 auto gs1 = energy_dict(ctx_, kset_, rho, pot, ewald_energy_);
152
153 auto calc_rms = [&](Field4D& a, Field4D& b) -> double
154 {
155 double rms{0};
156 for (int j = 0; j < ctx_.num_mag_dims() + 1; j++) {
157 for (int ig = 0; ig < ctx_.gvec().count(); ig++) {
158 rms += std::pow(std::abs(a.component(j).rg().f_pw_local(ig) - b.component(j).rg().f_pw_local(ig)), 2);
159 }
160 }
161 ctx_.comm().allreduce(&rms, 1);
162 return std::sqrt(rms / ctx_.gvec().num_gvec());
163 };
164
165 double rms = calc_rms(density_, rho);
166 double rms_veff = calc_rms(potential_, pot);
167
168 json dict;
169 dict["rms"] = rms;
170 dict["detot"] = gs0["energy"]["total"].get<double>() - gs1["energy"]["total"].get<double>();
171
172 rms_veff = std::sqrt(rms_veff / ctx_.gvec().num_gvec());
173
174 if (ctx_.verbosity() >= 1) {
175 RTE_OUT(ctx_.out()) << "RMS_rho: " << dict["rms"].get<double>() << std::endl
176 << "RMS_veff: " << rms_veff << std::endl
177 << "Eold: " << gs0["energy"]["total"].get<double>()
178 << " Enew: " << gs1["energy"]["total"].get<double>() << std::endl;
179
180 std::vector<std::string> labels({"total", "vha", "vxc", "exc", "bxc", "veff", "eval_sum", "kin", "ewald",
181 "vloc", "scf_correction", "entropy_sum"});
182
183 for (auto e: labels) {
184 RTE_OUT(ctx_.out()) << "energy component: " << e << ", diff: " <<
185 std::abs(gs0["energy"][e].get<double>() - gs1["energy"][e].get<double>()) << std::endl;
186 }
187 }
188
189 return dict;
190}
191
192json
193DFT_ground_state::find(double density_tol__, double energy_tol__, double iter_solver_tol__, int num_dft_iter__,
194 bool write_state__)
195{
196 PROFILE("sirius::DFT_ground_state::scf_loop");
197
198 auto tstart = std::chrono::high_resolution_clock::now();
199
200 double eold{0}, rms{0};
201
202 density_.mixer_init(ctx_.cfg().mixer());
203
204 int num_iter{-1};
205 std::vector<double> rms_hist;
206 std::vector<double> etot_hist;
207
208 Density rho1(ctx_);
209
210 std::stringstream s;
211 s << "density_tol : " << density_tol__ << std::endl
212 << "energy_tol : " << energy_tol__ << std::endl
213 << "iter_solver_tol (initial) : " << iter_solver_tol__ << std::endl
214 << "iter_solver_tol (target) : " << ctx_.cfg().settings().itsol_tol_min() << std::endl
215 << "num_dft_iter : " << num_dft_iter__;
216 ctx_.message(1, __func__, s);
217
218 for (int iter = 0; iter < num_dft_iter__; iter++) {
219 PROFILE("sirius::DFT_ground_state::scf_loop|iteration");
220 std::stringstream s;
221 s << std::endl;
222 s << "+------------------------------+" << std::endl
223 << "| SCF iteration " << std::setw(3) << iter << " out of " << std::setw(3) << num_dft_iter__ << '|' << std::endl
224 << "+------------------------------+" << std::endl;
225 ctx_.message(2, __func__, s);
226
228
229 if (ctx_.cfg().parameters().precision_wf() == "fp32") {
230#if defined(SIRIUS_USE_FP32)
232 /* find new wave-functions */
233 if (ctx_.cfg().parameters().precision_hs() == "fp32") {
234 result = sirius::diagonalize<float, float>(H0, kset_, iter_solver_tol__);
235 } else {
236 result = sirius::diagonalize<float, double>(H0, kset_, iter_solver_tol__);
237 }
238 /* find band occupancies */
240 /* generate new density from the occupied wave-functions */
241 density_.generate<float>(kset_, ctx_.use_symmetry(), true, true);
242#else
243 RTE_THROW("not compiled with FP32 support");
244#endif
245 } else {
247 /* find new wave-functions */
248 result = sirius::diagonalize<double, double>(H0, kset_, iter_solver_tol__);
249 /* find band occupancies */
251 /* generate new density from the occupied wave-functions */
252 density_.generate<double>(kset_, ctx_.use_symmetry(), true, true);
253 }
254
255 double e1 = energy_potential(density_, potential_);
256 copy(density_, rho1);
257
258 /* mix density */
259 rms = density_.mix();
260
261 double eha_res = density_residual_hartree_energy(density_, rho1);
262
263 /* estimate new tolerance of the iterative solver */
264 double tol = rms;
265 if (ctx_.cfg().mixer().use_hartree()) {
266 //tol = rms * rms / std::max(1.0, unit_cell_.num_electrons());
267 tol = eha_res / std::max(1.0, unit_cell_.num_electrons());
268 }
269 tol = std::min(ctx_.cfg().settings().itsol_tol_scale()[0] * tol,
270 ctx_.cfg().settings().itsol_tol_scale()[1] * iter_solver_tol__);
271 /* tolerance can't be too small */
272 iter_solver_tol__ = std::max(ctx_.cfg().settings().itsol_tol_min(), tol);
273
274 bool iter_solver_converged{true};
275 if (ctx_.cfg().iterative_solver().type() != "exact") {
276 iter_solver_converged = (tol <= ctx_.cfg().settings().itsol_tol_min());
277 }
278
279#if defined(SIRIUS_USE_FP32)
280 if (ctx_.cfg().parameters().precision_gs() != "auto") {
281 /* if the final precision is not equal to the current precision */
282 if (ctx_.cfg().parameters().precision_gs() == "fp64" && ctx_.cfg().parameters().precision_wf() == "fp32") {
283 /* if we reached the mimimum tolerance for fp32 */
284 if ((ctx_.cfg().settings().fp32_to_fp64_rms() == 0 && iter_solver_tol__ <= ctx_.cfg().settings().itsol_tol_min()) ||
285 (rms < ctx_.cfg().settings().fp32_to_fp64_rms())) {
286 std::cout << "switching to FP64" << std::endl;
287 ctx_.cfg().unlock();
288 ctx_.cfg().settings().itsol_tol_min(std::numeric_limits<double>::epsilon() * 10);
289 ctx_.cfg().parameters().precision_wf("fp64");
290 ctx_.cfg().parameters().precision_hs("fp64");
291 ctx_.cfg().lock();
292
293 for (auto it : kset_.spl_num_kpoints()) {
294 for (int ispn = 0; ispn < ctx_.num_spins(); ispn++) {
295 wf::copy(sddk::memory_t::host, kset_.get<float>(it.i)->spinor_wave_functions(),
297 kset_.get<double>(it.i)->spinor_wave_functions(), wf::spin_index(ispn),
299 }
300 }
301 for (int ik = 0; ik < kset_.num_kpoints(); ik++) {
302 for (int ispn = 0; ispn < ctx_.num_spinors(); ispn++) {
303 for (int j = 0; j < ctx_.num_bands(); j++) {
304 kset_.get<double>(ik)->band_energy(j, ispn, kset_.get<float>(ik)->band_energy(j, ispn));
305 kset_.get<double>(ik)->band_occupancy(j, ispn, kset_.get<float>(ik)->band_occupancy(j, ispn));
306 }
307 }
308 }
309 }
310 }
311 }
312#endif
313 if (ctx_.cfg().control().verification() >= 1) {
314 /* check number of electrons */
316 }
317
318 /* compute new potential */
320
321 if (!ctx_.full_potential() && ctx_.cfg().control().verification() >= 2) {
322 if (ctx_.verbosity() >= 1) {
323 RTE_OUT(ctx_.out()) << "checking functional derivative of Exc\n";
324 }
325 sirius::check_xc_potential(density_);
326 }
327
328 if (ctx_.cfg().parameters().use_scf_correction()) {
329 double e2 = energy_potential(rho1, potential_);
330 this->scf_correction_energy_ = e2 - e1;
331 }
332
333 /* compute new total energy for a new density */
334 double etot = total_energy();
335
336 etot_hist.push_back(etot);
337
338 rms_hist.push_back(rms);
339
340 /* write some information */
341 std::stringstream out;
342 out << std::endl;
343 print_info(out);
344 out << std::endl;
345 out << "iteration : " << iter << ", RMS : " << std::setprecision(12) << std::scientific << rms
346 << ", energy difference : " << std::setprecision(12) << std::scientific << etot - eold;
347 if (!ctx_.full_potential()) {
348 out << std::endl
349 << "Hartree energy of density residual : " << eha_res << std::endl
350 << "bands are converged : " << boolstr(result.converged);
351 }
352 if (ctx_.cfg().iterative_solver().type() != "exact") {
353 out << std::endl
354 << "iterative solver converged : " << boolstr(iter_solver_converged);
355 }
356
357 ctx_.message(2, __func__, out);
358 /* check if the calculation has converged */
359 bool converged{true};
360 converged = (std::abs(eold - etot) < energy_tol__) && result.converged && iter_solver_converged;
361 if (ctx_.cfg().mixer().use_hartree()) {
362 converged = converged && (eha_res < density_tol__);
363 } else {
364 converged = converged && (rms < density_tol__);
365 }
366 if (converged) {
367 std::stringstream out;
368 out << std::endl;
369 out << "converged after " << iter + 1 << " SCF iterations!" << std::endl;
370 ctx_.message(1, __func__, out);
372 num_iter = iter;
373 break;
374 }
375
376 eold = etot;
377 }
378 std::stringstream out;
379 out << std::endl;
380 print_info(out);
381 ctx_.message(1, __func__, out);
382
383 if (write_state__) {
384 ctx_.create_storage_file(storage_file_name);
385 if (ctx_.full_potential()) { // TODO: why this is necessary?
386 density_.rho().rg().fft_transform(-1);
387 for (int j = 0; j < ctx_.num_mag_dims(); j++) {
388 density_.mag(j).rg().fft_transform(-1);
389 }
390 }
391 potential_.save(storage_file_name);
392 density_.save(storage_file_name);
393 // kset_.save(storage_file_name);
394 }
395
396 auto tstop = std::chrono::high_resolution_clock::now();
397
398 auto dict = serialize();
399
400 /* check density */
401 if (num_iter >= 0) {
402 density_.rho().rg().fft_transform(1);
403 double rho_min{1e100};
404 for (int ir = 0; ir < density_.rho().rg().spfft().local_slice_size(); ir++) {
405 rho_min = std::min(rho_min, density_.rho().rg().value(ir));
406 }
407 dict["rho_min"] = rho_min;
408 ctx_.comm().allreduce<double, mpi::op_t::min>(&rho_min, 1);
409 }
410
411 dict["scf_time"] = std::chrono::duration_cast<std::chrono::duration<double>>(tstop - tstart).count();
412 dict["etot_history"] = etot_hist;
413 if (num_iter >= 0) {
414 dict["converged"] = true;
415 dict["num_scf_iterations"] = num_iter;
416 dict["rms_history"] = rms_hist;
417 } else {
418 dict["converged"] = false;
419 }
420
421 if (env::check_scf_density()) {
423 }
424
425 // dict["volume"] = ctx.unit_cell().omega() * std::pow(bohr_radius, 3);
426 // dict["volume_units"] = "angstrom^3";
427 // dict["energy"] = dft.total_energy() * ha2ev;
428 // dict["energy_units"] = "eV";
429
430 return dict;
431}
432
433void
434DFT_ground_state::print_info(std::ostream& out__) const
435{
436 double evalsum1 = kset_.valence_eval_sum();
437 double evalsum2 = core_eval_sum(ctx_.unit_cell());
438 double s_sum = kset_.entropy_sum();
439 double ekin = energy_kin(ctx_, kset_, density_, potential_);
440 double evxc = energy_vxc(density_, potential_);
441 double eexc = energy_exc(density_, potential_);
442 double ebxc = energy_bxc(density_, potential_);
443 double evha = energy_vha(potential_);
444 double hub_one_elec = one_electron_energy_hubbard(density_, potential_);
445 double etot = total_energy();
446 double gap = kset_.band_gap() * ha2ev;
447 double ef = kset_.energy_fermi();
448 double enuc = energy_enuc(ctx_, potential_);
449
450 double one_elec_en = evalsum1 - (evxc + evha + ebxc);
451
452 if (ctx_.electronic_structure_method() == electronic_structure_method_t::pseudopotential) {
453 one_elec_en -= potential_.PAW_one_elec_energy(density_);
454 one_elec_en -= hub_one_elec;
455 }
456
457 density_.print_info(out__);
458
459 out__ << std::endl;
460 out__ << "Energy" << std::endl
461 << hbar(80, '-') << std::endl;
462
463 auto write_energy = [&](std::string label__, double value__) {
464 out__ << std::left << std::setw(30) << label__ << " : " << std::right << std::setw(16) << std::setprecision(8)
465 << std::fixed << value__ << std::endl;
466 };
467
468 auto write_energy2 = [&](std::string label__, double value__)
469 {
470 out__ << std::left << std::setw(30) << label__ << " : " << std::right
471 << std::setw(16) << std::setprecision(8) << std::fixed << value__ << " (Ha), "
472 << std::setw(16) << std::setprecision(8) << std::fixed << value__ * 2 << " (Ry)" << std::endl;
473 };
474
475 write_energy("valence_eval_sum", evalsum1);
476 if (ctx_.full_potential()) {
477 write_energy("core_eval_sum", evalsum2);
478 write_energy("kinetic energy", ekin);
479 write_energy("enuc", enuc);
480 }
481 write_energy("<rho|V^{XC}>", evxc);
482 write_energy("<rho|E^{XC}>", eexc);
483 write_energy("<mag|B^{XC}>", ebxc);
484 write_energy("<rho|V^{H}>", evha);
485 if (!ctx_.full_potential()) {
486 write_energy2("one-electron contribution", one_elec_en); // eband + deband in QE
487 write_energy("hartree contribution", 0.5 * evha);
488 write_energy("xc contribution", eexc);
489 write_energy("ewald contribution", ewald_energy_);
490 write_energy("PAW contribution", potential_.PAW_total_energy(density_));
491 }
492 write_energy("smearing (-TS)", s_sum);
493 write_energy("SCF correction", this->scf_correction_energy_);
494 if (ctx_.hubbard_correction()) {
495 auto e = ::sirius::energy(density_.occupation_matrix());
496 write_energy2("Hubbard energy", e);
497 write_energy2("Hubbard one-el contribution", hub_one_elec);
498 }
499 write_energy2("Total energy", etot);
500 out__ << std::endl;
501 write_energy("band gap (eV)", gap);
502 write_energy("Efermi", ef);
503}
504
505} // namespace sirius
namespace for Niels Lohmann
Potential potential_
Instance of the Potential class.
double scf_correction_energy_
Correction to total energy from the SCF density minimisation.
void initial_state()
Generate initial density, potential and a subspace of wave-functions.
void print_info(std::ostream &out__) const
Print the basic information (total energy, charges, moments, etc.).
Density density_
Instance of the Density class.
Simulation_context & ctx_
Context of simulation.
K_point_set & kset_
Set of k-points that are used to generate density.
double ewald_energy_
Store Ewald energy which is computed once and which doesn't change during the run.
json find(double density_tol, double energy_tol, double initial_tolerance, int num_dft_iter, bool write_state)
Run the SCF ground state calculation and find a total energy minimum.
std::shared_ptr< Hamiltonian0< double > > H0_
k-point independent part of the Hamiltonian.
void update()
Update the parameters after the change of lattice vectors or atomic positions.
Unit_cell & unit_cell_
Alias of the unit cell.
json check_scf_density()
A quick check of self-constent density in case of pseudopotential.
Generate charge density and magnetization from occupied spinor wave-functions.
Definition: density.hpp:214
void mixer_init(config_t::mixer_t const &mixer_cfg__)
Initialize density mixer.
Definition: density.cpp:1830
void update()
Update internal parameters after a change of lattice vectors or atomic coordinates.
Definition: density.cpp:120
bool check_num_electrons() const
Check total density for the correct number of electrons.
Definition: density.cpp:1057
double mix()
Mix new density.
Definition: density.cpp:1933
auto const & rho() const
Return const reference to charge density (scalar functions).
Definition: density.hpp:416
void initial_density()
Generate initial charge density and magnetization.
Definition: density.cpp:148
void generate(K_point_set const &ks__, bool symmetrize__, bool add_core__, bool transform_to_rg__)
Generate full charge density (valence + core) and magnetization from the wave functions.
Definition: density.cpp:1088
Four-component function consisting of scalar and vector parts.
Definition: field4d.hpp:40
Represent the k-point independent part of Hamiltonian.
Definition: hamiltonian.hpp:75
int num_kpoints() const
Return total number of k-points.
double entropy_sum() const
Return entropy contribution from smearing store in Kpoint<T>.
double valence_eval_sum() const
Return sum of valence eigen-values store in Kpoint<T>.
void find_band_occupancies()
Find Fermi energy and band occupation numbers.
void update()
Update k-points after moving atoms or changing the lattice vectors.
Generate effective potential from charge density and magnetization.
Definition: potential.hpp:46
void generate(Density const &density__, bool use_sym__, bool transform_to_rg__)
Generate effective potential and magnetic field from charge density and magnetization.
Definition: potential.cpp:241
void update()
Recompute some variables that depend on atomic positions or the muffin-tin radius.
Definition: potential.cpp:152
auto const & gvec() const
Return const reference to Gvec object.
void update()
Update context after setting new lattice vectors or atomic coordinates.
std::ostream & out() const
Return output stream.
mpi::Communicator const & comm() const
Total communicator of the simulation.
void message(int level__, char const *label__, std::stringstream const &s) const
Print message from the stringstream.
int num_spins() const
Number of spin components.
int num_spinors() const
Number of spinor wave-functions labeled by a sinlge band index.
int num_bands(int num_bands__)
Set the number of bands.
int num_mag_dims() const
Number of dimensions in the magnetization vector.
bool use_symmetry() const
Get a use_symmetry flag.
double num_electrons() const
Total number of electrons (core + valence).
Definition: unit_cell.hpp:369
Horisontal bar.
void allreduce(T *buffer__, int count__) const
Perform the in-place (the output buffer is used as the input buffer) all-to-all reduction.
Describe a range of bands.
Contains definition and partial implementation of sirius::DFT_ground_state class.
Entry point for Hamiltonain diagonalization.
Create intial subspace from atomic-like wave-functions.
void copy(sddk::memory_t mem__, Wave_functions< T > const &in__, wf::spin_index s_in__, wf::band_range br_in__, Wave_functions< F > &out__, wf::spin_index s_out__, wf::band_range br_out__)
Copy wave-functions.
Namespace of the SIRIUS library.
Definition: sirius.f90:5
double core_eval_sum(Unit_cell const &unit_cell)
Return eigen-value sum of core states.
Definition: energy.cpp:125
double total_energy(Simulation_context const &ctx, K_point_set const &kset, Density const &density, Potential const &potential, double ewald_energy)
Total energy of the electronic subsystem.
Definition: energy.cpp:186
@ pseudopotential
Pseudopotential (ultrasoft, norm-conserving, PAW).
double energy_vxc(Density const &density, Potential const &potential)
Returns exchange correlation potential.
Definition: energy.cpp:76
double energy_vha(Potential const &potential)
Returns Hatree potential.
Definition: energy.cpp:88
const double ha2ev
Hartree in electron-volt units.
Definition: constants.hpp:54
double energy_kin(Simulation_context const &ctx, K_point_set const &kset, Density const &density, Potential const &potential)
Return kinetic energy.
Definition: energy.cpp:147
double energy_enuc(Simulation_context const &ctx, Potential const &potential)
Return nucleus energy in the electrostatic field.
Definition: energy.cpp:104
void initialize_subspace(Hamiltonian_k< T > const &Hk__, K_point< T > &kp__, int num_ao__)
Initialize the wave-functions subspace at a given k-point.
double energy_bxc(const Density &density, const Potential &potential)
TODO doc.
Definition: energy.cpp:94
double energy_exc(Density const &density, Potential const &potential)
Returns exchange correlation energy.
Definition: energy.cpp:82
double ewald_energy(const Simulation_context &ctx, const fft::Gvec &gvec, const Unit_cell &unit_cell)
Compute the ion-ion electrostatic energy using Ewald method.
Definition: energy.cpp:30
A time-based profiler.