SCF diagnostics use stale `phi` scratch data after `multilevel_solve_for_new_phi()`

[WeiqunZhang]

Summary

do_hscf_solve() copies PhiGrav_Type into local scratch phi[lev] early in each iteration, then later calls gravity->multilevel_solve_for_new_phi(...), but the diagnostics still read from the old scratch phi[lev] instead of refreshed level data.

Location

• Source/scf/scf_relax.cpp:203
• Source/scf/scf_relax.cpp:605
• Source/scf/scf_relax.cpp:639

Problem Details

The sequence in one SCF iteration is:

1. Copy get_new_data(PhiGrav_Type) into phi[lev].
2. Modify density and solve gravity with multilevel_solve_for_new_phi.
3. Compute pot_eng from phi_arr = (*phi[lev])[mfi].array();.

So the printed potential energy and virial diagnostic are based on pre-solve phi, not the just-updated gravitational potential.

Impact

• Reported potential energy can lag the actual solved state.
• Virial error output can be misleading during relaxation monitoring.
• SCF convergence diagnostics become harder to interpret.

Suggested Patch

Use current level data in the diagnostic reduce loop (or refresh phi[lev] from PhiGrav_Type immediately after the gravity solve).

diff --git a/Source/scf/scf_relax.cpp b/Source/scf/scf_relax.cpp
--- a/Source/scf/scf_relax.cpp
+++ b/Source/scf/scf_relax.cpp
@@
         gravity->multilevel_solve_for_new_phi(0, finest_level);

+        // Refresh local phi scratch with the newly solved potential.
+        for (int lev = 0; lev <= finest_level; ++lev) {
+            MultiFab::Copy((*phi[lev]), getLevel(lev).get_new_data(PhiGrav_Type), 0, 0, 1, 0);
+            if (lev < finest_level) {
+                const MultiFab& mask = getLevel(lev+1).build_fine_mask();
+                MultiFab::Multiply((*phi[lev]), mask, 0, 0, 1, 0);
+            }
+        }
+
         // Update diagnostic quantities.

Prepared by Codex