cannot recompute rho after loading a restart file.

[dreamer2368]

For debugging purpose, I'm trying to check if I'm properly loading the density from a restart file. Following snippet attempts to load a restart file, recompute the density and compare it from the density in the restart file:

    MGmol_MPI& mmpi      = *(MGmol_MPI::instance());
    const int rank = mmpi.mypeGlobal();
    const int nprocs = mmpi.size();

    Control& ct = *(Control::instance());
    MGmol<OrbitalsType> *mgmol = static_cast<MGmol<OrbitalsType> *>(mgmol_);
    std::shared_ptr<Rho<OrbitalsType>> rho = mgmol->getRho();

    /* load a restart file */
    mgmol->loadRestartFile(ct.restart_file);

    /* save density from the restart file to elsewhere */
    std::vector<RHODTYPE> rho0(rho->rho_[0].size());
    rho0 = rho->rho_[0];

    /* recompute rho from the orbital */
    rho->update(*mgmol->getOrbitals());

    /* check if the recomputed density is the same */
    for (int d = 0; d < rho0.size(); d++)
    {
        if (rho0[d] != rho->rho_[0][d])
            printf("rank %d, rho[%d]=%.15e, rho0[%d]=%.15e\n", rank, d, rho->rho_[0][d], d, rho0[d]);
    }

This is also implemented in rho-bug branch, as an executable test_rho_restart. This does not need libROM library, same installation as release branch. It can be run just as mgmol_opt:

test_rho_restart -c carbyne.cfg

However, when I try this on any restart file of Carbyne example, I receive the following error:

 Rho<OrbitalsType>::rescaleTotalCharge(), charge = 1.96e+02
 Rescaling factor: 2.55e-01
 Num. electrons: 5.00e+01
Error on density too large to continue. Abort.

@jeanlucf22 , can you take a look at this minimal example and check if I'm doing something wrong here? This test is necessary before I get into debugging other parts of ROM. I can send you a config and restart file as well if needed.

[dreamer2368]

After PR #291, #292, #293, ROMFPMD and pot-rom branches are merged/rebased.

PR #291 - #293 did fix the inconsistency issue between density matrix and electronic density in the same restart file. Test on Poisson ROM solver verified the ROM for electronic density up to machine precision.

However, the resulting potential is still inconsistent, having about 6~7% error on a reproductive case. Just as for electronic density, it is expected to have a machine-precision error (<1e-10).

I rebased rho-bug branch from the up-to-date release. A new test routine is added to verify the consistency of potential with rho/rhoc in the same restart file. When the potential is computed with rho/rhoc from a restart file, it was different from the potential from the same restart file (having ~5% error).

@jeanlucf22 , can you look into this branch rho-bug, the same executable test_rho_restart? Perhaps something is missing when a restart file is loaded?

The executable can be run just as the main mgmol driver, with a config file specifying a restart file:

srun -n 144 test_rho_restart -c carbyne.cfg

[jeanlucf22]

@dreamer2368 I cannot reproduce the error. Actually, as is, your test code fails already in testing Rho to to a very tight tolerance. But if I loosen that tolerance, the error I see is very small.
Note that for the potential, there may be a small error due to the fact that the Hartree potential is computed iteratively. But not 5%.

[dreamer2368]

@jeanlucf22 , can you confirm that the test fails with a newly saved restart file (that has Density_Matrix_WF dataset)?

If the test_rho_restart still fails, it means PR #291 still doesn't fix the issue. With this inconsistency, we cannot discern whether an error comes from this inconsistency, ROM approximation, or a bug.

[jeanlucf22]

Actually, I was using a restart in single precision. I now switched to double and the test on rho passes.
For the test on the Hartree potential, I get an error

FOM potential relative error: 2.289e-05

which may be explained by the iterative solver not converged to machine precision.
I am running your test, but with my own input. If you share your input, I could try to reproduce the issue.