Skip to main content

powerio_matrix/io/
mtx.rs

1//! Matrix Market I/O.
2//!
3//! `sprs::io::write_matrix_market_sym` writes the *upper* triangle, but the
4//! Matrix Market spec calls for the *lower* triangle (i ≥ j). To stay
5//! compatible with strict readers (e.g. `fast_matrix_market`), we hand roll
6//! the symmetric writer. We delegate to `sprs` for general (non symmetric)
7//! output and for reading.
8
9use std::io::{BufWriter, Write};
10use std::path::Path;
11
12use sprs::CsMat;
13
14use crate::{Error, Result};
15
16pub fn write_mtx(matrix: &CsMat<f64>, path: impl AsRef<Path>) -> Result<()> {
17    let path = path.as_ref();
18    if is_structurally_symmetric(matrix) {
19        write_symmetric_mtx(matrix, path)
20    } else {
21        sprs::io::write_matrix_market(path, matrix.view()).map_err(|e| Error::Mtx(e.to_string()))
22    }
23}
24
25fn write_symmetric_mtx(matrix: &CsMat<f64>, path: &Path) -> Result<()> {
26    let f = std::fs::File::create(path)?;
27    let mut w = BufWriter::new(f);
28    writeln!(w, "%%MatrixMarket matrix coordinate real symmetric")?;
29    writeln!(w, "% written by powerio")?;
30
31    // Two-pass: count entries first so the header can carry the exact nnz.
32    let nnz = matrix
33        .iter()
34        .filter(|&(_, (i, j))| i >= j)
35        .filter(|&(&v, _)| v != 0.0)
36        .count();
37    writeln!(w, "{} {} {}", matrix.rows(), matrix.cols(), nnz)?;
38
39    for (&v, (i, j)) in matrix {
40        if i < j || v == 0.0 {
41            continue;
42        }
43        writeln!(w, "{} {} {:.16e}", i + 1, j + 1, v)?;
44    }
45    Ok(())
46}
47
48/// Read a Matrix Market file into a CSR matrix.
49pub fn read_mtx(path: impl AsRef<Path>) -> Result<CsMat<f64>> {
50    let tri: sprs::TriMat<f64> =
51        sprs::io::read_matrix_market(path).map_err(|e| Error::Mtx(e.to_string()))?;
52    Ok(tri.to_csr())
53}
54
55/// Read a dense vector written by [`write_vector_mtx`] (`array real general`):
56/// `%`-comment lines, a `<len> 1` dimensions line, then one value per line.
57pub fn read_vector_mtx(path: impl AsRef<Path>) -> Result<Vec<f64>> {
58    let text = std::fs::read_to_string(path)?;
59    let mut lines = text.lines().filter(|l| !l.starts_with('%'));
60    let header = lines
61        .next()
62        .ok_or_else(|| Error::Mtx("empty vector file".into()))?;
63    let len: usize = header
64        .split_whitespace()
65        .next()
66        .and_then(|t| t.parse().ok())
67        .ok_or_else(|| Error::Mtx(format!("bad vector dimensions line: {header:?}")))?;
68    let values = lines
69        .take(len)
70        .map(|l| {
71            l.trim()
72                .parse::<f64>()
73                .map_err(|_| Error::Mtx(format!("bad vector entry: {l:?}")))
74        })
75        .collect::<Result<Vec<_>>>()?;
76    if values.len() != len {
77        return Err(Error::Mtx(format!(
78            "expected {len} entries, got {}",
79            values.len()
80        )));
81    }
82    Ok(values)
83}
84
85/// Write a dense vector as Matrix Market `array real general`.
86pub fn write_vector_mtx(values: &[f64], path: impl AsRef<Path>) -> Result<()> {
87    let f = std::fs::File::create(path)?;
88    let mut w = BufWriter::new(f);
89    writeln!(w, "%%MatrixMarket matrix array real general")?;
90    writeln!(w, "% written by powerio")?;
91    writeln!(w, "{} 1", values.len())?;
92    for v in values {
93        writeln!(w, "{v:.16e}")?;
94    }
95    Ok(())
96}
97
98fn is_structurally_symmetric(a: &CsMat<f64>) -> bool {
99    if a.rows() != a.cols() {
100        return false;
101    }
102    for (i, row) in a.outer_iterator().enumerate() {
103        for (j, &v) in row.iter() {
104            let mirror = a.get(j, i).copied().unwrap_or(0.0);
105            if (v - mirror).abs() > 1e-12 {
106                return false;
107            }
108        }
109    }
110    true
111}