Feature: Support calculation of multi-k point band-decomposed charge density under LCAO basis

[AsTonyshment]

Linked Issue

Fix #2905, fix #3572, fix #4412, fix #4697.

What's changed?

• Added support for calculating band-decomposed charge density for multiple $\boldsymbol{k}$-points.
• Refactored istate_charge.cpp by extracting large segments of reused code into separate functions.
• Fixed the issue where charge density calculation was incorrect when nspin=2.
• Corrected the naming of output .cube files for the spin-down component when nspin=2.
• Standardized the $\boldsymbol{k}$-point and spin indices in output file names to start from 1.
• Fixed the issue where bands_to_print was not working properly due to recent INPUT refactoring.
• Added a new INPUT parameter if_separate_k to specify whether to write the partial charge densities for all $\boldsymbol{k}$-points to individual files or merge them.
• Updated the online documentation for ABACUS.

Any changes of core modules? (ignore if not applicable)

• Added a new DensityMatrix function to obtain DMR for certain $\boldsymbol{k}$-point.

[AsTonyshment]

@Qianruipku Maybe you can have a look at whether I used the new INPUT formalism correctly? I added a system parameter out_band_kb, which is obtained by parsing the direct input bands_to_print.

[mohanchen]

Why the docs/readthedocs.org:abacus-rtd has stopped? I think it doesn't related to the codes.