Skip to main content

powerio_dist/cim/
mod.rs

1//! Distribution CIM: IEC 61968 and 61970 CIM100 RDF/XML for unbalanced
2//! feeders, read into and written from the [`MulticonductorNetwork`]
3//! (format token `cim-dist`).
4//!
5//! The profile is the one the GridAPPS-D and CIMHub tools exchange: one
6//! RDF/XML document in the `http://iec.ch/TC57/CIM100#` namespace (the
7//! iec61970cim17v40 and iec61968cim13v13 release CGMES 3.0 also uses), with
8//! `ConnectivityNode` connectivity and per phase equipment
9//! (`ACLineSegmentPhase`, `EnergyConsumerPhase`, `SwitchPhase`, …). The
10//! reader also accepts the older versioned `CIM-schema-cimNN` namespaces.
11//! A directory of `.xml` documents reads as one case, with `rdf:about`
12//! objects extending the ones another document defines.
13//!
14//! CIM16 documents may state an enumeration as text rather than as a
15//! reference, and draw the network in a Diagram Layout; the reader reads
16//! both.
17//!
18//! Same format emission returns the retained source bytes of a single
19//! document. Canonical output is one CIM100 document with deterministic
20//! mRIDs: imported ids pass through the `cim_mrid` extras, and everything
21//! else takes a UUIDv5, so write, read, write is byte stable.
22//!
23//! # Phase convention
24//!
25//! CIM `SinglePhaseKind` values map to the wire coordinate terminal names
26//! (OpenDSS node numbers) the rest of the crate uses: `A`→`1`, `B`→`2`,
27//! `C`→`3`, `N`→`4`, and the split phase `s1`/`s2`→`1`/`2`. A conductor
28//! beyond those four has no `SinglePhaseKind` and is reported on output.
29//!
30//! # Transformer forms
31//!
32//! A `PowerTransformer` states its windings either as `PowerTransformerEnd`
33//! records with their own ratings and impedances, or as `TransformerTank`s,
34//! each with `TransformerTankEnd`s and catalog data (`TransformerTankInfo`,
35//! `TransformerEndInfo`, `ShortCircuitTest`, `NoLoadTest`) that an `Asset`
36//! links. Each tank reads as one transformer, with its parent and catalog
37//! entry in extras, and is written back as a tank; every other transformer is
38//! written in the `PowerTransformerEnd` form. No-load data from either form
39//! lands in the `%noloadloss` and `%imag` extras.
40//!
41//! [`MulticonductorNetwork`]: crate::MulticonductorNetwork
42
43mod read;
44mod write;
45mod xml;
46
47use std::f64::consts::PI;
48
49use crate::model::ConductorMatrix;
50
51pub(crate) use read::parse_cim_collecting;
52pub(crate) use write::emit_cim_dist_text;
53
54/// Per phase and distribution only classes: CGMES never uses them, so one of
55/// them in a header-less CIM document marks distribution CIM. The IEC 61968
56/// customer and metering classes (`ServiceLocation`, `UsagePoint`) mark a
57/// distribution feeder stated with sequence impedances and no per phase
58/// class.
59const DISTRIBUTION_MARKERS: [&[u8]; 12] = [
60    b"ACLineSegmentPhase",
61    b"EnergyConsumerPhase",
62    b"EnergySourcePhase",
63    b"SwitchPhase",
64    b"ShuntCompensatorPhase",
65    b"PerLengthPhaseImpedance",
66    b"PhaseImpedanceData",
67    b"TransformerTank",
68    b"PowerElectronicsConnectionPhase",
69    b"ConcentricNeutralCableInfo",
70    b"ServiceLocation",
71    b"UsagePoint",
72];
73
74/// How much of a document the header checks read.
75pub const CIM_DIST_HEAD_BYTES: usize = 16 * 1024;
76
77/// Whether a document opens like distribution CIM: RDF/XML in a CIM namespace
78/// with no `md:FullModel` header (which marks a CGMES profile) in its first
79/// [`CIM_DIST_HEAD_BYTES`]. A caller holding only a file's head can rule it
80/// out with this before reading the whole file for [`looks_like_cim_dist`].
81#[must_use]
82pub fn looks_like_cim_dist_header(head: &[u8]) -> bool {
83    let head = &head[..head.len().min(CIM_DIST_HEAD_BYTES)];
84    has(head, b"22-rdf-syntax-ns#")
85        && (has(head, b"CIM100#") || has(head, b"CIM-schema-cim"))
86        && !has(head, b"FullModel")
87}
88
89/// Whether a document reads as distribution CIM: its header passes
90/// [`looks_like_cim_dist_header`], and it holds at least one per phase or
91/// distribution only class anywhere.
92///
93/// The check scans bytes, so it never requires the document to be UTF-8.
94#[must_use]
95pub fn looks_like_cim_dist(bytes: &[u8]) -> bool {
96    looks_like_cim_dist_header(bytes)
97        && DISTRIBUTION_MARKERS.iter().any(|marker| has(bytes, marker))
98}
99
100fn has(haystack: &[u8], needle: &[u8]) -> bool {
101    memchr::memmem::find(haystack, needle).is_some()
102}
103
104/// Whether a directory entry is a document of the case: an `.xml` file at
105/// the directory's top level. Subdirectories (archives, schemas, other
106/// cases) are not part of it.
107fn is_case_document(entry: &str) -> bool {
108    !entry.contains('/')
109        && std::path::Path::new(entry)
110            .extension()
111            .and_then(|extension| extension.to_str())
112            .is_some_and(|extension| extension.eq_ignore_ascii_case("xml"))
113}
114
115/// Whether a source holds distribution CIM: its own buffer reads as
116/// distribution CIM ([`looks_like_cim_dist`]), or, for a directory, one of the
117/// `.xml` documents at its top level does. The reads go through the source,
118/// so they count against its acquisition budget and are cached for the parse
119/// that follows.
120#[must_use]
121pub fn source_looks_like_cim_dist(source: &powerio_core::Source) -> bool {
122    if !source.is_directory() {
123        return source
124            .primary_buffer()
125            .is_ok_and(|buffer| looks_like_cim_dist(buffer.content_bytes()));
126    }
127    let Ok(entries) = source.entry_names() else {
128        return false;
129    };
130    entries.iter().any(|entry| {
131        is_case_document(entry.as_str())
132            && source
133                .root_buffer(entry.as_str())
134                .is_ok_and(|buffer| looks_like_cim_dist(buffer.content_bytes()))
135    })
136}
137
138// The reader and the writer share these so the writer computes exactly what
139// the reader does: write, read, write stays byte stable.
140
141/// Impedance base `V²/S` (ohm), or None when degenerate.
142fn z_base(v_ref: f64, s_rating: f64) -> Option<f64> {
143    (v_ref > 0.0 && s_rating > 0.0).then_some(v_ref * v_ref / s_rating)
144}
145
146/// Expand zero and positive sequence scalars into an `n×n` phase matrix:
147/// self `(z0 + 2·z1)/3` on the diagonal, mutual `(z0 − z1)/3` off it.
148fn sequence_to_phase(z0: f64, z1: f64, n: usize) -> ConductorMatrix {
149    let self_z = (z0 + 2.0 * z1) / 3.0;
150    let mutual = (z0 - z1) / 3.0;
151    let mut m = vec![vec![mutual; n]; n];
152    for (i, row) in m.iter_mut().enumerate() {
153        row[i] = self_z;
154    }
155    m
156}
157
158/// An angle in radians wrapped to `(-π, π]`.
159fn wrap_angle(angle: f64) -> f64 {
160    let mut a = angle % (2.0 * PI);
161    if a <= -PI {
162        a += 2.0 * PI;
163    } else if a > PI {
164        a -= 2.0 * PI;
165    }
166    a
167}
168
169/// The standard CIM leakage magnitude. The alternative is test-only evidence
170/// of why producer-specific X values cannot be used as standard magnitudes.
171#[derive(Clone, Copy, Debug, PartialEq, Eq)]
172enum LeakageReading {
173    Magnitude,
174    #[cfg(test)]
175    Reactance,
176}
177
178const LEAKAGE_READING: LeakageReading = LeakageReading::Magnitude;
179
180/// The phase clock represented by the existing OpenDSS winding convention:
181/// equal three-phase connections are clock 0; mixed D/Y connections use the
182/// default lagging clock 1. A reversed single-phase neutral pair is clock 6.
183/// Other vector groups need an explicit orientation mapping and are refused.
184fn tank_phase_clock(
185    primary: &crate::DistWinding,
186    winding: &crate::DistWinding,
187    three_phase: bool,
188) -> f64 {
189    if three_phase {
190        if primary.conn == winding.conn {
191            0.0
192        } else {
193            1.0
194        }
195    } else if (primary.terminal_map.first().is_some_and(|t| t == "4"))
196        != (winding.terminal_map.first().is_some_and(|t| t == "4"))
197    {
198        6.0
199    } else {
200        0.0
201    }
202}
203
204/// One `ShortCircuitTest` as percentages on winding 1's power base.
205#[derive(Clone, Copy, Debug, PartialEq)]
206struct ShortCircuitPct {
207    /// The pair reactance.
208    x: f64,
209    /// The pair resistance from the load loss.
210    r: f64,
211}
212
213/// A short circuit test's leakage impedance (ohm, referred to the energised
214/// end of rated voltage `rated_u`) and load loss (kW), on a test base of
215/// `base_power` VA, as percentages on winding 1's rating `s_1`. `None` when
216/// a voltage or power is not positive.
217fn short_circuit_test_pct(
218    leakage: f64,
219    loss_kw: f64,
220    rated_u: f64,
221    base_power: f64,
222    s_1: f64,
223    reading: LeakageReading,
224) -> Option<ShortCircuitPct> {
225    let base = z_base(rated_u, base_power)?;
226    if s_1 <= 0.0 {
227        return None;
228    }
229    let z = 100.0 * leakage / base;
230    let r = 100.0 * 1000.0 * loss_kw / base_power;
231    let x = match reading {
232        LeakageReading::Magnitude => (z * z - r * r).max(0.0).sqrt(),
233        #[cfg(test)]
234        LeakageReading::Reactance => z,
235    };
236    let scale = s_1 / base_power;
237    Some(ShortCircuitPct {
238        x: x * scale,
239        r: r * scale,
240    })
241}
242
243/// A no-load test's loss (kW) and exciting current (percent of rated
244/// current at `base_power`) as the OpenDSS `(%noloadloss, %imag)` on winding
245/// 1's rating `s_1`. `None` when a power is not positive.
246fn no_load_test_pct(loss_kw: f64, exciting: f64, base_power: f64, s_1: f64) -> Option<(f64, f64)> {
247    if base_power <= 0.0 || s_1 <= 0.0 {
248        return None;
249    }
250    let loss = 100.0 * 1000.0 * loss_kw / s_1;
251    let iexc = exciting * base_power / s_1;
252    let imag = if loss > 0.0 && loss <= iexc {
253        (iexc * iexc - loss * loss).sqrt()
254    } else {
255        iexc
256    };
257    Some((loss, imag))
258}
259
260#[cfg(test)]
261mod tests {
262    use std::f64::consts::PI;
263
264    use super::{
265        LeakageReading, is_case_document, looks_like_cim_dist, no_load_test_pct, sequence_to_phase,
266        short_circuit_test_pct, wrap_angle,
267    };
268
269    fn close(a: f64, b: f64) -> bool {
270        (a - b).abs() <= 1e-9 * a.abs().max(b.abs()).max(1.0)
271    }
272
273    /// The IEEE 13 center tap tank `ct5` (2400/120/120 V, 5 kVA): its tests
274    /// state 23.5008, 23.5008 and 0.039168 ohm, exactly 2.04%, 2.04% and
275    /// 1.36% of the energised end's base, and only the 1 to 2 test states a
276    /// loss (0.09 kW, 1.8%). Read as a magnitude the 1 to 2 pair loses its
277    /// resistive part; read as a reactance every pair keeps the source value.
278    #[test]
279    fn leakage_impedance_interpretations() {
280        let tests = [
281            (23.5008, 0.09, 2400.0),
282            (23.5008, 0.0, 2400.0),
283            (0.039_168, 0.0, 120.0),
284        ];
285        let read = |reading| -> Vec<f64> {
286            tests
287                .iter()
288                .map(|&(z, loss, u)| {
289                    short_circuit_test_pct(z, loss, u, 5000.0, 5000.0, reading)
290                        .unwrap()
291                        .x
292                })
293                .collect()
294        };
295        let magnitude = read(LeakageReading::Magnitude);
296        assert!(close(
297            magnitude[0],
298            (2.04_f64.powi(2) - 1.8_f64.powi(2)).sqrt()
299        ));
300        assert!(close(magnitude[0], 0.96));
301        assert!(close(magnitude[1], 2.04));
302        assert!(close(magnitude[2], 1.36));
303        let reactance = read(LeakageReading::Reactance);
304        assert!(close(reactance[0], 2.04));
305        assert!(close(reactance[1], 2.04));
306        assert!(close(reactance[2], 1.36));
307        // The resistance comes from the loss either way.
308        let pair = short_circuit_test_pct(
309            23.5008,
310            0.09,
311            2400.0,
312            5000.0,
313            5000.0,
314            LeakageReading::Reactance,
315        )
316        .unwrap();
317        assert!(close(pair.r, 1.8));
318    }
319
320    /// The test base rescales to winding 1's rating, and a magnitude smaller
321    /// than its resistive part reads as no reactance.
322    #[test]
323    fn short_circuit_tests_rescale_and_flag_short_magnitudes() {
324        // 2.04% on a 5 kVA test base is 4.08% on a 10 kVA winding 1.
325        let pct = short_circuit_test_pct(
326            23.5008,
327            0.0,
328            2400.0,
329            5000.0,
330            10_000.0,
331            LeakageReading::Magnitude,
332        )
333        .unwrap();
334        assert!(close(pct.x, 4.08));
335        let short =
336            short_circuit_test_pct(1.0, 0.09, 2400.0, 5000.0, 5000.0, LeakageReading::Magnitude)
337                .unwrap();
338        assert!(close(short.x, 0.0));
339        assert!(
340            short_circuit_test_pct(1.0, 0.0, 0.0, 5000.0, 5000.0, LeakageReading::Magnitude)
341                .is_none()
342        );
343    }
344
345    /// 25 kVA with 0.05 kW of no-load loss and 0.538516% exciting current is
346    /// 0.2% loss and 0.5% magnetizing current.
347    #[test]
348    fn no_load_tests_become_percentages() {
349        let (loss, imag) = no_load_test_pct(0.05, 0.538_516, 25_000.0, 25_000.0).unwrap();
350        assert!(close(loss, 0.2));
351        assert!((imag - 0.5).abs() < 1e-6);
352        // A loss above the exciting current keeps the exciting current.
353        let (_, imag) = no_load_test_pct(0.05, 0.1, 25_000.0, 25_000.0).unwrap();
354        assert!(close(imag, 0.1));
355        // A test base other than winding 1's rating rescales the current.
356        let (_, imag) = no_load_test_pct(0.0, 0.5, 5000.0, 10_000.0).unwrap();
357        assert!(close(imag, 0.25));
358    }
359
360    #[test]
361    #[allow(clippy::float_cmp)]
362    fn sequence_expands_to_a_symmetric_phase_matrix() {
363        let m = sequence_to_phase(6.0, 3.0, 3);
364        // self = (6 + 6)/3 = 4, mutual = (6 - 3)/3 = 1.
365        assert_eq!(m[0][0], 4.0);
366        assert_eq!(m[0][1], 1.0);
367        assert_eq!(m[1][0], 1.0);
368    }
369
370    #[test]
371    fn angles_wrap_into_the_half_open_interval() {
372        assert!((wrap_angle(-4.0 * PI / 3.0) - 2.0 * PI / 3.0).abs() < 1e-12);
373        assert!((wrap_angle(PI) - PI).abs() < 1e-12);
374        assert!((wrap_angle(-PI) - PI).abs() < 1e-12);
375    }
376
377    #[test]
378    fn a_case_reads_only_the_top_level_xml_documents() {
379        assert!(is_case_document("feeder.xml"));
380        assert!(is_case_document("feeder.XML"));
381        assert!(!is_case_document("archive/feeder.xml"));
382        assert!(!is_case_document("notes.txt"));
383    }
384
385    const RDF: &str = r#"<rdf:RDF xmlns:rdf="http://www.w3.org/1999/02/22-rdf-syntax-ns#" xmlns:cim="http://iec.ch/TC57/CIM100#">"#;
386
387    #[test]
388    fn a_distribution_document_is_recognized_by_its_per_phase_classes() {
389        let doc = format!("{RDF}<cim:ACLineSegmentPhase rdf:ID=\"_p\"/></rdf:RDF>");
390        assert!(looks_like_cim_dist(doc.as_bytes()));
391    }
392
393    #[test]
394    fn a_feeder_with_usage_points_is_distribution_cim() {
395        let cim16 = r#"<rdf:RDF xmlns:rdf="http://www.w3.org/1999/02/22-rdf-syntax-ns#" xmlns:cim="http://iec.ch/TC57/2013/CIM-schema-cim16#">"#;
396        let doc = format!(
397            "{cim16}<cim:ACLineSegment rdf:ID=\"_l\"/><cim:UsagePoint rdf:ID=\"_u\"/></rdf:RDF>"
398        );
399        assert!(looks_like_cim_dist(doc.as_bytes()));
400    }
401
402    #[test]
403    fn a_cgmes_profile_or_a_document_without_markers_is_not() {
404        let header = format!(
405            "{RDF}<md:FullModel rdf:about=\"urn:uuid:1\"/><cim:ACLineSegmentPhase rdf:ID=\"_p\"/></rdf:RDF>"
406        );
407        assert!(!looks_like_cim_dist(header.as_bytes()));
408        let balanced = format!("{RDF}<cim:ACLineSegment rdf:ID=\"_l\"/></rdf:RDF>");
409        assert!(!looks_like_cim_dist(balanced.as_bytes()));
410        assert!(!looks_like_cim_dist(
411            b"<network xmlns=\"http://www.powsybl.org/schema/iidm/1_0\"/>"
412        ));
413    }
414}