plot3d/
differencing.rs

1//! Forward and backward differencing for structured grid edges.
2
3use crate::block::Block;
4use crate::Float;
5
6/// Backward and forward difference pair: `(backward, forward)`.
7/// Each is a displacement vector `[dx, dy, dz]`.
8pub type DiffPair = ([Float; 3], [Float; 3]);
9
10/// Differencing data at a single node on a 2D face.
11#[derive(Clone, Debug)]
12pub struct FaceDiff {
13    pub p: usize,
14    pub q: usize,
15    /// Backward and forward differences along p.
16    pub dp: DiffPair,
17    /// Backward and forward differences along q.
18    pub dq: DiffPair,
19}
20
21/// Differencing data at a single node in a 3D block.
22#[derive(Clone, Debug)]
23pub struct BlockDiff {
24    pub i: usize,
25    pub j: usize,
26    pub k: usize,
27    /// Backward and forward differences along i.
28    pub di: DiffPair,
29    /// Backward and forward differences along j.
30    pub dj: DiffPair,
31    /// Backward and forward differences along k.
32    pub dk: DiffPair,
33}
34
35/// Compute forward and backward differences along each direction for a 2D face.
36///
37/// `x`, `y`, `z` are flat arrays of length `pmax * qmax`, stored row-major
38/// (p varies fastest within each row of q).
39pub fn find_face_edges(
40    x: &[Float],
41    y: &[Float],
42    z: &[Float],
43    pmax: usize,
44    qmax: usize,
45) -> Vec<FaceDiff> {
46    let idx = |p: usize, q: usize| -> usize { q * pmax + p };
47    let mut result = Vec::with_capacity(pmax * qmax);
48
49    for p in 0..pmax {
50        for q in 0..qmax {
51            let id = idx(p, q);
52
53            // dp: backward and forward along p
54            let dp_b = if p > 0 {
55                let prev = idx(p - 1, q);
56                [x[prev] - x[id], y[prev] - y[id], z[prev] - z[id]]
57            } else {
58                [0.0, 0.0, 0.0]
59            };
60            let dp_f = if p < pmax - 1 {
61                let next = idx(p + 1, q);
62                [x[next] - x[id], y[next] - y[id], z[next] - z[id]]
63            } else {
64                [0.0, 0.0, 0.0]
65            };
66
67            // dq: backward and forward along q
68            let dq_b = if q > 0 {
69                let prev = idx(p, q - 1);
70                [x[prev] - x[id], y[prev] - y[id], z[prev] - z[id]]
71            } else {
72                [0.0, 0.0, 0.0]
73            };
74            let dq_f = if q < qmax - 1 {
75                let next = idx(p, q + 1);
76                [x[next] - x[id], y[next] - y[id], z[next] - z[id]]
77            } else {
78                [0.0, 0.0, 0.0]
79            };
80
81            result.push(FaceDiff {
82                p,
83                q,
84                dp: (dp_b, dp_f),
85                dq: (dq_b, dq_f),
86            });
87        }
88    }
89    result
90}
91
92/// Compute forward and backward differences along each direction for a 3D block.
93pub fn find_edges(block: &Block) -> Vec<BlockDiff> {
94    let (ni, nj, nk) = (block.imax, block.jmax, block.kmax);
95    let mut result = Vec::with_capacity(ni * nj * nk);
96
97    for i in 0..ni {
98        for j in 0..nj {
99            for k in 0..nk {
100                let (cx, cy, cz) = block.xyz(i, j, k);
101
102                // di: backward and forward along i
103                let di_b = if i > 0 {
104                    let (px, py, pz) = block.xyz(i - 1, j, k);
105                    [px - cx, py - cy, pz - cz]
106                } else {
107                    [0.0, 0.0, 0.0]
108                };
109                let di_f = if i < ni - 1 {
110                    let (px, py, pz) = block.xyz(i + 1, j, k);
111                    [px - cx, py - cy, pz - cz]
112                } else {
113                    [0.0, 0.0, 0.0]
114                };
115
116                // dj: backward and forward along j
117                let dj_b = if j > 0 {
118                    let (px, py, pz) = block.xyz(i, j - 1, k);
119                    [px - cx, py - cy, pz - cz]
120                } else {
121                    [0.0, 0.0, 0.0]
122                };
123                let dj_f = if j < nj - 1 {
124                    let (px, py, pz) = block.xyz(i, j + 1, k);
125                    [px - cx, py - cy, pz - cz]
126                } else {
127                    [0.0, 0.0, 0.0]
128                };
129
130                // dk: backward and forward along k
131                let dk_b = if k > 0 {
132                    let (px, py, pz) = block.xyz(i, j, k - 1);
133                    [px - cx, py - cy, pz - cz]
134                } else {
135                    [0.0, 0.0, 0.0]
136                };
137                let dk_f = if k < nk - 1 {
138                    let (px, py, pz) = block.xyz(i, j, k + 1);
139                    [px - cx, py - cy, pz - cz]
140                } else {
141                    [0.0, 0.0, 0.0]
142                };
143
144                result.push(BlockDiff {
145                    i,
146                    j,
147                    k,
148                    di: (di_b, di_f),
149                    dj: (dj_b, dj_f),
150                    dk: (dk_b, dk_f),
151                });
152            }
153        }
154    }
155    result
156}
plot3d/differencing.rs

plot3d/
differencing.rs
