1use super::{ElectronicStopping, Ion, StoppingError};
15use crate::material::Material;
16
17pub trait CompoundCorrection {
19 fn factor(&self, material: &Material, ion: &Ion, energy_ev: f64) -> f64;
23}
24
25#[derive(Debug, Clone, Copy, Default)]
27pub struct NoCorrection;
28
29impl CompoundCorrection for NoCorrection {
30 fn factor(&self, _: &Material, _: &Ion, _: f64) -> f64 {
31 1.0
32 }
33}
34
35#[derive(Debug, Clone, Copy)]
37pub struct ConstantCorrection(pub f64);
38
39impl CompoundCorrection for ConstantCorrection {
40 fn factor(&self, _: &Material, _: &Ion, _: f64) -> f64 {
41 self.0
42 }
43}
44
45pub fn bragg_cross_section_per_atom(
48 model: &dyn ElectronicStopping,
49 correction: &dyn CompoundCorrection,
50 ion: &Ion,
51 material: &Material,
52 energy_ev: f64,
53) -> Result<f64, StoppingError> {
54 let mut sum = 0.0;
55 for c in material.components() {
56 sum += c.atom_fraction() * model.stopping(ion, c.z(), energy_ev)?;
57 }
58 let factor = correction.factor(material, ion, energy_ev);
59 if !(factor.is_finite() && factor >= 0.0) {
60 return Err(StoppingError::InvalidParameter {
61 name: "compound_correction",
62 value: factor,
63 });
64 }
65 Ok(factor * sum)
66}
67
68pub fn bragg_stopping_power(
71 model: &dyn ElectronicStopping,
72 correction: &dyn CompoundCorrection,
73 ion: &Ion,
74 material: &Material,
75 energy_ev: f64,
76) -> Result<f64, StoppingError> {
77 Ok(material.atom_number_density()
78 * bragg_cross_section_per_atom(model, correction, ion, material, energy_ev)?)
79}
80
81#[cfg(test)]
82mod tests {
83 use super::super::lindhard_scharff::LindhardScharff;
84 use super::*;
85
86 #[test]
87 fn single_element_equals_element_stopping_exactly() {
88 let m = LindhardScharff::new();
89 let ion = Ion::new(15).unwrap();
90 let si = Material::from_atom_fractions(&[(14, 1.0)], None).unwrap();
91 let want = m.stopping(&ion, 14, 3.0e4).unwrap();
92 let got = bragg_cross_section_per_atom(&m, &NoCorrection, &ion, &si, 3.0e4).unwrap();
93 assert_eq!(got.to_bits(), want.to_bits());
94 }
95
96 #[test]
97 fn sio2_is_hand_weighted_sum() {
98 let m = LindhardScharff::new();
99 let ion = Ion::new(5).unwrap();
100 let sio2 = Material::from_atom_fractions(&[(14, 1.0), (8, 2.0)], Some(2200.0)).unwrap();
101 let e = 2.0e4;
102 let (s_si, s_o) = (
103 m.stopping(&ion, 14, e).unwrap(),
104 m.stopping(&ion, 8, e).unwrap(),
105 );
106 let hand = (s_si + 2.0 * s_o) / 3.0;
107 let got = bragg_cross_section_per_atom(&m, &NoCorrection, &ion, &sio2, e).unwrap();
108 assert!((got / hand - 1.0).abs() < 1e-14);
109 let p = bragg_stopping_power(&m, &NoCorrection, &ion, &sio2, e).unwrap();
110 assert!((p / (sio2.atom_number_density() * hand) - 1.0).abs() < 1e-14);
111 }
112
113 #[test]
114 fn correction_hook_scales() {
115 let m = LindhardScharff::new();
116 let ion = Ion::new(5).unwrap();
117 let sio2 = Material::from_atom_fractions(&[(14, 1.0), (8, 2.0)], Some(2200.0)).unwrap();
118 let a = bragg_cross_section_per_atom(&m, &NoCorrection, &ion, &sio2, 1e4).unwrap();
119 let b =
120 bragg_cross_section_per_atom(&m, &ConstantCorrection(0.9), &ion, &sio2, 1e4).unwrap();
121 assert!((b / a - 0.9).abs() < 1e-14);
122 }
123
124 #[test]
125 fn rejects_bad_compound_correction() {
126 let m = LindhardScharff::new();
127 let ion = Ion::new(5).unwrap();
128 let sio2 = Material::from_atom_fractions(&[(14, 1.0), (8, 2.0)], Some(2200.0)).unwrap();
129 for bad in [-1.0, f64::NAN, f64::INFINITY] {
130 let r = bragg_cross_section_per_atom(&m, &ConstantCorrection(bad), &ion, &sio2, 1e4);
131 assert!(matches!(
132 r,
133 Err(StoppingError::InvalidParameter {
134 name: "compound_correction",
135 ..
136 })
137 ));
138 let r = bragg_stopping_power(&m, &ConstantCorrection(bad), &ion, &sio2, 1e4);
139 assert!(r.is_err());
140 }
141 let z = bragg_cross_section_per_atom(&m, &ConstantCorrection(0.0), &ion, &sio2, 1e4);
142 assert_eq!(z.unwrap(), 0.0);
143 }
144
145 #[test]
146 fn bad_element_correction_propagates_through_bragg() {
147 let m = LindhardScharff::new().with_correction(8, -1.0);
148 let ion = Ion::new(5).unwrap();
149 let sio2 = Material::from_atom_fractions(&[(14, 1.0), (8, 2.0)], Some(2200.0)).unwrap();
150 assert!(matches!(
151 bragg_cross_section_per_atom(&m, &NoCorrection, &ion, &sio2, 1e4),
152 Err(StoppingError::InvalidParameter { .. })
153 ));
154 }
155}