small fixes

src/Potential/poisson.hpp

@@ -221,7 +221,7 @@ inline void Potential::poisson(Periodic_function<double> const& rho)
         /* add contribution from the pseudo-charge */
         poisson_add_pseudo_pw(qmt, qit, const_cast<double_complex*>(&rho.f_pw_local(0)));
 
-        if (check_pseudo_charge) {
+        if (ctx_.control().verification_ >= 2) {
             poisson_sum_G(ctx_.lmmax_rho(), &rho.f_pw_local(0), sbessel_mom_, qit);
 
             double d = 0.0;

src/SDDK/wave_functions.hpp

@@ -19,7 +19,7 @@
 
 /** \file wave_functions.hpp
  *
- *  \brief Contains declaration and implementation of wave_functions class.
+ *  \brief Contains declaration and implementation of Wave_functions class.
  */
 
 #ifndef __WAVE_FUNCTIONS_HPP__

src/SDDK/wf_inner.hpp

@@ -57,7 +57,7 @@ inline void inner(device_t        pu__,
                   int             irow0__,
                   int             jcol0__)
 {
-    PROFILE("sddk::wave_functions::inner");
+    PROFILE("sddk::Wave_functions::inner");
 
     static_assert(std::is_same<T, double>::value || std::is_same<T, double_complex>::value, "wrong type");
 

src/SDDK/wf_trans.hpp

@@ -1,5 +1,9 @@
 /// Linear transformation of the wave-functions.
-/** The transformation matrix is expected in the CPU memory. */
+/** The transformation matrix is expected in the CPU memory. The following operation is performed:
+ *  \f[
+ *     \psi^{out}_{j} = \alpha \sum_{i} \psi^{in}_{i} Z_{ij} + \beta \psi^{out}_{j}
+ *  \f]
+ */
 template <typename T>
 inline void transform(device_t                     pu__,
                       int                          ispn__,
@@ -15,7 +19,7 @@ inline void transform(device_t                     pu__,
                       int                          j0__,
                       int                          n__)
 {
-    PROFILE("sddk::wave_functions::transform");
+    PROFILE("sddk::Wave_functions::transform");
 
     static_assert(std::is_same<T, double>::value || std::is_same<T, double_complex>::value, "wrong type");
 
@@ -140,7 +144,7 @@ inline void transform(device_t                     pu__,
         }
     };
 
-    sddk::timer t1("sddk::wave_functions::transform|init");
+    sddk::timer t1("sddk::Wave_functions::transform|init");
     /* initial values for the resulting wave-functions */
     for (int iv = 0; iv < nwf; iv++) {
         if (beta__ == 0) {
@@ -370,7 +374,6 @@ inline void transform(device_t        pu__,
                       int             j0__,
                       int             n__)
 {
-
     transform<T>(pu__, ispn__, alpha__, {&wf_in__}, i0__, m__, mtrx__, irow0__, jcol0__, beta__, {&wf_out__}, j0__, n__);
 }
 
@@ -387,75 +390,6 @@ inline void transform(device_t        pu__,
                       int             j0__,
                       int             n__)
 {
-
     transform<T>(pu__, ispn__, 1.0, {&wf_in__}, i0__, m__, mtrx__, irow0__, jcol0__, 0.0, {&wf_out__}, j0__, n__);
 }
 
-//== /// Linear transformation of wave-functions.
-//== /** The following operation is performed:
-//==  *  \f[
-//==  *     \psi^{out}_{j} = \alpha \sum_{i} \psi^{in}_{i} Z_{ij} + \beta \psi^{out}_{j}
-//==  *  \f]
-//==  */
-//== template <typename T>
-//== inline void transform(device_t pu__,
-//==                       double          alpha__,
-//==                       wave_functions& wf_in__,
-//==                       int             i0__,
-//==                       int             m__,
-//==                       dmatrix<T>&     mtrx__,
-//==                       int             irow0__,
-//==                       int             jcol0__,
-//==                       double          beta__,
-//==                       wave_functions& wf_out__,
-//==                       int             j0__,
-//==                       int             n__)
-//== {
-//==     transform<T>(pu__, alpha__, {&wf_in__}, i0__, m__, mtrx__, irow0__, jcol0__, beta__, {&wf_out__}, j0__, n__);
-//== }
-//==
-//== template <typename T>
-//== inline void transform(device_t        pu__,
-//==                       wave_functions& wf_in__,
-//==                       int             i0__,
-//==                       int             m__,
-//==                       dmatrix<T>&     mtrx__,
-//==                       int             irow0__,
-//==                       int             jcol0__,
-//==                       wave_functions& wf_out__,
-//==                       int             j0__,
-//==                       int             n__)
-//== {
-//==     transform<T>(pu__, 1.0, {&wf_in__}, i0__, m__, mtrx__, irow0__, jcol0__, 0.0, {&wf_out__}, j0__, n__);
-//== }
-//==
-//== template <typename T>
-//== inline void transform(device_t                     pu__,
-//==                       double                       alpha__,
-//==                       std::vector<Wave_functions*> wf_in__,
-//==                       int                          i0__,
-//==                       int                          m__,
-//==                       dmatrix<T>&                  mtrx__,
-//==                       int                          irow0__,
-//==                       int                          jcol0__,
-//==                       double                       beta__,
-//==                       std::vector<Wave_functions*> wf_out__,
-//==                       int                          j0__,
-//==                       int                          n__)
-//== {
-//==     assert(wf_in__.size() == wf_out__.size());
-//==     for (size_t i = 0; i < wf_in__.size(); i++) {
-//==         assert(wf_in__[i]->num_components() == wf_in__[0]->num_components());
-//==         assert(wf_in__[i]->num_components() == wf_out__[i]->num_components());
-//==     }
-//==     int num_sc = wf_in__[0]->num_components();
-//==     for (int is = 0; is < num_sc; is++) {
-//==         std::vector<wave_functions*> wf_in;
-//==         std::vector<wave_functions*> wf_out;
-//==         for (size_t i = 0; i < wf_in__.size(); i++) {
-//==             wf_in.push_back(&wf_in__[i]->component(is));
-//==             wf_out.push_back(&wf_out__[i]->component(is));
-//==         }
-//==         transform(pu__, alpha__, wf_in, i0__, m__, mtrx__, irow0__, jcol0__, beta__, wf_out, j0__, n__);
-//==     }
-//== }

src/Unit_cell/free_atom.hpp

@@ -38,7 +38,7 @@ class Free_atom : public Atom_type
     /// Radial wave-functions.
     mdarray<double, 2> free_atom_wave_functions_;
     /// Atomic potential.
-    Spline<double> free_atom_potential_;
+    Spline<double> free_atom_potential_; // TODO: rename. this is potential of electronic density
     /// NIST total energy for LDA calculation.
     double NIST_LDA_Etot_{0};
     /// NIST total energy for scalar-relativistic LDA calculation.
@@ -211,7 +211,7 @@ class Free_atom : public Atom_type
             }
 
             for (int i = 0; i < np; i++) {
-                f(i) = (veff[i] - vnuc[i]) * rho(i);
+                f(i) = vrho[i] * rho(i);
             }
             /* kinetic energy */
             energy_kin = eval_sum - fourpi * (f.interpolate().integrate(2) - zn() * rho.integrate(1));
@@ -329,6 +329,11 @@ class Free_atom : public Atom_type
     {
         return free_atom_potential_.at_point(x);
     }
+
+    inline Spline<double>& free_atom_potential()
+    {
+        return free_atom_potential_;
+    }
 };
 } // namespace sirius
 

src/config.h

@@ -73,10 +73,7 @@ const bool test_spinor_wf = false;
 
 const bool hdf5_trace_errors = false;
 
-/// True if pseudocharge moments in LAPW Poisson solver have to be checked and compared with true multipole moments.
-const bool check_pseudo_charge = false;
-
-/// Tru if LAPW Hamiltonian is diagonalized with second-variational approach.
+/// True if LAPW Hamiltonian is diagonalized with second-variational approach.
 const bool use_second_variation = true;
 
 /// True if wave-functions are permanently kept on a GPU and not allocated/deallocated.

src/matching_coefficients.h

@@ -88,7 +88,7 @@ class Matching_coefficients
             /* invert matrix of radial derivatives */
             switch (N) {
                 case 1: {
-                    if (unit_cell_.parameters().control().verification_ > 0) {
+                    if (unit_cell_.parameters().control().verification_ >= 1) {
                         if (std::abs(A(0, 0)) < 1.0 / std::sqrt(unit_cell_.omega())) {   
                             std::stringstream s;
                             s << "Ill defined plane wave matching problem for atom type " << iat << ", l = " << l << std::endl
@@ -103,7 +103,7 @@ class Matching_coefficients
                 case 2: {
                     double det = A(0, 0) * A(1, 1) - A(0, 1) * A(1, 0);
 
-                    if (unit_cell_.parameters().control().verification_ > 0) {
+                    if (unit_cell_.parameters().control().verification_ >= 1) {
                         if (std::abs(det) < 1.0 / std::sqrt(unit_cell_.omega())) {   
                             std::stringstream s;
                             s << "Ill defined plane wave matching problem for atom type " << iat << ", l = " << l << std::endl