[timchem1508]

Summary

The cache_container type is reconstructed into one multicharge cache type. The linear system solver type is added to speed-up the charge calculations and adpat the multicharge for futher improvements. Next to the existing direct solver (LAPACK Bunch–Kaufman), it adds a Jacobi-preconditioned Conjugate Gradient (CG) solver. For large systems, CG brings the cost of the linear solve down from O(N³) to O(N²).

Scaling

Application of the PCG solver vs. the direct one on the methane molecules sphere calculations:

The direct solver factorizes the dense (N+1)×(N+1) Coulomb matrix with sytrf/sytrs, which costs O(N³). Each CG iteration costs one dense matrix–vector product, which is O(N²). With the Jacobi preconditioner the number of iterations stays about the same as the system grows, so the whole solve scales as O(N²). Use the (Nid)×(Nid) van-der-Waals radii matrix instead of the (N)×(N) one to save time and memory.

• At 5,000 atoms the direct solve takes 0.79 s and CG takes about 0.07 s, roughly 12× faster.
• CG solves 12,500 atoms in 0.40 s, about what the direct solver needs for 4,000 atoms. Extending the cubic fit puts the direct solver at about 12 s for 12,500 atoms.

Changes

Model and cache

• mchrg_model_type separates setup from solution. update() fills mchrg_cache with the CN, local charges, capacitance and Coulomb matrices, and their derivatives. solve() takes a solver and returns the charges, energy, gradient and virial, and the charge derivatives.
• mchrg_cache replaces cache_container.

Solvers (multicharge_solver)

• Adds the abstract types mchrg_solver_type and mchrg_solver_input. new_mchrg_solver builds a solver from its input.
• direct_solver: Bunch–Kaufman factorization of the system that includes the Lagrange multiplier for the total-charge constraint. It computes the inverse (sytri) only when charge derivatives are requested.
• cg_solver: Jacobi-preconditioned CG. The system with the Lagrange multiplier is indefinite, so CG cannot solve it directly. Instead, CG solves two positive-definite N×N systems, A u = 1 and A v = −x, and gets the multiplier from the total-charge constraint. Defaults are 1000 iterations and a relative residual tolerance of 1e-15.
• CG returns an error when it does not converge within the iteration limit, and when it breaks down because the matrix is not positive definite (non-positive curvature).

Other

• Generic dot, scal, axpy and symv wrappers in multicharge_blas.
• New command line options: --solver cg|direct, --maxiter, --tolerance, -qg/--qgrad, -v/--verbose and -s/--silent. The man page is updated.
• Fix the CLI multiple charges error.
• Use the elementwise rvdw instead of the atomwise to save time and memory.
• The get_charges and add_charge_hessian APIs uses the direct solver at the moment.

Behavior changes and limitations

• The charge derivatives dqdr and dqdL need the inverse matrix, so only the direct solver provides them. Requesting them with CG returns an error. The energy gradient and virial work with both solvers.
• -v now means --verbose. Use --version to print the version.
• API changes: cache_container is now mchrg_cache, and solve takes a solver argument.

Tests

• New test_solver suite for CG: identity, diagonal, SPD, ill-conditioned and random SPD matrices, zero right-hand side, invalid input, non-convergence, breakdown and verbosity levels. For the direct solver: inverse, singular matrix and verbosity.

[codecov]

Codecov Report

❌ Patch coverage is 23.07692% with 90 lines in your changes missing coverage. Please review.
✅ Project coverage is 24.81%. Comparing base (282626e) to head (22fd268).
⚠️ Report is 30 commits behind head on main.

Files with missing lines	Patch %	Lines
app/main.f90	23.07%	72 Missing and 18 partials ⚠️

Additional details and impacted files

@@             Coverage Diff             @@
##             main      #77       +/-   ##
===========================================
- Coverage   51.57%   24.81%   -26.77%     
===========================================
  Files          18        2       -16     
  Lines        1648      266     -1382     
  Branches      746      116      -630     
===========================================
- Hits          850       66      -784     
+ Misses        390      153      -237     
+ Partials      408       47      -361     

☔ View full report in Codecov by Harness.
📢 Have feedback on the report? Share it here.

🚀 New features to boost your workflow:

• ❄️ Test Analytics: Detect flaky tests, report on failures, and find test suite problems.