From a929bbdca458c7cbca52d0d6ba889544fbcfba44 Mon Sep 17 00:00:00 2001 From: Francois Gygi Date: Sun, 15 Jun 2008 20:44:58 +0000 Subject: [PATCH] Modified output for ground-state only calculations (niter=0) git-svn-id: http://qboxcode.org/svn/qb/trunk@627 cba15fb0-1239-40c8-b417-11db7ca47a34 --- src/BOSampleStepper.C | 240 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++-------------------------------------------------------------------------------------------------- 1 file changed, 142 insertions(+), 98 deletions(-) diff --git a/src/BOSampleStepper.C b/src/BOSampleStepper.C index a60ecd5..b8eb8b2 100644 --- a/src/BOSampleStepper.C +++ b/src/BOSampleStepper.C @@ -3,7 +3,7 @@ // BOSampleStepper.C // //////////////////////////////////////////////////////////////////////////////// -// $Id: BOSampleStepper.C,v 1.41 2008-04-15 01:36:15 fgygi Exp $ +// $Id: BOSampleStepper.C,v 1.42 2008-06-15 20:44:58 fgygi Exp $ #include "BOSampleStepper.h" #include "EnergyFunctional.h" @@ -98,6 +98,8 @@ void BOSampleStepper::step(int niter) const bool compute_stress = ( s_.ctrl.stress == "ON" ); const bool cell_moves = ( niter > 0 && compute_stress && cell_dyn != "LOCKED" ); + // GS-only calculation: + const bool gs_only = !atoms_move && !cell_moves; const bool use_confinement = ( s_.ctrl.ecuts > 0.0 ); Timer tm_iter; @@ -175,120 +177,123 @@ void BOSampleStepper::step(int niter) // compute energy and ionic forces using existing wavefunction - tmap["charge"].start(); - cd_.update_density(); - tmap["charge"].stop(); + if ( !gs_only ) + { + tmap["charge"].start(); + cd_.update_density(); + tmap["charge"].stop(); - ef_.update_vhxc(); - const bool compute_forces = true; - double energy = - ef_.energy(false,dwf,compute_forces,fion,compute_stress,sigma_eks); + ef_.update_vhxc(); + const bool compute_forces = true; + double energy = + ef_.energy(false,dwf,compute_forces,fion,compute_stress,sigma_eks); - if ( onpe0 ) - { - cout.setf(ios::fixed,ios::floatfield); - cout.setf(ios::right,ios::adjustfield); - cout << " " << setprecision(8) - << setw(15) << ef_.ekin() << " \n"; - if ( use_confinement ) - cout << " " << setw(15) << ef_.econf() << " \n"; - cout << " " << setw(15) << ef_.eps() << " \n" - << " " << setw(15) << ef_.enl() << " \n" - << " " << setw(15) << ef_.ecoul() << " \n" - << " " << setw(15) << ef_.exc() << " \n" - << " " << setw(15) << ef_.esr() << " \n" - << " " << setw(15) << ef_.eself() << " \n" - << " " << setw(15) << ef_.ets() << " \n" - << " " << setw(15) << ef_.etotal() << " \n"; - if ( compute_stress ) + if ( onpe0 ) { - const double pext = (sigma_ext[0]+sigma_ext[1]+sigma_ext[2])/3.0; - const double enthalpy = ef_.etotal() + pext * cell.volume(); - cout << " " << setw(15) << pext * cell.volume() - << " " << endl; - cout << " " << setw(15) << enthalpy << " \n" - << flush; + cout.setf(ios::fixed,ios::floatfield); + cout.setf(ios::right,ios::adjustfield); + cout << " " << setprecision(8) + << setw(15) << ef_.ekin() << " \n"; + if ( use_confinement ) + cout << " " << setw(15) << ef_.econf() << " \n"; + cout << " " << setw(15) << ef_.eps() << " \n" + << " " << setw(15) << ef_.enl() << " \n" + << " " << setw(15) << ef_.ecoul() << " \n" + << " " << setw(15) << ef_.exc() << " \n" + << " " << setw(15) << ef_.esr() << " \n" + << " " << setw(15) << ef_.eself() << " \n" + << " " << setw(15) << ef_.ets() << " \n" + << " " << setw(15) << ef_.etotal() << " \n"; + if ( compute_stress ) + { + const double pext = (sigma_ext[0]+sigma_ext[1]+sigma_ext[2])/3.0; + const double enthalpy = ef_.etotal() + pext * cell.volume(); + cout << " " << setw(15) << pext * cell.volume() + << " " << endl; + cout << " " << setw(15) << enthalpy << " \n" + << flush; + } } - } - if ( iter > 0 && ionic_stepper ) - { - ionic_stepper->compute_v(energy,fion); - } - // at this point, positions r0, velocities v0 and forces fion are - // consistent - double ekin_ion = 0.0, temp_ion = 0.0; - if ( ionic_stepper ) - { - ekin_ion = ionic_stepper->ekin(); - temp_ion = ionic_stepper->temp(); - } + if ( iter > 0 && ionic_stepper ) + { + ionic_stepper->compute_v(energy,fion); + } + // at this point, positions r0, velocities v0 and forces fion are + // consistent + double ekin_ion = 0.0, temp_ion = 0.0; + if ( ionic_stepper ) + { + ekin_ion = ionic_stepper->ekin(); + temp_ion = ionic_stepper->temp(); + } - // print positions, velocities and forces at time t0 - if ( onpe0 ) - { - cout << "" << endl; - cout << atoms.cell(); - for ( int is = 0; is < atoms.atom_list.size(); is++ ) + // print positions, velocities and forces at time t0 + if ( onpe0 ) { - int i = 0; - for ( int ia = 0; ia < atoms.atom_list[is].size(); ia++ ) + cout << "" << endl; + cout << atoms.cell(); + for ( int is = 0; is < atoms.atom_list.size(); is++ ) { - Atom* pa = atoms.atom_list[is][ia]; - cout << " name() << "\"" - << " species=\"" << pa->species() - << "\">\n" - << " " << pa->position() << " \n" - << " " << pa->velocity() << " \n" - << " " - << fion[is][i] << " " - << fion[is][i+1] << " " - << fion[is][i+2] - << " \n"; - cout << " " << endl; - i += 3; + int i = 0; + for ( int ia = 0; ia < atoms.atom_list[is].size(); ia++ ) + { + Atom* pa = atoms.atom_list[is][ia]; + cout << " name() << "\"" + << " species=\"" << pa->species() + << "\">\n" + << " " << pa->position() << " \n" + << " " << pa->velocity() << " \n" + << " " + << fion[is][i] << " " + << fion[is][i+1] << " " + << fion[is][i+2] + << " \n"; + cout << " " << endl; + i += 3; + } } + cout << "" << endl; + cout << " " << energy+ekin_ion << " \n"; + cout << " " << ekin_ion << " \n"; + cout << " " << temp_ion << " \n"; } - cout << "" << endl; - cout << " " << energy+ekin_ion << " \n"; - cout << " " << ekin_ion << " \n"; - cout << " " << temp_ion << " \n"; - } - if ( atoms_move ) - { - if ( s_.constraints.size() > 0 ) + if ( atoms_move ) { - s_.constraints.compute_forces(ionic_stepper->r0(), fion); - if ( onpe0 ) + if ( s_.constraints.size() > 0 ) { - s_.constraints.list_constraints(cout); + s_.constraints.compute_forces(ionic_stepper->r0(), fion); + if ( onpe0 ) + { + s_.constraints.list_constraints(cout); + } } + // move atoms to new position: r0 <- r0 + v0*dt + dt2/m * fion + ionic_stepper->compute_r(energy,fion); + ef_.atoms_moved(); } - // move atoms to new position: r0 <- r0 + v0*dt + dt2/m * fion - ionic_stepper->compute_r(energy,fion); - ef_.atoms_moved(); - } - if ( compute_stress ) - { - compute_sigma(); - print_stress(); - - if ( cell_moves ) + if ( compute_stress ) { - cell_stepper->compute_new_cell(sigma); + compute_sigma(); + print_stress(); + + if ( cell_moves ) + { + cell_stepper->compute_new_cell(sigma); - // Update cell - cell_stepper->update_cell(); + // Update cell + cell_stepper->update_cell(); - ef_.cell_moved(); - ef_.atoms_moved(); // modifications of the cell also move ions + ef_.cell_moved(); + ef_.atoms_moved(); // modifications of the cell also move ions - if ( use_preconditioner ) - preconditioner->update(); + if ( use_preconditioner ) + preconditioner->update(); + } } - } + } // if !gs_only // Recalculate ground state wavefunctions @@ -609,7 +614,7 @@ void BOSampleStepper::step(int niter) // subspace diagonalization if ( compute_eigvec || s_.ctrl.wf_diag == "EIGVAL" ) { - energy = ef_.energy(true,dwf,false,fion,false,sigma_eks); + ef_.energy(true,dwf,false,fion,false,sigma_eks); s_.wf.diag(dwf,compute_eigvec); if ( onpe0 ) { @@ -654,11 +659,50 @@ void BOSampleStepper::step(int niter) cout << " BOSampleStepper: end scf iteration" << endl; } // for itscf - // If GS calculation only, print energy at end of iterations - if ( !atoms_move && !cell_moves ) + // If GS calculation only, print energy and atomset at end of iterations + if ( gs_only ) + { + tmap["charge"].start(); + cd_.update_density(); + tmap["charge"].stop(); + + ef_.update_vhxc(); + const bool compute_forces = true; + ef_.energy(false,dwf,compute_forces,fion,compute_stress,sigma_eks); + if ( onpe0 ) + { cout << ef_; - + cout << "" << endl; + cout << atoms.cell(); + for ( int is = 0; is < atoms.atom_list.size(); is++ ) + { + int i = 0; + for ( int ia = 0; ia < atoms.atom_list[is].size(); ia++ ) + { + Atom* pa = atoms.atom_list[is][ia]; + cout << " name() << "\"" + << " species=\"" << pa->species() + << "\">\n" + << " " << pa->position() << " \n" + << " " << pa->velocity() << " \n" + << " " + << fion[is][i] << " " + << fion[is][i+1] << " " + << fion[is][i+2] + << " \n"; + cout << " " << endl; + i += 3; + } + } + cout << "" << endl; + if ( compute_stress ) + { + compute_sigma(); + print_stress(); + } + } + } wf_stepper->postprocess(); } else -- libgit2 0.26.0