Skip to content

Commit

Permalink
Merge pull request #75 from electronic-structure/develop
Browse files Browse the repository at this point in the history
Develop
  • Loading branch information
electronic-structure authored Aug 16, 2017
2 parents 1b01d37 + a9f9835 commit f1df939
Show file tree
Hide file tree
Showing 92 changed files with 2,483 additions and 12,882 deletions.
2 changes: 1 addition & 1 deletion apps/dft_loop/sirius.scf.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -114,7 +114,7 @@ double ground_state(Simulation_context& ctx,
double e1 = dict["ground_state"]["energy"]["total"];
double e2 = dict_ref["ground_state"]["energy"]["total"];

if (std::abs(e1 - e2) > 1e-7) {
if (std::abs(e1 - e2) > 1e-6) {
printf("total energy is different\n");
sirius::terminate(1);
}
Expand Down
24 changes: 15 additions & 9 deletions apps/upf/plot_uspp.py
Original file line number Diff line number Diff line change
Expand Up @@ -11,23 +11,29 @@
fin.close()

rgrid = jin['pseudo_potential']['radial_grid']
print("number of radial grid points: %i"%len(rgrid))

z = jin['pseudo_potential']['header']['z_valence']

vloc = jin['pseudo_potential']['local_potential']
ir0=0
for i in range(len(rgrid)):
vloc[i] = (vloc[i] * rgrid[i] + z) * rgrid[i]
if rgrid[i] < 10: ir0 = i
#vloc = jin['pseudo_potential']['local_potential']
#ir0=0
#for i in range(len(rgrid)):
# vloc[i] = (vloc[i] * rgrid[i] + z) * rgrid[i]
# if rgrid[i] < 10: ir0 = i


#for b in jin['pseudo_potential']['beta_projectors']:
for b in jin['pseudo_potential']['beta_projectors']:
np = len(b['radial_function'])
print("number of data points: %i"%np)
#for b in jin['pseudo_potential']['atomic_wave_functions']:
# plt.plot(rgrid[0:len(b['radial_function'])], b['radial_function'], "-", linewidth = 2.0)
#plt.plot(rgrid[0:len(b['radial_function'])], b['radial_function'], "-", linewidth = 2.0)
plt.plot(rgrid[0:np], b['radial_function'][0:np], "-", linewidth = 2.0)

#plt.plot(rgrid, jin['pseudo_potential']['core_charge_density'], '-', linewidth=2.0)
#plt.plot(rgrid, jin['pseudo_potential']['total_charge_density'], '-', linewidth=2.0)
plt.plot(rgrid[ir0:], vloc[ir0:], '-', linewidth=2.0)
#plt.plot(rgrid[ir0:], vloc[ir0:], '-', linewidth=2.0)
plt.grid(which = "major", axis = "both")
plt.xlim(xmax=rgrid[-1])
#plt.xlim(xmax=rgrid[-1])

fname = os.path.splitext(os.path.basename(sys.argv[1]))[0]
plt.savefig(fname+".pdf", format="pdf")
Expand Down
2 changes: 1 addition & 1 deletion doc/doxygen.cfg
Original file line number Diff line number Diff line change
Expand Up @@ -38,7 +38,7 @@ PROJECT_NAME = "SIRIUS"
# could be handy for archiving the generated documentation or if some version
# control system is used.

PROJECT_NUMBER = "4.7"
PROJECT_NUMBER = "4.8"

# Using the PROJECT_BRIEF tag one can provide an optional one line description
# for a project that appears at the top of each page and should give viewer a
Expand Down
6 changes: 4 additions & 2 deletions examples/pp-pw/STRESS/Cu_paw/sirius.json
Original file line number Diff line number Diff line change
Expand Up @@ -5,7 +5,9 @@
"!std_evp_solver_type" : "lapack",
"!gen_evp_solver_type" : "lapack",
"verbosity": 1,
"reduce_gvec": true
"reduce_gvec": true,
"print_stress" : true,
"print_forces" : true
},

"parameters" : {
Expand All @@ -16,7 +18,7 @@

"smearing_width" : 0.02,

"use_symmetry" : 1,
"use_symmetry" : true,

"num_mag_dims" : 0,

Expand Down
9 changes: 5 additions & 4 deletions examples/pp-pw/STRESS/LiF_paw/sirius.json
Original file line number Diff line number Diff line change
Expand Up @@ -6,7 +6,7 @@
"!std_evp_solver_type" : "lapack",
"!gen_evp_solver_type" : "lapack",
"reduce_gvec": true,
"verbosity": 1,
"verbosity": 0,
"print_stress": true,
"print_forces": true
},
Expand All @@ -32,7 +32,8 @@

"num_dft_iter" : 100,

"ngridk" : [4,4,4]
"ngridk" : [4,4,4],
"reduce_gvec": 1
},

"iterative_solver" : {
Expand All @@ -55,8 +56,8 @@
"atom_types" : ["Li", "F"],

"atom_files" : {
"Li" : "Li.pz-s-kjpaw_psl.0.2.1.UPF.json",
"F" : "F.pz-n-kjpaw_psl.0.1.UPF.json"
"Li" : "li_lda_v1.4.uspp.F.UPF.json",
"F" : "f_lda_v1.4.uspp.F.UPF.json"
},

"atoms" : {
Expand Down
4 changes: 3 additions & 1 deletion examples/pp-pw/SrVO3_paw/sirius.json
Original file line number Diff line number Diff line change
Expand Up @@ -3,7 +3,9 @@
"cyclic_block_size" : 16,
"processing_unit" : "cpu",
"std_evp_solver_type" : "lapack",
"gen_evp_solver_type" : "lapack"
"gen_evp_solver_type" : "lapack",
"print_stress" : true,
"print_forces" : true
},

"parameters" : {
Expand Down
63 changes: 44 additions & 19 deletions src/Band/diagonalize.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -85,7 +85,7 @@ inline void Band::diag_fv_full_potential_exact(K_point* kp, Potential const& pot
t.stop();
kp->set_fv_eigen_values(&eval[0]);

if (ctx_.control().verbosity_ >= 3 && kp->comm().rank() == 0) {
if (ctx_.control().verbosity_ >= 4 && kp->comm().rank() == 0) {
for (int i = 0; i < ctx_.num_fv_states(); i++) {
DUMP("eval[%i]=%20.16f", i, eval[i]);
}
Expand Down Expand Up @@ -346,7 +346,7 @@ inline void Band::get_singular_components(K_point* kp__) const
/* setup eigen-value problem
* N is the number of previous basis functions
* n is the number of new basis functions */
set_subspace_mtrx(N, n, phi, ophi, ovlp, ovlp_old);
set_subspace_mtrx(1, N, n, phi, ophi, ovlp, ovlp_old);

/* increase size of the variation space */
N += n;
Expand All @@ -362,10 +362,12 @@ inline void Band::get_singular_components(K_point* kp__) const
TERMINATE(s);
}

if (ctx_.control().verbosity_ > 2 && kp__->comm().rank() == 0) {
if (ctx_.control().verbosity_ >= 3 && kp__->comm().rank() == 0) {
DUMP("step: %i, current subspace size: %i, maximum subspace size: %i", k, N, num_phi);
for (int i = 0; i < ncomp; i++) {
DUMP("eval[%i]=%20.16f, diff=%20.16f", i, eval[i], std::abs(eval[i] - eval_old[i]));
if (ctx_.control().verbosity_ >= 4) {
for (int i = 0; i < ncomp; i++) {
DUMP("eval[%i]=%20.16f, diff=%20.16f", i, eval[i], std::abs(eval[i] - eval_old[i]));
}
}
}

Expand All @@ -380,14 +382,14 @@ inline void Band::get_singular_components(K_point* kp__) const
sddk::timer t1("sirius::Band::diag_fv_full_potential_davidson|update_phi");
/* recompute wave-functions */
/* \Psi_{i} = \sum_{mu} \phi_{mu} * Z_{mu, i} */
transform(phi, 0, N, evec, 0, 0, psi, 0, ncomp);
transform(ctx_.processing_unit(), phi, 0, N, evec, 0, 0, psi, 0, ncomp);

/* exit the loop if the eigen-vectors are converged or this is a last iteration */
if (n <= itso.min_num_res_ || k == (itso.num_steps_ - 1)) {
break;
}
else { /* otherwise, set Psi as a new trial basis */
if (ctx_.control().verbosity_ > 2 && kp__->comm().rank() == 0) {
if (ctx_.control().verbosity_ >= 3 && kp__->comm().rank() == 0) {
DUMP("subspace size limit reached");
}

Expand Down Expand Up @@ -537,7 +539,7 @@ inline void Band::diag_fv_full_potential_davidson(K_point* kp) const
/* number of newly added basis functions */
int n = nlo + ncomp + num_bands;

if (ctx_.control().verbosity_ > 2 && kp->comm().rank() == 0) {
if (ctx_.control().verbosity_ >= 3 && kp->comm().rank() == 0) {
DUMP("iterative solver tolerance: %18.12f", ctx_.iterative_solver_tolerance());
}

Expand All @@ -555,7 +557,7 @@ inline void Band::diag_fv_full_potential_davidson(K_point* kp) const
/* setup eigen-value problem
* N is the number of previous basis functions
* n is the number of new basis functions */
set_subspace_mtrx(N, n, phi, hphi, hmlt, hmlt_old);
set_subspace_mtrx(1, N, n, phi, hphi, hmlt, hmlt_old);

/* increase size of the variation space */
N += n;
Expand All @@ -571,10 +573,12 @@ inline void Band::diag_fv_full_potential_davidson(K_point* kp) const
TERMINATE(s);
}

if (ctx_.control().verbosity_ > 2 && kp->comm().rank() == 0) {
if (ctx_.control().verbosity_ >= 3 && kp->comm().rank() == 0) {
DUMP("step: %i, current subspace size: %i, maximum subspace size: %i", k, N, num_phi);
for (int i = 0; i < num_bands; i++) {
DUMP("eval[%i]=%20.16f, diff=%20.16f", i, eval[i], std::abs(eval[i] - eval_old[i]));
if (ctx_.control().verbosity_ >= 4) {
for (int i = 0; i < num_bands; i++) {
DUMP("eval[%i]=%20.16f, diff=%20.16f", i, eval[i], std::abs(eval[i] - eval_old[i]));
}
}
}

Expand All @@ -589,14 +593,14 @@ inline void Band::diag_fv_full_potential_davidson(K_point* kp) const
sddk::timer t1("sirius::Band::diag_fv_full_potential_davidson|update_phi");
/* recompute wave-functions */
/* \Psi_{i} = \sum_{mu} \phi_{mu} * Z_{mu, i} */
transform(phi, 0, N, evec, 0, 0, psi, 0, num_bands);
transform(ctx_.processing_unit(), phi, 0, N, evec, 0, 0, psi, 0, num_bands);

/* exit the loop if the eigen-vectors are converged or this is a last iteration */
if (n <= itso.min_num_res_ || k == (itso.num_steps_ - 1)) {
break;
}
else { /* otherwise, set Psi as a new trial basis */
if (ctx_.control().verbosity_ > 2 && kp->comm().rank() == 0) {
if (ctx_.control().verbosity_ >= 3 && kp->comm().rank() == 0) {
DUMP("subspace size limit reached");
}

Expand Down Expand Up @@ -741,6 +745,17 @@ inline void Band::diag_pseudo_potential_davidson(K_point* kp__,
auto h_diag = get_h_diag(kp__, *local_op_, d_op__);
auto o_diag = get_o_diag(kp__, q_op__);

if (ctx_.control().print_checksum_) {
auto cs1 = h_diag.checksum();
auto cs2 = o_diag.checksum();
kp__->comm().allreduce(&cs1, 1);
kp__->comm().allreduce(&cs2, 1);
if (kp__->comm().rank() == 0) {
print_checksum("h_diag", cs1);
print_checksum("o_diag", cs2);
}
}

sddk::timer t3("sirius::Band::diag_pseudo_potential_davidson|iter");
for (int ispin_step = 0; ispin_step < num_spin_steps; ispin_step++) {

Expand Down Expand Up @@ -775,6 +790,16 @@ inline void Band::diag_pseudo_potential_davidson(K_point* kp__,
* N is the number of previous basis functions
* n is the number of new basis functions */
set_subspace_mtrx(num_sc, N, n, phi, hphi, hmlt, hmlt_old);

//== static int counter{0};
//== std::stringstream s;
//== if (ctx_.processing_unit() == CPU) {
//== s<<"hmlt_cpu"<<counter;
//== } else {
//== s<<"hmlt_gpu"<<counter;
//== }
//== hmlt.serialize(s.str(), N + n);
//== counter++;

if (ctx_.control().verification_ >= 1) {
double max_diff = Utils::check_hermitian(hmlt, N + n);
Expand Down Expand Up @@ -830,7 +855,7 @@ inline void Band::diag_pseudo_potential_davidson(K_point* kp__,

if (ctx_.control().verbosity_ >= 2 && kp__->comm().rank() == 0) {
DUMP("step: %i, current subspace size: %i, maximum subspace size: %i", k, N, num_phi);
if (ctx_.control().verbosity_ >= 3) {
if (ctx_.control().verbosity_ >= 4) {
for (int i = 0; i < num_bands; i++) {
DUMP("eval[%i]=%20.16f, diff=%20.16f, occ=%20.16f", i, eval[i], std::abs(eval[i] - eval_old[i]),
kp__->band_occupancy(i + ispin_step * ctx_.num_fv_states()));
Expand All @@ -851,16 +876,16 @@ inline void Band::diag_pseudo_potential_davidson(K_point* kp__,
/* recompute wave-functions */
/* \Psi_{i} = \sum_{mu} \phi_{mu} * Z_{mu, i} */
if (nc_mag) {
transform<T>(1.0, {&phi}, 0, N, evec, 0, 0, 0.0, {&psi}, 0, num_bands);
transform<T>(ctx_.processing_unit(), 1.0, {&phi}, 0, N, evec, 0, 0, 0.0, {&psi}, 0, num_bands);
} else {
transform<T>(phi.component(0), 0, N, evec, 0, 0, psi.component(ispin_step), 0, num_bands);
transform<T>(ctx_.processing_unit(), phi.component(0), 0, N, evec, 0, 0, psi.component(ispin_step), 0, num_bands);
}

/* exit the loop if the eigen-vectors are converged or this is a last iteration */
if (n <= itso.min_num_res_ || k == (itso.num_steps_ - 1)) {
break;
} else { /* otherwise, set Psi as a new trial basis */
if (ctx_.control().verbosity_ > 2 && kp__->comm().rank() == 0) {
if (ctx_.control().verbosity_ >= 3 && kp__->comm().rank() == 0) {
DUMP("subspace size limit reached");
}
hmlt_old.zero();
Expand All @@ -876,7 +901,7 @@ inline void Band::diag_pseudo_potential_davidson(K_point* kp__,

/* need to compute all hpsi and opsi states (not only unconverged) */
if (converge_by_energy) {
transform<T>(1.0, std::vector<Wave_functions*>({&hphi, &ophi}), 0, N, evec, 0, 0, 0.0, {&hpsi, &opsi}, 0, num_bands);
transform<T>(ctx_.processing_unit(), 1.0, std::vector<Wave_functions*>({&hphi, &ophi}), 0, N, evec, 0, 0, 0.0, {&hpsi, &opsi}, 0, num_bands);
}

/* update basis functions, hphi and ophi */
Expand Down
4 changes: 1 addition & 3 deletions src/Band/get_h_o_diag.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -218,11 +218,9 @@ Band::get_o_diag(K_point* kp__,
}
}
}
#ifdef __GPU
if (ctx_.processing_unit() == GPU) {
o_diag.allocate(memory_t::device);
o_diag.copy_to_device();
o_diag.copy<memory_t::host, memory_t::device>();
}
#endif
return std::move(o_diag);
}
Loading

0 comments on commit f1df939

Please sign in to comment.