Skip to main content

sciforge_lib/physics/electrodynamics/
fields.rs

1use crate::constants::{EPSILON_0, K_COULOMB, MU_0};
2
3pub fn electric_field_point_charge(q: f64, r: [f64; 3]) -> [f64; 3] {
4    let r_mag = (r[0] * r[0] + r[1] * r[1] + r[2] * r[2]).sqrt();
5    if r_mag < 1e-30 {
6        return [0.0; 3];
7    }
8    let e_mag = K_COULOMB * q / (r_mag * r_mag);
9    [
10        e_mag * r[0] / r_mag,
11        e_mag * r[1] / r_mag,
12        e_mag * r[2] / r_mag,
13    ]
14}
15
16pub fn electric_potential_point(q: f64, r: f64) -> f64 {
17    K_COULOMB * q / r
18}
19
20pub fn magnetic_field_wire(current: f64, r: f64) -> f64 {
21    MU_0 * current / (2.0 * std::f64::consts::PI * r)
22}
23
24pub fn magnetic_field_solenoid(n_per_length: f64, current: f64) -> f64 {
25    MU_0 * n_per_length * current
26}
27
28pub fn magnetic_field_loop(current: f64, radius: f64, z: f64) -> f64 {
29    MU_0 * current * radius * radius / (2.0 * (radius * radius + z * z).powf(1.5))
30}
31
32pub fn biot_savart_segment(current: f64, dl: [f64; 3], r: [f64; 3]) -> [f64; 3] {
33    let r_mag = (r[0] * r[0] + r[1] * r[1] + r[2] * r[2]).sqrt();
34    if r_mag < 1e-30 {
35        return [0.0; 3];
36    }
37    let cross = [
38        dl[1] * r[2] - dl[2] * r[1],
39        dl[2] * r[0] - dl[0] * r[2],
40        dl[0] * r[1] - dl[1] * r[0],
41    ];
42    let factor = MU_0 * current / (4.0 * std::f64::consts::PI * r_mag.powi(3));
43    [factor * cross[0], factor * cross[1], factor * cross[2]]
44}
45
46pub fn lorentz_force(q: f64, v: [f64; 3], e: [f64; 3], b: [f64; 3]) -> [f64; 3] {
47    let vxb = [
48        v[1] * b[2] - v[2] * b[1],
49        v[2] * b[0] - v[0] * b[2],
50        v[0] * b[1] - v[1] * b[0],
51    ];
52    [
53        q * (e[0] + vxb[0]),
54        q * (e[1] + vxb[1]),
55        q * (e[2] + vxb[2]),
56    ]
57}
58
59pub fn poynting_vector(e: [f64; 3], b: [f64; 3]) -> [f64; 3] {
60    let inv_mu0 = 1.0 / MU_0;
61    [
62        inv_mu0 * (e[1] * b[2] - e[2] * b[1]),
63        inv_mu0 * (e[2] * b[0] - e[0] * b[2]),
64        inv_mu0 * (e[0] * b[1] - e[1] * b[0]),
65    ]
66}
67
68pub fn energy_density_em(e: [f64; 3], b: [f64; 3]) -> f64 {
69    let e_sq = e[0] * e[0] + e[1] * e[1] + e[2] * e[2];
70    let b_sq = b[0] * b[0] + b[1] * b[1] + b[2] * b[2];
71    0.5 * EPSILON_0 * e_sq + 0.5 * b_sq / MU_0
72}
73
74pub fn electric_dipole_field(p: [f64; 3], r: [f64; 3]) -> [f64; 3] {
75    let r_mag = (r[0] * r[0] + r[1] * r[1] + r[2] * r[2]).sqrt();
76    if r_mag < 1e-30 {
77        return [0.0; 3];
78    }
79    let r5 = r_mag.powi(5);
80    let r3 = r_mag.powi(3);
81    let p_dot_r = p[0] * r[0] + p[1] * r[1] + p[2] * r[2];
82    let factor = 1.0 / (4.0 * std::f64::consts::PI * EPSILON_0);
83    [
84        factor * (3.0 * p_dot_r * r[0] / r5 - p[0] / r3),
85        factor * (3.0 * p_dot_r * r[1] / r5 - p[1] / r3),
86        factor * (3.0 * p_dot_r * r[2] / r5 - p[2] / r3),
87    ]
88}
89
90pub fn magnetic_dipole_field(m: [f64; 3], r: [f64; 3]) -> [f64; 3] {
91    let r_mag = (r[0] * r[0] + r[1] * r[1] + r[2] * r[2]).sqrt();
92    if r_mag < 1e-30 {
93        return [0.0; 3];
94    }
95    let r5 = r_mag.powi(5);
96    let r3 = r_mag.powi(3);
97    let m_dot_r = m[0] * r[0] + m[1] * r[1] + m[2] * r[2];
98    let factor = MU_0 / (4.0 * std::f64::consts::PI);
99    [
100        factor * (3.0 * m_dot_r * r[0] / r5 - m[0] / r3),
101        factor * (3.0 * m_dot_r * r[1] / r5 - m[1] / r3),
102        factor * (3.0 * m_dot_r * r[2] / r5 - m[2] / r3),
103    ]
104}
105
106pub fn capacitance_parallel_plate(area: f64, distance: f64, epsilon_r: f64) -> f64 {
107    epsilon_r * EPSILON_0 * area / distance
108}
109
110pub fn inductance_solenoid(n_turns: f64, length: f64, area: f64) -> f64 {
111    MU_0 * n_turns * n_turns * area / length
112}
113
114pub fn cyclotron_frequency(charge: f64, mass: f64, b: f64) -> f64 {
115    charge.abs() * b / mass
116}
117
118pub fn larmor_radius(mass: f64, v_perp: f64, charge: f64, b: f64) -> f64 {
119    mass * v_perp / (charge.abs() * b)
120}
121
122pub fn plasma_frequency(number_density: f64, mass: f64, charge: f64) -> f64 {
123    (number_density * charge * charge / (EPSILON_0 * mass)).sqrt()
124}
125
126pub fn debye_length(temperature: f64, number_density: f64, charge: f64) -> f64 {
127    use crate::constants::K_B;
128    (EPSILON_0 * K_B * temperature / (number_density * charge * charge)).sqrt()
129}