1use nalgebra::DMatrix;
7
8#[derive(Debug, Clone)]
10pub struct ScfOutput {
11 pub energy: f64,
13 pub orbital_energies: Vec<f64>,
15 pub iterations: usize,
17 pub converged: bool,
19 pub density: Option<DMatrix<f64>>,
21 pub mulliken_charges: Option<Vec<f64>>,
23}
24
25#[derive(Debug, Clone)]
27pub struct ScfConvergenceConfig {
28 pub max_iter: usize,
30 pub energy_threshold: f64,
32 pub density_threshold: f64,
34 pub diis_size: usize,
36 pub level_shift: f64,
38 pub use_adiis: bool,
40 pub adiis_switch_iter: usize,
42}
43
44impl Default for ScfConvergenceConfig {
45 fn default() -> Self {
46 ScfConvergenceConfig {
47 max_iter: 100,
48 energy_threshold: 1e-8,
49 density_threshold: 1e-6,
50 diis_size: 8,
51 level_shift: 0.0,
52 use_adiis: false,
53 adiis_switch_iter: 5,
54 }
55 }
56}
57
58pub trait ScfSolver {
60 fn solve(&self, config: &ScfConvergenceConfig) -> Result<ScfOutput, String>;
62
63 fn method_name(&self) -> &str;
65
66 fn n_basis(&self) -> usize;
68
69 fn n_electrons(&self) -> usize;
71}
72
73pub struct HfScfSolver {
75 pub h_core: DMatrix<f64>,
76 pub s_mat: DMatrix<f64>,
77 pub eris: Vec<f64>,
78 pub n_elec: usize,
79}
80
81impl ScfSolver for HfScfSolver {
82 fn solve(&self, config: &ScfConvergenceConfig) -> Result<ScfOutput, String> {
83 let hf_config = super::scf::ScfConfig {
84 max_iter: config.max_iter,
85 energy_threshold: config.energy_threshold,
86 density_threshold: config.density_threshold,
87 diis_size: config.diis_size,
88 level_shift: config.level_shift,
89 };
90
91 let result = super::scf::solve_scf(
92 &self.h_core,
93 &self.s_mat,
94 &self.eris,
95 self.n_elec,
96 &hf_config,
97 );
98
99 Ok(ScfOutput {
100 energy: result.energy,
101 orbital_energies: result.orbital_energies,
102 iterations: result.iterations,
103 converged: result.converged,
104 density: Some(result.density),
105 mulliken_charges: None,
106 })
107 }
108
109 fn method_name(&self) -> &str {
110 "HF"
111 }
112 fn n_basis(&self) -> usize {
113 self.h_core.nrows()
114 }
115 fn n_electrons(&self) -> usize {
116 self.n_elec
117 }
118}
119
120pub struct Pm3ScfSolver {
122 pub elements: Vec<u8>,
123 pub positions: Vec<[f64; 3]>,
124}
125
126impl ScfSolver for Pm3ScfSolver {
127 fn solve(&self, _config: &ScfConvergenceConfig) -> Result<ScfOutput, String> {
128 let result = crate::pm3::solver::solve_pm3(&self.elements, &self.positions)?;
129 Ok(ScfOutput {
130 energy: result.total_energy,
131 orbital_energies: result.orbital_energies,
132 iterations: result.scf_iterations,
133 converged: result.converged,
134 density: None,
135 mulliken_charges: Some(result.mulliken_charges),
136 })
137 }
138
139 fn method_name(&self) -> &str {
140 "PM3"
141 }
142 fn n_basis(&self) -> usize {
143 0
144 } fn n_electrons(&self) -> usize {
146 0
147 }
148}
149
150pub struct XtbScfSolver {
152 pub elements: Vec<u8>,
153 pub positions: Vec<[f64; 3]>,
154}
155
156impl ScfSolver for XtbScfSolver {
157 fn solve(&self, _config: &ScfConvergenceConfig) -> Result<ScfOutput, String> {
158 let result = crate::xtb::solver::solve_xtb(&self.elements, &self.positions)?;
159 Ok(ScfOutput {
160 energy: result.total_energy,
161 orbital_energies: result.orbital_energies,
162 iterations: result.scc_iterations,
163 converged: result.converged,
164 density: None,
165 mulliken_charges: Some(result.mulliken_charges),
166 })
167 }
168
169 fn method_name(&self) -> &str {
170 "xTB"
171 }
172 fn n_basis(&self) -> usize {
173 0
174 }
175 fn n_electrons(&self) -> usize {
176 0
177 }
178}