diff --git a/.github/workflows/nalgebra-ci-build.yml b/.github/workflows/nalgebra-ci-build.yml index 37536518..3a56df13 100644 --- a/.github/workflows/nalgebra-ci-build.yml +++ b/.github/workflows/nalgebra-ci-build.yml @@ -71,10 +71,10 @@ jobs: - name: test nalgebra-sparse # Manifest-path is necessary because cargo otherwise won't correctly forward features # We increase number of proptest cases to hopefully catch more potential bugs - run: PROPTEST_CASES=10000 cargo test --manifest-path=nalgebra-sparse/Cargo.toml --features compare,proptest-support + run: PROPTEST_CASES=10000 cargo test --manifest-path=nalgebra-sparse/Cargo.toml --features compare,proptest-support,io - name: test nalgebra-sparse (slow tests) # Unfortunately, the "slow-tests" take so much time that we need to run them with --release - run: PROPTEST_CASES=10000 cargo test --release --manifest-path=nalgebra-sparse/Cargo.toml --features compare,proptest-support,slow-tests slow + run: PROPTEST_CASES=10000 cargo test --release --manifest-path=nalgebra-sparse/Cargo.toml --features compare,proptest-support,io,slow-tests slow test-nalgebra-macros: runs-on: ubuntu-latest steps: diff --git a/nalgebra-sparse/Cargo.toml b/nalgebra-sparse/Cargo.toml index 76fb408a..7692984e 100644 --- a/nalgebra-sparse/Cargo.toml +++ b/nalgebra-sparse/Cargo.toml @@ -16,6 +16,9 @@ license = "Apache-2.0" proptest-support = ["proptest", "nalgebra/proptest-support"] compare = [ "matrixcompare-core" ] +# Enable matrix market I/O +io = [ "pest", "pest_derive" ] + # Enable to enable running some tests that take a lot of time to run slow-tests = [] @@ -24,6 +27,8 @@ nalgebra = { version="0.29", path = "../" } num-traits = { version = "0.2", default-features = false } proptest = { version = "1.0", optional = true } matrixcompare-core = { version = "0.1.0", optional = true } +pest = { version = "2", optional = true } +pest_derive = { version = "2", optional = true } [dev-dependencies] itertools = "0.10" diff --git a/nalgebra-sparse/src/io/matrix_market.pest b/nalgebra-sparse/src/io/matrix_market.pest new file mode 100644 index 00000000..aeeba406 --- /dev/null +++ b/nalgebra-sparse/src/io/matrix_market.pest @@ -0,0 +1,53 @@ +WHITESPACE = _{ " "|"\t" } + +// + +Sparsity = {^"coordinate" | ^"array"} +DataType = {^"real" | ^"complex" | ^"pattern" | ^"integer" } +StorageScheme = {^"symmetric" | ^"general" | ^"skew-symmetric" | ^"hermitian"} +// Only consider matrices here. +Header = { ^"%%matrixmarket matrix" ~ Sparsity ~ DataType ~ StorageScheme } + +// + +Comments = _{ "%" ~ (!NEWLINE ~ ANY)* } + +// + +Dimension = @{ ASCII_DIGIT+ } +SparseShape = { Dimension ~ Dimension ~ Dimension} +DenseShape = { Dimension ~ Dimension} +Shape = {SparseShape | DenseShape } + +// + +// grammar from https://doc.rust-lang.org/std/primitive.f64.html#grammar + +Sign = {("+" | "-")} +Exp = @{ ^"e" ~ Sign? ~ ASCII_DIGIT+} +Number = @{ ((ASCII_DIGIT+ ~ "." ~ ASCII_DIGIT*) | (ASCII_DIGIT* ~ "." ~ASCII_DIGIT+) | ASCII_DIGIT+ ) ~ Exp? } +Real = @{ Sign? ~ ("inf" | "NaN" | Number) } + + +SparseReal = {Dimension~ Dimension~ Real } +SparseComplex = {Dimension ~ Dimension ~ Real ~ Real} +SparsePattern = {Dimension ~ Dimension} + +DenseReal = {Real} +DenseComplex = {Real ~ Real} + + +Entry = { SparseComplex | SparseReal | SparsePattern | DenseComplex | DenseReal } + +// end of file, a silent way, see https://github.com/pest-parser/pest/issues/304#issuecomment-427198507 +eoi = _{ !ANY } + +Document = { + SOI ~ + NEWLINE* ~ + Header ~ + (NEWLINE ~ Comments)* ~ + (NEWLINE ~ Shape) ~ + (NEWLINE ~ Entry?)* ~ + eoi +} \ No newline at end of file diff --git a/nalgebra-sparse/src/io/matrix_market.rs b/nalgebra-sparse/src/io/matrix_market.rs new file mode 100644 index 00000000..dea284ee --- /dev/null +++ b/nalgebra-sparse/src/io/matrix_market.rs @@ -0,0 +1,1331 @@ +//! Implementation of matrix market io code. +//! +//! See the [website](https://math.nist.gov/MatrixMarket/formats.html) or the [paper](https://www.researchgate.net/publication/2630533_The_Matrix_Market_Exchange_Formats_Initial_Design) for more details about matrix market. +use crate::coo::CooMatrix; +use crate::SparseFormatError; +use crate::SparseFormatErrorKind; +use nalgebra::Complex; +use pest::iterators::Pairs; +use pest::Parser; +use std::cmp::PartialEq; +use std::convert::Infallible; +use std::convert::TryFrom; +use std::fmt; +use std::fmt::Formatter; +use std::fs; +use std::num::ParseIntError; +use std::num::TryFromIntError; +use std::path::Path; +use std::str::FromStr; + +/// A description of the error that occurred during importing a matrix from a matrix market format data. +#[derive(Debug)] +pub struct MatrixMarketError { + error_kind: MatrixMarketErrorKind, + message: String, +} + +/// Errors produced by functions that expect well-formed matrix market format data. +#[non_exhaustive] +#[derive(Copy, Clone, Debug, PartialEq)] +pub enum MatrixMarketErrorKind { + /// Parsing failure. + /// + /// Indicates that the parser failed, for example due to an unexpected string. + /// + /// Examples + /// -------- + /// ```rust + /// # use nalgebra_sparse::io::load_coo_from_matrix_market_str; + /// # use nalgebra_sparse::io::MatrixMarketErrorKind; + /// let str = r#" + /// %%MatrixMarket invalid invalid invalid invalid + /// 1 1 1 + /// 1 1 5 + /// "#; + /// let matrix_result = load_coo_from_matrix_market_str::(str); + /// assert_eq!(matrix_result.is_err(), true); + /// assert_eq!(matrix_result.unwrap_err().kind(), MatrixMarketErrorKind::ParsingError); + /// ``` + ParsingError, + + /// Indicates that the matrix market header is invalid. + /// + /// Examples + /// -------- + /// ```rust + /// # use nalgebra_sparse::io::load_coo_from_matrix_market_str; + /// # use nalgebra_sparse::io::MatrixMarketErrorKind; + /// let str = r#" + /// %%MatrixMarket matrix coordinate real hermitian + /// % a real matrix can't be hermitian + /// 1 1 1 + /// 1 1 5 + /// "#; + /// let matrix_result = load_coo_from_matrix_market_str::(str); + /// assert_eq!(matrix_result.is_err(), true); + /// assert_eq!(matrix_result.unwrap_err().kind(),MatrixMarketErrorKind::InvalidHeader); + /// ``` + InvalidHeader, + + /// Indicates that the number of data entries in the matrix market file does not match the header. + /// + /// Examples + /// -------- + /// ```rust + /// # use nalgebra_sparse::io::load_coo_from_matrix_market_str; + /// # use nalgebra_sparse::io::MatrixMarketErrorKind; + /// let str = r#" + /// %%matrixmarket matrix coordinate real general + /// % it has one more data entry than specified. + /// 3 3 1 + /// 2 2 2 + /// 2 3 2 + /// "#; + /// let matrix_result = load_coo_from_matrix_market_str::(str); + /// assert_eq!(matrix_result.is_err(), true); + /// assert_eq!(matrix_result.unwrap_err().kind(),MatrixMarketErrorKind::EntryMismatch); + /// ``` + EntryMismatch, + + /// Indicates that the scalar type requested is not compatible with the scalar type stored + /// in the matrix market file. + /// + /// Examples + /// -------- + /// ```rust + /// # use nalgebra_sparse::io::load_coo_from_matrix_market_str; + /// # use nalgebra_sparse::io::MatrixMarketErrorKind; + /// let str = r#" + /// %%matrixmarket matrix coordinate real general + /// % it should be loaded with load_coo_from_matrix_market_str::(str) (or f32) + /// 3 3 2 + /// 2 2 2.22 + /// 2 3 2.22 + /// "#; + /// let matrix_result = load_coo_from_matrix_market_str::(str); + /// assert_eq!(matrix_result.is_err(), true); + /// assert_eq!(matrix_result.unwrap_err().kind(),MatrixMarketErrorKind::TypeMismatch); + /// ``` + TypeMismatch, + + /// Indicates that zero has been used as an index in the data. + /// + /// **Note**: The matrix market format uses 1-based indexing. + /// + /// Examples + /// -------- + /// ```rust + /// # use nalgebra_sparse::io::load_coo_from_matrix_market_str; + /// # use nalgebra_sparse::io::MatrixMarketErrorKind; + /// let str = r#" + /// %%matrixmarket matrix coordinate real general + /// 1 1 1 + /// 0 0 10 + /// "#; + /// let matrix_result = load_coo_from_matrix_market_str::(str); + /// assert_eq!(matrix_result.is_err(), true); + /// assert_eq!(matrix_result.unwrap_err().kind(),MatrixMarketErrorKind::ZeroError); + /// ``` + ZeroError, + + /// Indicates [SparseFormatError] while creating the sparse matrix. + /// + /// + /// Examples + /// -------- + /// ```rust + /// # use nalgebra_sparse::io::load_coo_from_matrix_market_str; + /// # use nalgebra_sparse::io::MatrixMarketErrorKind; + /// # use nalgebra_sparse::SparseFormatErrorKind; + /// let str = r#" + /// %%matrixmarket matrix coordinate real general + /// 1 1 1 + /// 4 2 10 + /// "#; + /// let matrix_result = load_coo_from_matrix_market_str::(str); + /// assert_eq!(matrix_result.is_err(), true); + /// assert_eq!(matrix_result.unwrap_err().kind(), + /// MatrixMarketErrorKind::SparseFormatError(SparseFormatErrorKind::IndexOutOfBounds)); + /// ``` + SparseFormatError(SparseFormatErrorKind), + + /// Indicates that a wrong diagonal element has been provided to the matrix. + /// + /// Examples + /// -------- + /// ```rust + /// # use nalgebra_sparse::io::load_coo_from_matrix_market_str; + /// # use nalgebra_sparse::io::MatrixMarketErrorKind; + /// # use nalgebra::Complex; + /// let str = r#" + /// %%matrixmarket matrix coordinate real skew-symmetric + /// % skew-symmetric matrix can't have element on diagonal + /// 5 5 2 + /// 1 1 10 + /// 2 1 5 + /// "#; + /// let matrix_result = load_coo_from_matrix_market_str::(str); + /// assert_eq!(matrix_result.is_err(), true); + /// assert_eq!(matrix_result.unwrap_err().kind(),MatrixMarketErrorKind::DiagonalError); + /// + /// let str = r#" + /// %%matrixmarket matrix coordinate complex hermitian + /// % hermitian matrix diagonal element must be a real number + /// 5 5 2 + /// 1 1 10 2 + /// 2 1 5 2 + /// "#; + /// let matrix_result = load_coo_from_matrix_market_str::>(str); + /// assert_eq!(matrix_result.is_err(), true); + /// assert_eq!(matrix_result.unwrap_err().kind(),MatrixMarketErrorKind::DiagonalError); + /// ``` + /// Here the skew matrix shouldn't have an element on the diagonal. + DiagonalError, + + /// Indicates an [IO error](`std::io::Error`) while reading the data from file. + /// + /// Examples + /// -------- + /// ```rust + /// # use nalgebra_sparse::io::load_coo_from_matrix_market_file; + /// # use nalgebra_sparse::io::MatrixMarketErrorKind; + /// let matrix_result = load_coo_from_matrix_market_file::("matrix.mtx"); + /// assert_eq!(matrix_result.is_err(), true); + /// assert_eq!(matrix_result.unwrap_err().kind(),MatrixMarketErrorKind::IOError(std::io::ErrorKind::NotFound)); + /// ``` + IOError(std::io::ErrorKind), + + /// Indicates that a (skew-)symmetric (or hermitian) matrix is not a lower triangular matrix. + /// + /// Examples + /// -------- + /// ```rust + /// # use nalgebra_sparse::io::load_coo_from_matrix_market_str; + /// # use nalgebra_sparse::io::MatrixMarketErrorKind; + /// let str = r#" + /// %%matrixmarket matrix coordinate integer symmetric + /// 5 5 2 + /// 1 1 10 + /// 2 3 5 + /// "#; + /// let matrix_result = load_coo_from_matrix_market_str::(str); + /// assert_eq!(matrix_result.is_err(), true); + /// assert_eq!(matrix_result.unwrap_err().kind(),MatrixMarketErrorKind::NotLowerTriangle); + /// ``` + NotLowerTriangle, + + /// Indicates that a (skew-)symmetric (or hermitian) matrix is not a square matrix. + /// + /// Examples + /// -------- + /// ```rust + /// # use nalgebra_sparse::io::load_coo_from_matrix_market_str; + /// # use nalgebra_sparse::io::MatrixMarketErrorKind; + /// let str = r#" + /// %%matrixmarket matrix coordinate integer symmetric + /// 5 4 2 + /// 1 1 10 + /// 3 2 5 + /// "#; + /// let matrix_result = load_coo_from_matrix_market_str::(str); + /// assert_eq!(matrix_result.is_err(), true); + /// assert_eq!(matrix_result.unwrap_err().kind(),MatrixMarketErrorKind::NonSquare); + /// ``` + NonSquare, +} + +impl MatrixMarketError { + fn from_kind_and_message(error_type: MatrixMarketErrorKind, message: String) -> Self { + Self { + error_kind: error_type, + message, + } + } + + /// The matrix market error kind. + #[must_use] + pub fn kind(&self) -> MatrixMarketErrorKind { + self.error_kind + } + + /// The underlying error message. + #[must_use] + pub fn message(&self) -> &str { + self.message.as_str() + } +} + +impl fmt::Display for MatrixMarketError { + fn fmt(&self, f: &mut Formatter<'_>) -> fmt::Result { + write!(f, "Matrix Market error: ")?; + match self.kind() { + MatrixMarketErrorKind::ParsingError => { + write!(f, "ParsingError,")?; + } + MatrixMarketErrorKind::InvalidHeader => { + write!(f, "InvalidHeader,")?; + } + MatrixMarketErrorKind::EntryMismatch => { + write!(f, "EntryNumUnmatched,")?; + } + MatrixMarketErrorKind::TypeMismatch => { + write!(f, "TypeMismatch,")?; + } + MatrixMarketErrorKind::SparseFormatError(_) => { + write!(f, "SparseFormatError,")?; + } + MatrixMarketErrorKind::ZeroError => { + write!(f, "ZeroError,")?; + } + MatrixMarketErrorKind::IOError(_) => { + write!(f, "IOError,")?; + } + MatrixMarketErrorKind::DiagonalError => { + write!(f, "DiagonalError,")?; + } + MatrixMarketErrorKind::NotLowerTriangle => { + write!(f, "NotLowerTriangle,")?; + } + MatrixMarketErrorKind::NonSquare => { + write!(f, "NotSquareMatrix,")?; + } + } + write!(f, " message: {}", self.message) + } +} + +impl std::error::Error for MatrixMarketError {} + +impl MatrixMarketError { + fn from_pest_error(error: pest::error::Error) -> Self + where + T: fmt::Debug + std::hash::Hash + std::marker::Copy + Ord, + { + Self::from_kind_and_message( + MatrixMarketErrorKind::ParsingError, + format!("Can't parse the data.\n Error: {}", error), + ) + } +} + +impl From for MatrixMarketError { + fn from(err: ParseIntError) -> Self { + Self::from_kind_and_message( + MatrixMarketErrorKind::ParsingError, + format!("Can't parse data as i128.\n Error: {}", err), + ) + } +} + +impl From for MatrixMarketError { + fn from(err: SparseFormatError) -> Self { + Self::from_kind_and_message( + MatrixMarketErrorKind::SparseFormatError(*err.kind()), + format!("{}", &err), + ) + } +} + +impl From for MatrixMarketError { + fn from(err: std::io::Error) -> Self { + Self::from_kind_and_message( + MatrixMarketErrorKind::IOError(err.kind()), + format!("{}", &err), + ) + } +} + +impl From for MatrixMarketError { + fn from(err: TryFromIntError) -> Self { + Self::from_kind_and_message( + MatrixMarketErrorKind::TypeMismatch, + format!( + "Please consider using a larger integer type. Error message: {}", + &err + ), + ) + } +} + +// This is needed when calling `i128::try_from(i: i128)` +// but it won't happen +impl From for MatrixMarketError { + fn from(_err: Infallible) -> Self { + Self::from_kind_and_message( + MatrixMarketErrorKind::TypeMismatch, + format!("This won't happen"), + ) + } +} + +#[derive(Debug, PartialEq)] +enum Sparsity { + Sparse, + Dense, +} +#[derive(Debug, PartialEq)] +enum DataType { + Real, + Complex, + Pattern, + Integer, +} +#[derive(Debug, PartialEq)] +enum StorageScheme { + Symmetric, + General, + Skew, + Hermitian, +} +#[derive(Debug, PartialEq)] +struct Typecode { + sparsity: Sparsity, + datatype: DataType, + storagescheme: StorageScheme, +} + +impl FromStr for Sparsity { + type Err = MatrixMarketError; + /// Assumes that `word` is already lower case. + fn from_str(word: &str) -> Result { + match word { + "coordinate" => Ok(Sparsity::Sparse), + "array" => Ok(Sparsity::Dense), + _ => Err(MatrixMarketError::from_kind_and_message( + MatrixMarketErrorKind::ParsingError, + format!("keyword {} is unknown", word), + )), + } + } +} + +impl FromStr for DataType { + type Err = MatrixMarketError; + /// Assumes that `word` is already lower case. + fn from_str(word: &str) -> Result { + match word { + "real" => Ok(DataType::Real), + "complex" => Ok(DataType::Complex), + "integer" => Ok(DataType::Integer), + "pattern" => Ok(DataType::Pattern), + _ => Err(MatrixMarketError::from_kind_and_message( + MatrixMarketErrorKind::ParsingError, + format!("keyword {} is unknown", word), + )), + } + } +} + +impl FromStr for StorageScheme { + type Err = MatrixMarketError; + /// Assumes that `word` is already lower case. + fn from_str(word: &str) -> Result { + match word { + "skew-symmetric" => Ok(StorageScheme::Skew), + "general" => Ok(StorageScheme::General), + "symmetric" => Ok(StorageScheme::Symmetric), + "hermitian" => Ok(StorageScheme::Hermitian), + _ => Err(MatrixMarketError::from_kind_and_message( + MatrixMarketErrorKind::ParsingError, + format!("keyword {} is unknown", word), + )), + } + } +} + +/// Precheck if it's a valid header. +/// +/// For more details, please check +/// +/// Boisvert, Ronald F., Roldan Pozo, and Karin A. Remington. +/// The matrix market formats: Initial design. +/// Technical report, Applied and Computational Mathematics Division, NIST, 1996. Section 3. +fn typecode_precheck(tc: &Typecode) -> Result<(), MatrixMarketError> { + match tc { + Typecode { + datatype: DataType::Real, + storagescheme: StorageScheme::Hermitian, + .. + } => Err(MatrixMarketError::from_kind_and_message( + MatrixMarketErrorKind::InvalidHeader, + String::from("Real matrix can't be hermitian."), + )), + Typecode { + datatype: DataType::Integer, + storagescheme: StorageScheme::Hermitian, + .. + } => Err(MatrixMarketError::from_kind_and_message( + MatrixMarketErrorKind::InvalidHeader, + String::from("Integer matrix can't be hermitian."), + )), + Typecode { + datatype: DataType::Pattern, + storagescheme: StorageScheme::Hermitian, + .. + } => Err(MatrixMarketError::from_kind_and_message( + MatrixMarketErrorKind::InvalidHeader, + String::from("Pattern matrix can't be hermitian."), + )), + Typecode { + datatype: DataType::Pattern, + storagescheme: StorageScheme::Skew, + .. + } => Err(MatrixMarketError::from_kind_and_message( + MatrixMarketErrorKind::InvalidHeader, + String::from("Pattern matrix can't be skew-symmetric."), + )), + Typecode { + datatype: DataType::Pattern, + sparsity: Sparsity::Dense, + .. + } => Err(MatrixMarketError::from_kind_and_message( + MatrixMarketErrorKind::InvalidHeader, + String::from("Dense matrix can't be pattern matrix."), + )), + // precheck success + _ => Ok(()), + } +} + +/// Scalar types supported by the matrix market parser. +mod internal { + use crate::io::MatrixMarketError; + use na::{Complex, Scalar}; + + pub trait SupportedMatrixMarketScalar: Scalar { + /// When the matrix is an integer matrix, it will convert a [i128] number to this type. + fn from_i128(i: i128) -> Result; + /// When matrix is a Real matrix, it will convert a [f64] number to this type. + fn from_f64(f: f64) -> Result; + /// When matrix is a Complex matrix, it will convert a [Complex] number to this type. + fn from_c64(c: Complex) -> Result; + /// When matrix is a Pattern matrix, it will convert a unit type [unit] to this type. + fn from_pattern(p: ()) -> Result; + /// When matrix is a Skew-symmetric matrix, it will convert itself to its negative. + fn negative(self) -> Result; + /// When matrix is a Hermitian matrix, it will convert itself to its conjugate. + fn conjugate(self) -> Result; + } +} + +/// A marker trait for supported matrix market scalars. +/// +/// This is a sealed trait; it cannot be implemented by external crates. This is done in order to prevent leaking +/// some of the implementation details we currently rely on. We may relax this restriction in the future. +pub trait MatrixMarketScalar: internal::SupportedMatrixMarketScalar {} + +/// Implement MatrixMarketScalar for primitive integer types. +macro_rules! mm_int_impl { + ($T:ty) => { + impl MatrixMarketScalar for $T {} + + impl internal::SupportedMatrixMarketScalar for $T { + #[inline] + fn from_i128(i: i128) -> Result { + Ok(Self::try_from(i)?) + } + #[inline] + fn from_f64(_f: f64) -> Result { + Err(MatrixMarketError::from_kind_and_message( + MatrixMarketErrorKind::TypeMismatch, + format!("Int type can't be parsed from f64"), + )) + } + #[inline] + fn from_c64(_c: Complex) -> Result { + Err(MatrixMarketError::from_kind_and_message( + MatrixMarketErrorKind::TypeMismatch, + format!("Int type can't be parsed from Complex"), + )) + } + #[inline] + fn from_pattern(_p: ()) -> Result { + Err(MatrixMarketError::from_kind_and_message( + MatrixMarketErrorKind::TypeMismatch, + format!("Int type can't be parsed from ()"), + )) + } + #[inline] + fn conjugate(self) -> Result { + Err(MatrixMarketError::from_kind_and_message( + MatrixMarketErrorKind::TypeMismatch, + format!("Int type has no conjugate"), + )) + } + #[inline] + fn negative(self) -> Result { + Ok(-self) + } + } + }; +} +/// Implement MatrixMarketScalar for primitive real types. +macro_rules! mm_real_impl { + ($T:ty) => { + impl MatrixMarketScalar for $T {} + + impl internal::SupportedMatrixMarketScalar for $T { + #[inline] + fn from_i128(_i: i128) -> Result { + Err(MatrixMarketError::from_kind_and_message( + MatrixMarketErrorKind::TypeMismatch, + format!("real type can't be parsed from i128"), + )) + } + #[inline] + fn from_f64(f: f64) -> Result { + Ok(f as Self) + } + #[inline] + fn from_c64(_c: Complex) -> Result { + Err(MatrixMarketError::from_kind_and_message( + MatrixMarketErrorKind::TypeMismatch, + format!("real type can't be parsed from Complex"), + )) + } + #[inline] + fn from_pattern(_p: ()) -> Result { + Err(MatrixMarketError::from_kind_and_message( + MatrixMarketErrorKind::TypeMismatch, + format!("real type can't be parsed from ()"), + )) + } + #[inline] + fn conjugate(self) -> Result { + Err(MatrixMarketError::from_kind_and_message( + MatrixMarketErrorKind::TypeMismatch, + format!("real type has no conjugate"), + )) + } + #[inline] + fn negative(self) -> Result { + Ok(-self) + } + } + }; +} + +/// Implement MatrixMarketScalar for primitive complex types. +macro_rules! mm_complex_impl { + ($T:ty) => { + impl MatrixMarketScalar for Complex<$T> {} + + impl internal::SupportedMatrixMarketScalar for Complex<$T> { + #[inline] + fn from_i128(_i: i128) -> Result { + Err(MatrixMarketError::from_kind_and_message( + MatrixMarketErrorKind::TypeMismatch, + format!("Complex type can't be parsed from i128"), + )) + } + #[inline] + fn from_f64(_f: f64) -> Result { + Err(MatrixMarketError::from_kind_and_message( + MatrixMarketErrorKind::TypeMismatch, + format!("Complex type can't be parsed from f64"), + )) + } + #[inline] + fn from_c64(c: Complex) -> Result { + Ok(Self { + re: c.re as $T, + im: c.im as $T, + }) + } + #[inline] + fn from_pattern(_p: ()) -> Result { + Err(MatrixMarketError::from_kind_and_message( + MatrixMarketErrorKind::TypeMismatch, + format!("Complex type can't be parsed from ()"), + )) + } + #[inline] + fn conjugate(self) -> Result { + Ok(self.conj()) + } + #[inline] + fn negative(self) -> Result { + Ok(-self) + } + } + }; +} +/// Implement MatrixMarketScalar for primitive unit types. +macro_rules! mm_pattern_impl { + ($T:ty) => { + impl MatrixMarketScalar for $T {} + + impl internal::SupportedMatrixMarketScalar for $T { + #[inline] + fn from_i128(_i: i128) -> Result { + Err(MatrixMarketError::from_kind_and_message( + MatrixMarketErrorKind::TypeMismatch, + format!("Pattern type can't be parsed from i128"), + )) + } + #[inline] + fn from_f64(_f: f64) -> Result { + Err(MatrixMarketError::from_kind_and_message( + MatrixMarketErrorKind::TypeMismatch, + format!("Pattern type can't be parsed from f64"), + )) + } + #[inline] + fn from_c64(_c: Complex) -> Result { + Err(MatrixMarketError::from_kind_and_message( + MatrixMarketErrorKind::TypeMismatch, + format!("Pattern type can't be parsed from Complex"), + )) + } + #[inline] + fn from_pattern(p: ()) -> Result { + Ok(p) + } + + #[inline] + fn conjugate(self) -> Result { + Err(MatrixMarketError::from_kind_and_message( + MatrixMarketErrorKind::TypeMismatch, + format!("Pattern type has no conjugate"), + )) + } + #[inline] + fn negative(self) -> Result { + Err(MatrixMarketError::from_kind_and_message( + MatrixMarketErrorKind::TypeMismatch, + format!("Pattern type has no negative"), + )) + } + } + }; +} + +mm_int_impl!(i8); +mm_int_impl!(i16); +mm_int_impl!(i32); +mm_int_impl!(i64); +mm_int_impl!(i128); + +mm_real_impl!(f32); +mm_real_impl!(f64); + +mm_complex_impl!(f32); +mm_complex_impl!(f64); + +mm_pattern_impl!(()); + +#[derive(Parser)] +#[grammar = "io/matrix_market.pest"] +struct MatrixMarketParser; + +/// Parses a Matrix Market file at the given path as a `CooMatrix`. +/// +/// The matrix market format specification does not clarify whether duplicate entries are allowed. Our importer +/// assumes that this is permitted and produces a `CooMatrix` with possibly duplicate entries. +/// +/// **Note**: A current restriction of the importer is that you must use a compatible scalar type when importing. +/// For example, in order to import a matrix stored as `integer` in the matrix market format, you must +/// import it as an integer matrix, otherwise a [TypeMismatch](MatrixMarketErrorKind::TypeMismatch) error +/// will be returned. This restriction may be lifted in the future, and is +/// tracked by issue [#1038](https://github.com/dimforge/nalgebra/issues/1038). +/// +/// Errors +/// -------- +/// +/// See [MatrixMarketErrorKind] for a list of possible error conditions. +/// +/// Examples +/// -------- +/// ```no_run +/// use nalgebra_sparse::io::load_coo_from_matrix_market_file; +/// // Use e.g. `i32` for integer matrices +/// let matrix = load_coo_from_matrix_market_file::("path/to/matrix.mtx").unwrap(); +/// ``` +pub fn load_coo_from_matrix_market_file>( + path: P, +) -> Result, MatrixMarketError> +where + T: MatrixMarketScalar, +{ + let file = fs::read_to_string(path)?; + load_coo_from_matrix_market_str(&file) +} + +/// Parses a Matrix Market file described by the given string as a `CooMatrix`. +/// +/// See [load_coo_from_matrix_market_file] for more information. +/// +/// Errors +/// -------- +/// +/// See [MatrixMarketErrorKind] for a list of possible error conditions. +/// +/// Examples +/// -------- +/// ``` +/// use nalgebra_sparse::io::load_coo_from_matrix_market_str; +/// let str = r#" +/// %%matrixmarket matrix coordinate integer general +/// 5 4 2 +/// 1 1 10 +/// 2 3 5 +/// "#; +/// // Use e.g. `i32` for integer matrices +/// let matrix = load_coo_from_matrix_market_str::(str).unwrap(); +/// ``` +pub fn load_coo_from_matrix_market_str(data: &str) -> Result, MatrixMarketError> +where + T: MatrixMarketScalar, +{ + // unwrap() in this function are guaranteed by parsing the data + let file = MatrixMarketParser::parse(Rule::Document, data) + .map_err(MatrixMarketError::from_pest_error)? + .next() + .unwrap(); + + let mut rows: Vec = Vec::new(); + let mut cols: Vec = Vec::new(); + let mut data: Vec = Vec::new(); + let mut lines = file.into_inner(); + + let header_line = lines.next().unwrap(); + let header_type = parse_header(&mut header_line.into_inner()); + typecode_precheck(&header_type)?; + + let shape_line = lines.next().unwrap(); + // shape here is number of rows, columns, non-zeros + let shape: (usize, usize, usize); + match header_type.sparsity { + Sparsity::Sparse => { + shape = parse_sparse_shape(&mut shape_line.into_inner(), &header_type.storagescheme)?; + } + Sparsity::Dense => { + shape = parse_dense_shape(&mut shape_line.into_inner(), &header_type.storagescheme)?; + } + } + + // used when constructing dense matrix. + // If it's sparse matrix, it has no effect. + let mut current_dense_coordinate: (usize, usize) = (0, 0); + if header_type.storagescheme == StorageScheme::Skew { + // for skew dense matrix, the first element starts from (1,0) + current_dense_coordinate = (1, 0); + } + // count how many entries in the matrix data + let count = lines.clone().count(); + if count != shape.2 { + return Err(MatrixMarketError::from_kind_and_message( + MatrixMarketErrorKind::EntryMismatch, + format!( + "{} entries required for the matrix, but {} was provided", + shape.2, count, + ), + )); + } + + for data_line in lines { + let entry: (usize, usize, T); + match header_type { + Typecode { + sparsity: Sparsity::Sparse, + datatype: DataType::Real, + .. + } => { + entry = parse_sparse_real::(&mut data_line.into_inner())?; + } + Typecode { + sparsity: Sparsity::Sparse, + datatype: DataType::Integer, + .. + } => { + entry = parse_sparse_int::(&mut data_line.into_inner())?; + } + Typecode { + sparsity: Sparsity::Sparse, + datatype: DataType::Pattern, + .. + } => { + entry = parse_sparse_pattern::(&mut data_line.into_inner())?; + } + Typecode { + sparsity: Sparsity::Sparse, + datatype: DataType::Complex, + .. + } => { + entry = parse_sparse_complex::(&mut data_line.into_inner())?; + } + Typecode { + sparsity: Sparsity::Dense, + datatype: DataType::Complex, + .. + } => { + entry = ( + current_dense_coordinate.0, + current_dense_coordinate.1, + parse_dense_complex::(&mut data_line.into_inner())?, + ); + next_dense_coordinate( + &mut current_dense_coordinate, + shape, + &header_type.storagescheme, + ); + } + Typecode { + sparsity: Sparsity::Dense, + datatype: DataType::Real, + .. + } => { + entry = ( + current_dense_coordinate.0, + current_dense_coordinate.1, + parse_dense_real::(&mut data_line.into_inner())?, + ); + + next_dense_coordinate( + &mut current_dense_coordinate, + shape, + &header_type.storagescheme, + ); + } + Typecode { + sparsity: Sparsity::Dense, + datatype: DataType::Integer, + .. + } => { + entry = ( + current_dense_coordinate.0, + current_dense_coordinate.1, + parse_dense_int::(&mut data_line.into_inner())?, + ); + next_dense_coordinate( + &mut current_dense_coordinate, + shape, + &header_type.storagescheme, + ); + } + _ => { + // it shouldn't happen here, because dense matrix can't be pattern. And it will give InvalidHeader error beforehand. + entry = (1, 1, T::from_i128(1)?) + } + } + + let (r, c, d) = entry; + + match header_type.storagescheme { + StorageScheme::General => { + rows.push(r); + cols.push(c); + data.push(d); + } + StorageScheme::Symmetric => { + check_lower_triangle(r, c)?; + rows.push(r); + cols.push(c); + data.push(d.clone()); + // don't need to add twice if the element in on diagonal + if r != c { + rows.push(c); + cols.push(r); + data.push(d); + } + } + StorageScheme::Skew => { + check_lower_triangle(r, c)?; + rows.push(r); + cols.push(c); + data.push(d.clone()); + // skew-symmetric matrix shouldn't have diagonal element + if r == c { + return Err(MatrixMarketError::from_kind_and_message( + MatrixMarketErrorKind::DiagonalError, + format!( + "There is a diagonal element in skew matrix, in row(and column) {}", + r + 1 + ), + )); + } + rows.push(c); + cols.push(r); + data.push(d.negative()?); + } + StorageScheme::Hermitian => { + check_lower_triangle(r, c)?; + rows.push(r); + cols.push(c); + data.push(d.clone()); + + if r == c && d != d.clone().conjugate()? { + return Err(MatrixMarketError::from_kind_and_message( + MatrixMarketErrorKind::DiagonalError, + format!( + "There is a diagonal element in hermitian matrix, which is not a real number, in row(and column) {}", + r + 1 + ), + )); + } + // don't need to add twice if the element in on diagonal + if r != c { + rows.push(c); + cols.push(r); + data.push(d.conjugate()?); + } + } + } + } + Ok(CooMatrix::try_from_triplets( + shape.0, shape.1, rows, cols, data, + )?) +} + +#[inline] +/// do a quick check it the entry is in the lower triangle part of the matrix +fn check_lower_triangle(r: usize, c: usize) -> Result<(), MatrixMarketError> { + if c > r { + return Err(MatrixMarketError::from_kind_and_message( + MatrixMarketErrorKind::NotLowerTriangle, + format!( + "Entry: row {} col {} should be put into lower triangle", + r, c + ), + )); + } + Ok(()) +} + +#[inline] +/// Parse a pest structure to a Typecode of the matrix. +fn parse_header(inner: &mut Pairs<'_, Rule>) -> Typecode { + // unwrap() in this function are guaranteed by parsing the data + Typecode { + sparsity: inner + .next() + .unwrap() + .as_str() + .to_ascii_lowercase() + .parse::() + .unwrap(), + datatype: inner + .next() + .unwrap() + .as_str() + .to_ascii_lowercase() + .parse::() + .unwrap(), + storagescheme: inner + .next() + .unwrap() + .as_str() + .to_ascii_lowercase() + .parse::() + .unwrap(), + } +} + +// Parse shape starts here------------------------------------------------- + +/// Parse a pest structure to sparse shape information, including 3 int, which are number of rows, cols and non-zeros. +fn parse_sparse_shape( + inner: &mut Pairs<'_, Rule>, + storagescheme: &StorageScheme, +) -> Result<(usize, usize, usize), MatrixMarketError> { + // unwrap() in this function are guaranteed by parsing the data + let shape_inner = inner.next().unwrap(); + if shape_inner.as_rule() != Rule::SparseShape { + return Err(MatrixMarketError::from_kind_and_message(MatrixMarketErrorKind::ParsingError,format!(" + Shape shape line requires 3 int numbers as number of rows, columns and non-zeros, but line {} was provided here. + ",shape_inner.as_str()))); + } + + let mut inner = shape_inner.into_inner(); + + let r = inner.next().unwrap().as_str().parse::().unwrap(); + let c = inner.next().unwrap().as_str().parse::().unwrap(); + let nnz = inner.next().unwrap().as_str().parse::().unwrap(); + + // check for square matrix, when it's not a general matrix + if *storagescheme != StorageScheme::General && r != c { + return Err(MatrixMarketError::from_kind_and_message(MatrixMarketErrorKind::NonSquare, format!("(Skew-)Symmetric or hermitian matrix should be square matrix, but it has dimension {} and {}", r, c))); + } + + Ok((r, c, nnz)) +} + +/// Parse a pest structure to dense shape information, including 2 int, which are number of rows, cols. +fn parse_dense_shape( + inner: &mut Pairs<'_, Rule>, + storagescheme: &StorageScheme, +) -> Result<(usize, usize, usize), MatrixMarketError> { + // unwrap() in this function are guaranteed by parsing the data + let shape_inner = inner.next().unwrap(); + if shape_inner.as_rule() != Rule::DenseShape { + return Err(MatrixMarketError::from_kind_and_message(MatrixMarketErrorKind::ParsingError,format!(" + Shape shape line requires 2 int numbers as number of rows, columns, but line {} was provided here. + ",shape_inner.as_str()))); + } + + let mut inner = shape_inner.into_inner(); + let r = inner.next().unwrap().as_str().parse::().unwrap(); + let c = inner.next().unwrap().as_str().parse::().unwrap(); + + // check for square matrix, when it's not a general matrix + if *storagescheme != StorageScheme::General && r != c { + return Err(MatrixMarketError::from_kind_and_message(MatrixMarketErrorKind::NonSquare, format!("(Skew-)Symmetric or hermitian matrix should be square matrix, but it has dimension {} and {}", r, c))); + } + + let n: usize; + // Calculate the number of entries in the dense matrix + match storagescheme { + StorageScheme::General => { + // general matrix should contain r*c entries + n = r * c; + } + StorageScheme::Symmetric | StorageScheme::Hermitian => { + // it must be square matrix, so r==c is true here + // Symmetric or Hermitian should contain 1+2...+r = r*(r+1)/2 entries + n = r * (r + 1) / 2; + } + StorageScheme::Skew => { + // it must be square matrix, so r==c is true here + // Skew-Symmetric should contain 1+2...+r-1 = r*(r-1)/2 entries + n = r * (r - 1) / 2; + } + } + + Ok((r, c, n)) +} + +// Parse shape ends here------------------------------------------------- + +// Parse entry starts here------------------------------------------------- + +/// Parse a pest structure to sparse real entry, including 2 int, which are number of rows, cols, and a real number as data +fn parse_sparse_real(inner: &mut Pairs<'_, Rule>) -> Result<(usize, usize, T), MatrixMarketError> +where + T: MatrixMarketScalar, +{ + // unwrap() in this function are guaranteed by parsing the data + let entry_inner = inner.next().unwrap(); + if entry_inner.as_rule() != Rule::SparseReal { + return Err(MatrixMarketError::from_kind_and_message(MatrixMarketErrorKind::ParsingError,format!(" + Spare real matrix requires 2 int number as coordinates and 1 real number as data, but line {} was provided. + ",entry_inner.as_str() ))); + } + + let mut inner = entry_inner.into_inner(); + let (r, c) = parse_sparse_coordinate(&mut inner)?; + let d = inner.next().unwrap().as_str().parse::().unwrap(); + Ok((r, c, T::from_f64(d)?)) +} + +/// Parse a pest structure to sparse integer entry, including 2 int, which are number of rows, cols, and a int number as data +fn parse_sparse_int(inner: &mut Pairs<'_, Rule>) -> Result<(usize, usize, T), MatrixMarketError> +where + T: MatrixMarketScalar, +{ + // unwrap() in this function are guaranteed by parsing the data + let entry_inner = inner.next().unwrap(); + // Because integer numbers can also be parsed as float numbers, it will be checked again in `parse::()?` + if entry_inner.as_rule() != Rule::SparseReal { + return Err(MatrixMarketError::from_kind_and_message( + MatrixMarketErrorKind::ParsingError, + format!( + " + Spare real matrix requires 3 int number as coordinates and data, but line {} was provided. + ", + entry_inner.as_str() + ), + )); + } + let mut inner = entry_inner.into_inner(); + let (r, c) = parse_sparse_coordinate(&mut inner)?; + // Here to guarantee it is an integer number + let d = inner.next().unwrap().as_str().parse::()?; + Ok((r, c, T::from_i128(d)?)) +} + +/// Parse a pest structure to sparse pattern entry, including 2 int, which are number of rows, cols +fn parse_sparse_pattern( + inner: &mut Pairs<'_, Rule>, +) -> Result<(usize, usize, T), MatrixMarketError> +where + T: MatrixMarketScalar, +{ + // unwrap() in this function are guaranteed by parsing the data + let entry_inner = inner.next().unwrap(); + if entry_inner.as_rule() != Rule::SparsePattern { + return Err(MatrixMarketError::from_kind_and_message( + MatrixMarketErrorKind::ParsingError, + format!( + " + Spare real matrix requires 2 int number as coordinates, but line {} was provided. + ", + entry_inner.as_str() + ), + )); + } + let mut inner = entry_inner.into_inner(); + let (r, c) = parse_sparse_coordinate(&mut inner)?; + Ok((r, c, T::from_pattern(())?)) +} + +/// Parse a pest structure to sparse complex entry, including 2 int, which are number of rows, cols, and 2 real number as complex data +fn parse_sparse_complex( + inner: &mut Pairs<'_, Rule>, +) -> Result<(usize, usize, T), MatrixMarketError> +where + T: MatrixMarketScalar, +{ + // unwrap() in this function are guaranteed by parsing the data + let entry_inner = inner.next().unwrap(); + if entry_inner.as_rule() != Rule::SparseComplex { + return Err(MatrixMarketError::from_kind_and_message(MatrixMarketErrorKind::ParsingError,format!(" + Spare real matrix requires 2 int number as coordinates and 2 real number as complex data, but line {} was provided. + ",entry_inner.as_str() ))); + } + let mut inner = entry_inner.into_inner(); + let (r, c) = parse_sparse_coordinate(&mut inner)?; + let real = inner.next().unwrap().as_str().parse::().unwrap(); + let imag = inner.next().unwrap().as_str().parse::().unwrap(); + let complex = Complex::::new(real, imag); + Ok((r, c, T::from_c64(complex)?)) +} + +/// Parse a pest structure to dense real entry, including a real number as data +fn parse_dense_real(inner: &mut Pairs<'_, Rule>) -> Result +where + T: MatrixMarketScalar, +{ + // unwrap() in this function are guaranteed by parsing the data + let entry_inner = inner.next().unwrap(); + if entry_inner.as_rule() != Rule::DenseReal { + return Err(MatrixMarketError::from_kind_and_message( + MatrixMarketErrorKind::ParsingError, + format!( + " + Dense real matrix requires 1 real number as data, but line {} was provided. + ", + entry_inner.as_str() + ), + )); + } + let mut inner = entry_inner.into_inner(); + let d = inner.next().unwrap().as_str().parse::().unwrap(); + Ok(T::from_f64(d)?) +} + +/// Parse a pest structure to dense integer entry, including a integer number as data +fn parse_dense_int(inner: &mut Pairs<'_, Rule>) -> Result +where + T: MatrixMarketScalar, +{ + // unwrap() in this function are guaranteed by parsing the data + let entry_inner = inner.next().unwrap(); + // Because integer numbers can also be parsed as float numbers, it will be checked again in `parse::()?` + if entry_inner.as_rule() != Rule::DenseReal { + return Err(MatrixMarketError::from_kind_and_message( + MatrixMarketErrorKind::ParsingError, + format!( + " + Dense real matrix requires 1 int number as data, but line {} was provided. + ", + entry_inner.as_str() + ), + )); + } + let mut inner = entry_inner.into_inner(); + // Here to guarantee it is an integer number + let d = inner.next().unwrap().as_str().parse::()?; + Ok(T::from_i128(d)?) +} + +/// Parse a pest structure to dense complex entry, including 2 real number as complex data +fn parse_dense_complex(inner: &mut Pairs<'_, Rule>) -> Result +where + T: MatrixMarketScalar, +{ + // unwrap() in this function are guaranteed by parsing the data + let entry_inner = inner.next().unwrap(); + // Note: theoretically, 2 positive integers could also become the complex number, + // but it would be parsed as SparsePattern, because SparsePattern has higher priority. + // But DenseComplex can't have higher priority, + // because, it's more often to deal with "normal" SparsePattern, rather than "unnormal" DenseComplex + if entry_inner.as_rule() != Rule::DenseComplex && entry_inner.as_rule() != Rule::SparsePattern { + return Err(MatrixMarketError::from_kind_and_message( + MatrixMarketErrorKind::ParsingError, + format!( + " + Dense real matrix requires 2 real number as complex data, but line {} was provided. + ", + entry_inner.as_str() + ), + )); + } + let mut inner = entry_inner.into_inner(); + let real = inner.next().unwrap().as_str().parse::().unwrap(); + let imag = inner.next().unwrap().as_str().parse::().unwrap(); + let complex = Complex::::new(real, imag); + Ok(T::from_c64(complex)?) +} + +// Parse entry ends here------------------------------------------------- + +/// Parse the coordinates information used for sparse matrix +fn parse_sparse_coordinate( + inner: &mut Pairs<'_, Rule>, +) -> Result<(usize, usize), MatrixMarketError> { + // unwrap() in this function are guaranteed by parsing the data + let r = inner.next().unwrap().as_str().parse::().unwrap(); + let c = inner.next().unwrap().as_str().parse::().unwrap(); + if r * c == 0 { + return Err(MatrixMarketError::from_kind_and_message( + MatrixMarketErrorKind::ZeroError, + String::from("The data has to be one-indexed"), + )); + } + // The coordinates in matrix market is one-based, but in CooMatrix is zero-based. + Ok((r - 1, c - 1)) +} + +/// Calculate the next coordinates used for dense matrix +fn next_dense_coordinate( + current_dense_coordinate: &mut (usize, usize), + shape: (usize, usize, usize), + storagescheme: &StorageScheme, +) { + // matrix market is column based format. + // so it follows the order (0,0) -> (1,0) -> ... -> (row, 0) -> (0,1) -> ... ->(row,col) + // current_dense_coordinate is (row, column) + match storagescheme { + StorageScheme::General => { + if current_dense_coordinate.0 < shape.0 - 1 { + current_dense_coordinate.0 += 1 + } else { + // jump to next column, reset row to 1, column add 1 + current_dense_coordinate.0 = 0; + current_dense_coordinate.1 += 1; + } + } + StorageScheme::Symmetric | StorageScheme::Hermitian => { + if current_dense_coordinate.0 < shape.0 - 1 { + current_dense_coordinate.0 += 1 + } else { + // jump to next column, column add 1, then set row equals to current column + // for example (0,0) -> (1,0) -> ... -> (row, 0) -> (1,1) -> ... + current_dense_coordinate.1 += 1; + current_dense_coordinate.0 = current_dense_coordinate.1; + } + } + StorageScheme::Skew => { + if current_dense_coordinate.0 < shape.0 - 1 { + current_dense_coordinate.0 += 1; + } else { + // jump to next column, set row equals to current column, then column add 1 + // skew matrix doesn't have element on diagonal + // for example (1,0) -> (2,0) -> ... -> (row, 0) -> (2,1) -> ... + current_dense_coordinate.1 += 1; + current_dense_coordinate.0 = current_dense_coordinate.1 + 1; + } + } + } +} diff --git a/nalgebra-sparse/src/io/mod.rs b/nalgebra-sparse/src/io/mod.rs new file mode 100644 index 00000000..89b21ffb --- /dev/null +++ b/nalgebra-sparse/src/io/mod.rs @@ -0,0 +1,38 @@ +//! Functionality for importing and exporting sparse matrices to and from files. +//! +//! **Available only when the `io` feature is enabled.** +//! +//! The following formats are currently supported: +//! +//! | Format | Import | Export | +//! | ------------------------------------------------|------------|------------| +//! | [Matrix market](#matrix-market-format) | Yes | No | +//! +//! [Matrix market]: https://math.nist.gov/MatrixMarket/formats.html +//! +//! ## Matrix Market format +//! +//! The Matrix Market format is a simple ASCII-based file format for sparse matrices, and was initially developed for +//! the [NIST Matrix Market](https://math.nist.gov/MatrixMarket/), a repository of example sparse matrices. +//! In later years it has largely been superseded by the +//! [SuiteSparse Matrix Collection](https://sparse.tamu.edu/) (formerly University of Florida Sparse Matrix Collection), +//! which also uses the Matrix Market file format. +//! +//! We currently offer functionality for importing a Matrix market file to an instance of a +//! [CooMatrix](crate::CooMatrix) through the function [load_coo_from_matrix_market_file]. It is also possible to load +//! a matrix stored in the matrix market format with the function [load_coo_from_matrix_market_str]. +//! +//! Export is currently not implemented, but [planned](https://github.com/dimforge/nalgebra/issues/1037). +//! +//! Our implementation is based on the [format description](https://math.nist.gov/MatrixMarket/formats.html) +//! on the Matrix Market website and the +//! [following NIST whitepaper](https://math.nist.gov/MatrixMarket/reports/MMformat.ps): +//! +//! > Boisvert, Ronald F., Roldan Pozo, and Karin A. Remington.
+//! > "*The Matrix Market Exchange Formats: Initial Design.*" (1996). + +pub use self::matrix_market::{ + load_coo_from_matrix_market_file, load_coo_from_matrix_market_str, MatrixMarketError, + MatrixMarketErrorKind, MatrixMarketScalar, +}; +mod matrix_market; diff --git a/nalgebra-sparse/src/lib.rs b/nalgebra-sparse/src/lib.rs index bf845757..edbf83bd 100644 --- a/nalgebra-sparse/src/lib.rs +++ b/nalgebra-sparse/src/lib.rs @@ -19,6 +19,7 @@ //! - Sparsity patterns in CSR and CSC matrices are explicitly represented by the //! [SparsityPattern](pattern::SparsityPattern) type, which encodes the invariants of the //! associated index data structures. +//! - [Matrix market format support](`io`) when the `io` feature is enabled. //! - [proptest strategies](`proptest`) for sparse matrices when the feature //! `proptest-support` is enabled. //! - [matrixcompare support](https://crates.io/crates/matrixcompare) for effortless @@ -142,11 +143,19 @@ )] pub extern crate nalgebra as na; +#[cfg(feature = "io")] +extern crate pest; +#[macro_use] +#[cfg(feature = "io")] +extern crate pest_derive; + pub mod convert; pub mod coo; pub mod csc; pub mod csr; pub mod factorization; +#[cfg(feature = "io")] +pub mod io; pub mod ops; pub mod pattern; diff --git a/nalgebra-sparse/tests/unit.rs b/nalgebra-sparse/tests/unit.rs index 73a95cd7..74d32a40 100644 --- a/nalgebra-sparse/tests/unit.rs +++ b/nalgebra-sparse/tests/unit.rs @@ -1,6 +1,9 @@ //! Unit tests -#[cfg(any(not(feature = "proptest-support"), not(feature = "compare")))] -compile_error!("Tests must be run with features `proptest-support` and `compare`"); +#[cfg(not(all(feature = "proptest-support", feature = "compare", feature = "io",)))] +compile_error!( + "Please enable the `proptest-support`, `compare` and `io` features in order to compile and run the tests. + Example: `cargo test -p nalgebra-sparse --features proptest-support,compare,io`" +); mod unit_tests; diff --git a/nalgebra-sparse/tests/unit_tests/matrix_market.rs b/nalgebra-sparse/tests/unit_tests/matrix_market.rs new file mode 100644 index 00000000..48ff1a78 --- /dev/null +++ b/nalgebra-sparse/tests/unit_tests/matrix_market.rs @@ -0,0 +1,345 @@ +use matrixcompare::assert_matrix_eq; +use nalgebra::dmatrix; +use nalgebra::Complex; +use nalgebra_sparse::io::load_coo_from_matrix_market_str; +use nalgebra_sparse::CooMatrix; + +#[test] +#[rustfmt::skip] +fn test_matrixmarket_sparse_real_general_empty() { + // Test several valid zero-shapes of a sparse matrix + let shapes = vec![ (0, 0), (1, 0), (0, 1) ]; + let strings: Vec = shapes + .iter() + .map(|(m, n)| format!("%%MatrixMarket matrix coordinate real general\n {} {} 0", m, n)) + .collect(); + + for (shape,string) in shapes.iter().zip(strings.iter()) { + let sparse_mat = load_coo_from_matrix_market_str::(string).unwrap(); + assert_eq!(sparse_mat.nrows(), shape.0); + assert_eq!(sparse_mat.ncols(), shape.1); + assert_eq!(sparse_mat.nnz(), 0); + } +} + +#[test] +#[rustfmt::skip] +fn test_matrixmarket_dense_real_general_empty() { + // Test several valid zero-shapes of a dense matrix + let shapes = vec![ (0, 0), (1, 0), (0, 1) ]; + let strings: Vec = shapes + .iter() + .map(|(m, n)| format!("%%MatrixMarket matrix array real general\n {} {}", m, n)) + .collect(); + + for (shape,string) in shapes.iter().zip(strings.iter()) { + let sparse_mat = load_coo_from_matrix_market_str::(string).unwrap(); + assert_eq!(sparse_mat.nrows(), shape.0); + assert_eq!(sparse_mat.ncols(), shape.1); + assert_eq!(sparse_mat.nnz(), 0); + } +} + +#[test] +#[rustfmt::skip] +fn test_matrixmarket_sparse_real_general() { + let file_str = r#" +%%MatrixMarket matrix CoOrdinate real general +% This is also an example of free-format features. +%================================================================================= +% +% This ASCII file represents a sparse MxN matrix with L +% nonzeros in the following Matrix Market format: +% +% +----------------------------------------------+ +% |%%MatrixMarket matrix coordinate real general | <--- header line +% |% | <--+ +% |% comments | |-- 0 or more comment lines +% |% | <--+ +% | M T L | <--- rows, columns, entries +% | I1 J1 A(I1, J1) | <--+ +% | I2 J2 A(I2, J2) | | +% | I3 J3 A(I3, J3) | |-- L lines +% | . . . | | +% | IL JL A(IL, JL) | <--+ +% +----------------------------------------------+ +% +% Indices are 1-based, i.e. A(1,1) is the first element. +% +%================================================================================= + 5 5 8 + 1 1 1 + + 2 2 1.050e+01 + 3 3 1.500e-02 + + + + + + + + + 1 4 6.000e+00 + + 4 2 2.505e+02 + +4 4 -2.800e+02 +4 5 3.332e+01 + 5 5 1.200e+01 +"#; + let sparse_mat = load_coo_from_matrix_market_str::(file_str).unwrap(); + let expected = dmatrix![ + 1.0, 0.0, 0.0, 6.0, 0.0; + 0.0, 10.5, 0.0, 0.0, 0.0; + 0.0, 0.0, 0.015, 0.0, 0.0; + 0.0, 250.5, 0.0, -280.0, 33.32; + 0.0, 0.0, 0.0, 0.0, 12.0; + ]; + assert_matrix_eq!(sparse_mat, expected); +} + +#[test] +#[rustfmt::skip] +fn test_matrixmarket_sparse_int_symmetric() { + let file_str = r#" +%%MatrixMarket matrix coordinate integer symmetric +% + 5 5 9 + 1 1 11 + 2 2 22 + 3 2 23 + 3 3 33 + 4 2 24 + 4 4 44 + 5 1 -15 + 5 3 35 + 5 5 55 +"#; + let sparse_mat = load_coo_from_matrix_market_str::(file_str).unwrap(); + let expected = dmatrix![ + 11, 0, 0, 0, -15; + 0, 22, 23, 24, 0; + 0, 23, 33, 0, 35; + 0, 24, 0, 44, 0; + -15, 0, 35, 0, 55; + ]; + assert_matrix_eq!(sparse_mat, expected); +} + +#[test] +#[rustfmt::skip] +fn test_matrixmarket_sparse_complex_hermitian() { + let file_str = r#" +%%MatrixMarket matrix coordinate complex hermitian +% + 5 5 7 + 1 1 1.0 0.0 + 2 2 10.5 0.0 + 4 2 250.5 22.22 + 3 3 0.015 0.0 + 4 4 -2.8e2 0.0 + 5 5 12.0 0.0 + 5 4 0.0 33.32 + +"#; + let sparse_mat = load_coo_from_matrix_market_str::>(file_str).unwrap(); + let expected = dmatrix![ + Complex::{re:1.0,im:0.0}, Complex::{re:0.0,im:0.0}, Complex::{re:0.0,im:0.0}, Complex::{re:0.0,im:0.0},Complex::{re:0.0,im:0.0}; + Complex::{re:0.0,im:0.0}, Complex::{re:10.5,im:0.0}, Complex::{re:0.0,im:0.0}, Complex::{re:250.5,im:-22.22},Complex::{re:0.0,im:0.0}; + Complex::{re:0.0,im:0.0}, Complex::{re:0.0,im:0.0}, Complex::{re:0.015,im:0.0}, Complex::{re:0.0,im:0.0},Complex::{re:0.0,im:0.0}; + Complex::{re:0.0,im:0.0}, Complex::{re:250.5,im:22.22}, Complex::{re:0.0,im:0.0}, Complex::{re:-280.0,im:0.0},Complex::{re:0.0,im:-33.32}; + Complex::{re:0.0,im:0.0}, Complex::{re:0.0,im:0.0}, Complex::{re:0.0,im:0.0}, Complex::{re:0.0,im:33.32},Complex::{re:12.0,im:0.0}; + ]; + assert_matrix_eq!(sparse_mat, expected); +} + +#[test] +#[rustfmt::skip] +fn test_matrixmarket_sparse_real_skew() { + let file_str = r#" +%%MatrixMarket matrix coordinate real skew-symmetric +% + 5 5 4 + 3 2 -23.0 + 4 2 -24.0 + 5 1 -15.0 + 5 3 -35.0 +"#; + let sparse_mat = load_coo_from_matrix_market_str::(file_str).unwrap(); + let expected = dmatrix![ + 0.0, 0.0, 0.0, 0.0, 15.0; + 0.0, 0.0, 23.0, 24.0, 0.0; + 0.0, -23.0, 0.0, 0.0, 35.0; + 0.0, -24.0, 0.0, 0.0, 0.0; + -15.0, 0.0, -35.0, 0.0, 0.0; + ]; + assert_matrix_eq!(sparse_mat, expected); +} + +#[test] +#[rustfmt::skip] +fn test_matrixmarket_sparse_pattern_general() { + let file_str = r#" +%%MatrixMarket matrix coordinate pattern general +% + 5 5 10 + 1 1 + 1 5 + 2 3 + 2 4 + 3 2 + 3 5 + 4 1 + 5 2 + 5 4 + 5 5 +"#; + let pattern_matrix = load_coo_from_matrix_market_str::<()>(file_str).unwrap(); + let nrows = pattern_matrix.nrows(); + let ncols = pattern_matrix.ncols(); + let (row_idx, col_idx, val) = pattern_matrix.disassemble(); + let values = vec![1; val.len()]; + let sparse_mat = CooMatrix::try_from_triplets(nrows, ncols, row_idx, col_idx, values).unwrap(); + let expected = dmatrix![ + 1, 0, 0, 0, 1; + 0, 0, 1, 1, 0; + 0, 1, 0, 0, 1; + 1, 0, 0, 0, 0; + 0, 1, 0, 1, 1; + ]; + assert_matrix_eq!(sparse_mat, expected); +} + +#[test] +#[rustfmt::skip] +fn test_matrixmarket_dense_real_general() { + let file_str = r#" +%%MatrixMarket matrix array real general +% +4 3 +1.0 +2.0 +3.0 +4.0 +5.0 +6.0 +7.0 +8.0 +9.0 +10.0 +11.0 +12.0 + +"#; + let sparse_mat = load_coo_from_matrix_market_str::(file_str).unwrap(); + let expected = dmatrix![ + 1.0, 5.0, 9.0; + 2.0, 6.0, 10.0; + 3.0, 7.0, 11.0; + 4.0, 8.0, 12.0; + ]; + assert_matrix_eq!(sparse_mat, expected); +} + +#[test] +#[rustfmt::skip] +fn test_matrixmarket_dense_real_symmetric() { + let file_str = r#" +%%MatrixMarket matrix array real symmetric +% +4 4 +1.0 +2.0 +3.0 +4.0 +5.0 +6.0 +7.0 +8.0 +9.0 +10.0 + +"#; + let sparse_mat = load_coo_from_matrix_market_str::(file_str).unwrap(); + let expected = dmatrix![ + 1.0, 2.0, 3.0, 4.0; + 2.0, 5.0, 6.0, 7.0; + 3.0, 6.0, 8.0, 9.0; + 4.0, 7.0, 9.0, 10.0; + ]; + assert_matrix_eq!(sparse_mat, expected); +} + +#[test] +#[rustfmt::skip] +fn test_matrixmarket_dense_complex_hermitian() { + let file_str = r#" +%%MatrixMarket matrix array complex hermitian +% +4 4 +1.0 0.0 +2.0 2.0 +3.0 3.0 +4.0 4.0 +5.0 0.0 +6.0 6.0 +7.0 7.0 +8.0 0.0 +9.0 9.0 +10.0 0.0 + +"#; + let sparse_mat = load_coo_from_matrix_market_str::>(file_str).unwrap(); + let expected = dmatrix![ + Complex::{re:1.0,im:0.0}, Complex::{re:2.0,im:-2.0} ,Complex::{re:3.0,im:-3.0} ,Complex::{re:4.0,im:-4.0}; + Complex::{re:2.0,im:2.0}, Complex::{re:5.0,im:0.0} ,Complex::{re:6.0,im:-6.0} ,Complex::{re:7.0,im:-7.0}; + Complex::{re:3.0,im:3.0}, Complex::{re:6.0,im:6.0} ,Complex::{re:8.0,im:0.0} ,Complex::{re:9.0,im:-9.0}; + Complex::{re:4.0,im:4.0}, Complex::{re:7.0,im:7.0} ,Complex::{re:9.0,im:9.0} ,Complex::{re:10.0,im:0.0}; + ]; + assert_matrix_eq!(sparse_mat, expected); +} + +#[test] +#[rustfmt::skip] +fn test_matrixmarket_dense_int_skew() { + let file_str = r#" +%%MatrixMarket matrix array integer skew-symmetric +% +4 4 +1 +2 +3 +4 +5 +6 +"#; + let sparse_mat = load_coo_from_matrix_market_str::(file_str).unwrap(); + let expected = dmatrix![ + 0,-1,-2,-3; + 1, 0,-4,-5; + 2, 4, 0,-6; + 3, 5, 6, 0; + ]; + assert_matrix_eq!(sparse_mat, expected); +} + +#[test] +#[rustfmt::skip] +fn test_matrixmarket_dense_complex_general() { + let file_str = r#" +%%MatrixMarket matrix array complex general +% +2 2 +1 0 +1 0 +1 0 +1 0 +"#; + let sparse_mat = load_coo_from_matrix_market_str::>(file_str).unwrap(); + let expected = dmatrix![ + Complex::{re:1.0,im:0.0},Complex::{re:1.0,im:0.0}; + Complex::{re:1.0,im:0.0},Complex::{re:1.0,im:0.0}; + ]; + assert_matrix_eq!(sparse_mat, expected); +} diff --git a/nalgebra-sparse/tests/unit_tests/mod.rs b/nalgebra-sparse/tests/unit_tests/mod.rs index 0099a246..7090a493 100644 --- a/nalgebra-sparse/tests/unit_tests/mod.rs +++ b/nalgebra-sparse/tests/unit_tests/mod.rs @@ -3,6 +3,7 @@ mod convert_serial; mod coo; mod csc; mod csr; +mod matrix_market; mod ops; mod pattern; mod proptest;