Testing ROM MVP MD (#338)

* Fix EnergyAndForces tests (#277)

* have them work in debug mode too

* Move factor 4pi out og linear solvers (#278)

* Move some code into PoissonSolverFactory (#279)

* Clean up class Potentials (#280)

* Clean up class Ions, add test for it (#281)

* Add test MD_MVP (#290)

* Clean up code related to DM restart data (#292)

* Write dm (#291)

* Update use of DM in restart

* Remove unused function in Control (#294)

* Change symlink to restart in tests (#295)

* enable tests when old link present

* Extract number empty orbitals from restart file (#296)

* Clean up MD_IonicStepper restart data write (#297)

* Clean up code based on compiler warnings (#299)

* Add getForces for evaluateDMandEnergyAndForces (#300)

* Add check for compatibility MVP/Mehrstellen (#301)

* Add support for Br atom (#302)

* Clean up some HDFrestart functions (#303)

* Fix and test restart single hdf5 file (#305)

* Save Hartree potential for write in restart file (#306)

* Speed-up recently added tests (#309)

* simply use a smaller domain/mesh

* Encapsulate some functions in Potentials (#310)

* some functions were called at random places and confusing

* Remove confusing 0 in naming restart files (#308)

* use added integer only in case of fail/retry

* Add functionalities for extra info in restart file (#312)

* MPI abort (#313)

* Fix error code used by MPI_Abort

* Rho and VH restart  (#311)

* enable restart with consistent rho and VHartree

---------

Co-authored-by: Seung Whan Chung <[email protected]>

* Strenghten testIons (#315)

* More clean up in class Potentials (#316)

* remove some unused functions
* make use of class MGmol_MPI
* change some function names for clarity
* fix some minor bug in iterative indexes

* Misc code fixes (#317)

* exposed by porting on new environment

* Use unique restart filenames in test HDF5single (#318)

* Swap ions in existing test (#319)

* strengthen testing by swapping ions in testRestartEnergyAndForces

* Code clean up (#320)

* use more const
* use initNuc() instead of moveVnuc()
* remove incorrect assert
* use MGmol_MPI in more cases
* setup g_kbpsi_ inside initNuc()
* rename initNuc() into setupPotentials()

* Added functionalities to set local forces (#321)

* Update mixed precision code (#322)

* enable cmake build with mixed precision
* enable use of BLIS
* isolate blas3 loop rewrite in separate files
* fix misc issues with mixed precision code

* Clean up and fixes Ions (#326)

* reset static index counter ion Ion when needed
* simplify some loops
* make Ions copy constructor private

* Fixes for build without HDF5P (#324)

* Fix testRhoVhRestart (#325)

* Fix a few more issues with class Ions (#328)

* Introduce new constructor for class Ions (#330)

* modify various functions to use argument of type Ions instead of MGmol::ions_

* Add cleanup and MGmolInterface (#331)

* Add function to set local forces (#333)

* based on matching coordinates to local ions

* Fix setLocalForces (#335)

* Fix test WFEnergyAndForces (#336)

* writing restart twice in same directory was failing

* Add mixing option in MVP (#337)

* Timer

* Update MVP with mixing (#339)

* use tolerance to terminate iterations
* tune verbosity

* Updates

* use setupPotentials

* Setup potentials again

* setupRhoc is done within setupPotentials

* Avoid setPositions

* Rescaling ROM ion density

* Relax threshold

* Increase tolerance for testROMIonDensity

---------

Co-authored-by: Jean-Luc Fattebert <[email protected]>
Co-authored-by: Seung Whan Chung <[email protected]>

Author: 3 people

examples/PinnedH2O/job.basis_1_50

@@ -10,7 +10,7 @@ setenv OMP_NUM_THREADS 1
 
 set ncpus = 8
 
-set maindir = /p/lustre2/cheung26/mgmol-20250219
+set maindir = /p/lustre2/cheung26/mgmol
 
 setenv LD_LIBRARY_PATH ${maindir}/build_quartz/libROM/build/lib:$LD_LIBRARY_PATH
 

examples/PinnedH2O/job.offline

@@ -10,7 +10,7 @@ setenv OMP_NUM_THREADS 1
 
 set ncpus = 8
 
-set maindir = /p/lustre2/cheung26/mgmol-20250219
+set maindir = /p/lustre2/cheung26/mgmol
 
 setenv LD_LIBRARY_PATH ${maindir}/build_quartz/libROM/build/lib:$LD_LIBRARY_PATH
 

examples/PinnedH2O/job.ref

@@ -11,7 +11,7 @@ setenv OMP_NUM_THREADS 1
 set ncpus = 8
 set case = 2
 
-set maindir = /p/lustre2/cheung26/mgmol-20250219
+set maindir = /p/lustre2/cheung26/mgmol
 
 setenv LD_LIBRARY_PATH ${maindir}/build_quartz/libROM/build/lib:$LD_LIBRARY_PATH
 

examples/PinnedH2O/job.rom

@@ -10,7 +10,7 @@ setenv OMP_NUM_THREADS 1
 
 set ncpus = 8
 
-set maindir = /p/lustre2/cheung26/mgmol-20250219
+set maindir = /p/lustre2/cheung26/mgmol
 
 setenv LD_LIBRARY_PATH ${maindir}/build_quartz/libROM/build/lib:$LD_LIBRARY_PATH
 

examples/PinnedH2O/job.rom_3DOF

@@ -11,7 +11,7 @@ setenv OMP_NUM_THREADS 1
 set ncpus = 8
 set case = 2
 
-set maindir = /p/lustre2/cheung26/mgmol-20250219
+set maindir = /p/lustre2/cheung26/mgmol
 
 setenv LD_LIBRARY_PATH ${maindir}/build_quartz/libROM/build/lib:$LD_LIBRARY_PATH
 

examples/PinnedH2O/mgmol_rom_3DOF_test1.cfg

@@ -1,4 +1,4 @@
-verbosity=1
+verbosity=2
 xcFunctional=PBE
 FDtype=4th
 [Mesh]
@@ -34,7 +34,9 @@ initial_width=2.
 nempty=30
 [DensityMatrix]
 solver=MVP
-nb_inner_it=100
+nb_inner_it=20
+mixing=0.8
+tol=1.e-8
 [Restart]
 output_level=4
 [ROM]

examples/PinnedH2O/mgmol_rom_3DOF_test2.cfg

@@ -1,4 +1,4 @@
-verbosity=1
+verbosity=2
 xcFunctional=PBE
 FDtype=4th
 [Mesh]
@@ -30,7 +30,9 @@ initial_width=2.
 nempty=30
 [DensityMatrix]
 solver=MVP
-nb_inner_it=100
+nb_inner_it=20
+mixing=0.8
+tol=1.e-8
 [Restart]
 output_level=4
 [ROM]

src/Control.cc

@@ -224,6 +224,7 @@ void Control::print(std::ostream& os)
        << conv_tol << std::endl;
     os << std::fixed;
     os << " Density matrix mixing = " << dm_mix << std::endl;
+    os << std::setprecision(4) << std::scientific << " Density matrix tol = " << dm_tol << std::endl;
     if (DMEigensolver() == DMEigensolverType::Eigensolver)
     {
         os << " Density matrix computation algorithm = "
@@ -439,7 +440,7 @@ void Control::sync(void)
         memset(&int_buffer[0], 0, size_int_buffer * sizeof(int));
     }
 
-    const short size_float_buffer = 43;
+    const short size_float_buffer = 44;
     float* float_buffer           = new float[size_float_buffer];
     if (mype_ == 0)
     {
@@ -485,6 +486,7 @@ void Control::sync(void)
         float_buffer[40] = threshold_eigenvalue_gram_quench_;
         float_buffer[41] = pair_mlwf_distance_threshold_;
         float_buffer[42] = e0_;
+        float_buffer[43] = dm_tol;
     }
     else
     {
@@ -680,6 +682,7 @@ void Control::sync(void)
     threshold_eigenvalue_gram_quench_ = float_buffer[40];
     pair_mlwf_distance_threshold_     = float_buffer[41];
     e0_                               = float_buffer[42];
+    dm_tol                            = float_buffer[43];
     max_electronic_steps_loose_       = max_electronic_steps;
 
     delete[] short_buffer;
@@ -699,8 +702,8 @@ void Control::setDefaultValues()
 
 void Control::adjust()
 {
-    // change dm_mix default to 1. if not using Davidson
-    if (it_algo_type_ != 2 && dm_mix < 0.) dm_mix = 1.;
+    // change dm_mix default to 1. if not using Davidson or MVP
+    if ((it_algo_type_ != 2 && DM_solver_ != 1) && dm_mix < 0.) dm_mix = 1.;
 
     if (nel_ - 2 * numst == 0)
     {
@@ -1720,6 +1723,7 @@ void Control::setOptions(const boost::program_options::variables_map& vm)
             lrs_extrapolation = 10;
 
         dm_mix         = vm["DensityMatrix.mixing"].as<float>();
+        dm_tol         = vm["DensityMatrix.tol"].as<float>();
         dm_inner_steps = vm["DensityMatrix.nb_inner_it"].as<short>();
         dm_use_old_    = vm["DensityMatrix.use_old"].as<bool>() ? 1 : 0;
         str            = vm["DensityMatrix.algo"].as<std::string>();
@@ -1739,8 +1743,6 @@ void Control::setOptions(const boost::program_options::variables_map& vm)
         else
             dm_algo_ = 2;
 
-        dm_tol = vm["DensityMatrix.tol"].as<float>();
-
         str = vm["DensityMatrix.solver"].as<std::string>();
         if (str.compare("Mixing") == 0) DM_solver_ = 0;
         if (str.compare("MVP") == 0) DM_solver_ = 1;

src/Ions.cc

@@ -1315,22 +1315,6 @@ void Ions::setLocalForces(
         {
             if (ion->compareName(*s))
             {
-                //std::cout << "Ion found: " << ion->name() << std::endl;
-                //std::cout << "Ion force: (" << *it << ", " << *(it + 1) << ", " << *(it + 2) << ")" << std::endl;
-                //std::cout << "names: ";
-                //for (int i = 0; i < names.size(); i++)
-                //{
-                //    std::cout << names[i];
-                //    if (i == forces.size() - 1) std::cout << std::endl;
-                //    else std::cout << ", ";
-                //}
-                //std::cout << "forces: ";
-                //for (int i = 0; i < forces.size(); i++)
-                //{
-                //    std::cout << forces[i];
-                //    if (i == forces.size() - 1) std::cout << ")" << std::endl;
-                //    else std::cout << ", ";
-                //}
                 ion->set_force(0, *it);
                 ion->set_force(1, *(it + 1));
                 ion->set_force(2, *(it + 2));
@@ -1364,9 +1348,6 @@ void Ions::setLocalForces(
             double d = std::sqrt(d2);
             if (d < tol)
             {
-                std::cout << "Ion found: " << ion->name() << std::endl;
-                std::cout << "Ion position:( " << *cit << ", " << *(cit + 1) << ", " << *(cit + 2) << ")" << std::endl;
-                std::cout << "Ion force: (" << *fit << ", " << *(fit + 1) << ", " << *(fit + 2) << ")" << std::endl;
                 ion->set_force(0, *fit);
                 ion->set_force(1, *(fit + 1));
                 ion->set_force(2, *(fit + 2));

src/MGmol.cc

@@ -102,6 +102,7 @@ extern Timer md_moveVnuc_tm;
 extern Timer md_updateMasks_tm;
 extern Timer md_extrapolateOrbitals_tm;
 extern Timer md_updateRhoAndPot_tm;
+extern Timer md_updateDMandEnergy_tm;
 extern Timer quench_tm;
 extern Timer ions_setupInteractingIons_tm;
 extern Timer ions_setup_tm;
@@ -928,6 +929,7 @@ void MGmol<OrbitalsType>::printTimers()
     init_nuc_tm_.print(os_);
     md_updateMasks_tm.print(os_);
     md_extrapolateOrbitals_tm.print(os_);
+    md_updateDMandEnergy_tm.print(os_);
     quench_tm.print(os_);
     evnl_tm_.print(os_);
     ions_setupInteractingIons_tm.print(os_);