From d0cf6cb77332a1e8320f84abedfe387357e6feb4 Mon Sep 17 00:00:00 2001 From: Timo Betcke Date: Fri, 8 Mar 2024 17:46:39 +0000 Subject: [PATCH] New rlst scalar (#63) * WIP: Remove rlst_common * Removed rlst_common * Removed cyclic dependency * Renamed Scalar to RlstScalar * Added pointer to original licese for Scalar package derived RlstScalar trait. * Removed cauchy from dependencies * Fixed serde issue --- Cargo.toml | 1 - blis/Cargo.toml | 1 - blis/examples/threads.rs | 10 - blis/src/interface.rs | 27 +- blis/src/interface/gemm.rs | 21 +- blis/src/interface/types.rs | 30 +- blis/src/raw.rs | 2 - common/Cargo.toml | 22 - common/src/lib.rs | 4 - common/src/types.rs | 47 --- dense/Cargo.toml | 7 +- dense/src/array.rs | 49 ++- dense/src/array/empty_axis.rs | 22 +- dense/src/array/iterators.rs | 42 +- dense/src/array/mult_into.rs | 16 +- dense/src/array/operations.rs | 48 +-- dense/src/array/operators/addition.rs | 17 +- .../src/array/operators/cmp_wise_division.rs | 17 +- dense/src/array/operators/cmp_wise_product.rs | 17 +- dense/src/array/operators/conjugate.rs | 12 +- dense/src/array/operators/other.rs | 2 +- dense/src/array/operators/scalar_mult.rs | 14 +- dense/src/array/operators/to_complex.rs | 32 +- dense/src/array/operators/transpose.rs | 16 +- dense/src/array/random.rs | 8 +- dense/src/array/slice.rs | 27 +- dense/src/array/views.rs | 81 ++-- dense/src/base_array.rs | 29 +- dense/src/data_container.rs | 38 +- dense/src/gemm.rs | 26 ++ dense/src/gemm/blis.rs | 102 +++++ dense/src/lib.rs | 2 + dense/src/linalg/inverse.rs | 2 +- dense/src/linalg/lu.rs | 6 +- dense/src/linalg/pseudo_inverse.rs | 16 +- dense/src/linalg/qr.rs | 18 +- dense/src/linalg/svd.rs | 34 +- dense/src/matrix_multiply.rs | 20 +- dense/src/tools.rs | 20 +- dense/src/traits.rs | 14 +- dense/src/traits/accessors.rs | 18 +- dense/src/types.rs | 379 ++++++++++++++++++ io/Cargo.toml | 2 - io/src/matrix_market.rs | 20 +- lapack/Cargo.toml | 1 - operator/Cargo.toml | 2 - operator/src/interface/array_vector_space.rs | 18 +- .../src/interface/dense_matrix_operator.rs | 15 +- operator/src/interface/sparse_operator.rs | 26 +- .../src/operations/conjugate_gradients.rs | 14 +- .../src/operations/modified_gram_schmidt.rs | 8 +- operator/src/operator.rs | 2 +- operator/src/space/dual_space.rs | 2 +- operator/src/space/element.rs | 4 +- operator/src/space/linear_space.rs | 4 +- operator/src/space/normed_space.rs | 6 +- rlst/Cargo.toml | 1 - rlst/examples/fixed_size_array.rs | 2 +- rlst/tests/array_operations.rs | 2 +- sparse/Cargo.toml | 2 - sparse/src/distributed_vector.rs | 6 +- sparse/src/ghost_communicator.rs | 6 +- .../index_layout/default_mpi_index_layout.rs | 4 +- .../default_serial_index_layout.rs | 2 +- sparse/src/sparse/csc_mat.rs | 18 +- sparse/src/sparse/csr_mat.rs | 18 +- sparse/src/sparse/mpi_csr_mat.rs | 8 +- sparse/src/sparse/tools.rs | 4 +- sparse/src/sparse/umfpack.rs | 14 +- sparse/src/traits/index_layout.rs | 2 +- sparse/src/traits/indexable_vector.rs | 22 +- 71 files changed, 964 insertions(+), 587 deletions(-) delete mode 100644 blis/examples/threads.rs delete mode 100644 common/Cargo.toml delete mode 100644 common/src/lib.rs delete mode 100644 common/src/types.rs create mode 100644 dense/src/gemm.rs create mode 100644 dense/src/gemm/blis.rs create mode 100644 dense/src/types.rs diff --git a/Cargo.toml b/Cargo.toml index 0c40b483..d34a6c9d 100644 --- a/Cargo.toml +++ b/Cargo.toml @@ -5,7 +5,6 @@ members = [ "dense", "operator", "sparse", - "common", "io", "suitesparse-src", "blis-src", diff --git a/blis/Cargo.toml b/blis/Cargo.toml index 8b1ae173..2522b8bf 100644 --- a/blis/Cargo.toml +++ b/blis/Cargo.toml @@ -12,7 +12,6 @@ edition = "2021" rlst-blis-src = { path = "../blis-src" } num_cpus = "1" num = "0.4" -cauchy = "0.4" [lib] name = "rlst_blis" diff --git a/blis/examples/threads.rs b/blis/examples/threads.rs deleted file mode 100644 index 9ff9700d..00000000 --- a/blis/examples/threads.rs +++ /dev/null @@ -1,10 +0,0 @@ -//! Modify the threading behaviour in BLIS - -use rlst_blis::interface::threading; - -fn main() { - println!("Num threads: {}", threading::get_num_threads()); - println!("Setting threads to maximum cpu threads."); - threading::enable_threading(); - println!("Num threads: {}", threading::get_num_threads()); -} diff --git a/blis/src/interface.rs b/blis/src/interface.rs index 96446c59..93dd3a34 100644 --- a/blis/src/interface.rs +++ b/blis/src/interface.rs @@ -1,29 +1,4 @@ //! Public interface to Blis routines. -pub mod gemm; +// pub mod gemm; pub mod threading; -pub mod types; - -use cauchy::Scalar; - -/// Compute expected size of a data slice from stride and shape. -pub fn get_expected_data_size(stride: [usize; 2], shape: [usize; 2]) -> usize { - if shape[0] == 0 || shape[1] == 0 { - return 0; - } - - 1 + (shape[0] - 1) * stride[0] + (shape[1] - 1) * stride[1] -} - -/// Panic if expected data size is not identical to actual data size. -pub fn assert_data_size(data: &[T], stride: [usize; 2], shape: [usize; 2]) { - let expected = get_expected_data_size(stride, shape); - - assert_eq!( - expected, - data.len(), - "Wrong size for data slice. Actual size {}. Expected size {}.", - data.len(), - expected - ); -} diff --git a/blis/src/interface/gemm.rs b/blis/src/interface/gemm.rs index 21769ceb..b52f6f66 100644 --- a/blis/src/interface/gemm.rs +++ b/blis/src/interface/gemm.rs @@ -2,10 +2,23 @@ use super::types::TransMode; use crate::interface::assert_data_size; use crate::raw; -use cauchy::{c32, c64, Scalar}; +use cauchy::{c32, c64}; + +/// Transposition Mode. +#[derive(Clone, Copy, PartialEq)] +#[repr(u32)] + +pub enum TransMode { + /// Complex conjugate of matrix. + ConjNoTrans, + /// No modification of matrix. + NoTrans, + /// Transposition of matrix. + Trans, + /// Conjugate transpose of matrix. + ConjTrans, +} -/// Safe interface to using the Blis Gemm routine. -/// /// Performs the matrix multiplication `C = alpha * A * B + beta * C`. /// # Arguments /// @@ -25,7 +38,7 @@ use cauchy::{c32, c64, Scalar}; /// - `c` - Reference to data of C. /// - `rsc` - Row stride of C. /// - `csc` - Column stride of C. -pub trait Gemm: Scalar { +pub trait Gemm: Sized { #[allow(clippy::too_many_arguments)] fn gemm( transa: TransMode, diff --git a/blis/src/interface/types.rs b/blis/src/interface/types.rs index cfe3849f..ec2dabdc 100644 --- a/blis/src/interface/types.rs +++ b/blis/src/interface/types.rs @@ -1,16 +1,16 @@ -//! Interface to Blis types -use crate::raw; +// //! Interface to Blis types +// use crate::raw; -/// Transposition Mode. -#[derive(Clone, Copy, PartialEq)] -#[repr(u32)] -pub enum TransMode { - /// Complex conjugate of matrix. - ConjNoTrans = raw::trans_t_BLIS_CONJ_NO_TRANSPOSE, - /// No modification of matrix. - NoTrans = raw::trans_t_BLIS_NO_TRANSPOSE, - /// Transposition of matrix. - Trans = raw::trans_t_BLIS_TRANSPOSE, - /// Conjugate transpose of matrix. - ConjTrans = raw::trans_t_BLIS_CONJ_TRANSPOSE, -} +// /// Transposition Mode. +// #[derive(Clone, Copy, PartialEq)] +// #[repr(u32)] +// pub enum TransMode { +// /// Complex conjugate of matrix. +// ConjNoTrans = raw::trans_t_BLIS_CONJ_NO_TRANSPOSE, +// /// No modification of matrix. +// NoTrans = raw::trans_t_BLIS_NO_TRANSPOSE, +// /// Transposition of matrix. +// Trans = raw::trans_t_BLIS_TRANSPOSE, +// /// Conjugate transpose of matrix. +// ConjTrans = raw::trans_t_BLIS_CONJ_TRANSPOSE, +// } diff --git a/blis/src/raw.rs b/blis/src/raw.rs index 663a8042..ab5174be 100644 --- a/blis/src/raw.rs +++ b/blis/src/raw.rs @@ -4,6 +4,4 @@ #![allow(non_camel_case_types)] #![allow(non_snake_case)] -extern crate rlst_blis_src; - include!(concat!(env!("OUT_DIR"), "/bindings.rs")); diff --git a/common/Cargo.toml b/common/Cargo.toml deleted file mode 100644 index cd785191..00000000 --- a/common/Cargo.toml +++ /dev/null @@ -1,22 +0,0 @@ -[features] -strict = [] - -[package] -name = "rlst-common" -version = "0.1.0" -edition = "2021" - -# See more keys and their definitions at https://doc.rust-lang.org/cargo/reference/manifest.html - -[dependencies] -num = "0.4" -cauchy = "0.4" -thiserror = "1.0.40" -rand = "0.8" -rand_distr = "0.4" -approx = { version = "0.5", features = ["num-complex"] } -rlst-blis = { path = "../blis" } - - -[lib] -name = "rlst_common" diff --git a/common/src/lib.rs b/common/src/lib.rs deleted file mode 100644 index 2d5cb50d..00000000 --- a/common/src/lib.rs +++ /dev/null @@ -1,4 +0,0 @@ -//! Common Rlst data structures -#![cfg_attr(feature = "strict", deny(warnings))] - -pub mod types; diff --git a/common/src/types.rs b/common/src/types.rs deleted file mode 100644 index 2b65e0e7..00000000 --- a/common/src/types.rs +++ /dev/null @@ -1,47 +0,0 @@ -//! Basic types. - -pub use cauchy::{c32, c64, Scalar}; -use thiserror::Error; - -/// The Rlst error type. -#[derive(Error, Debug)] -pub enum RlstError { - #[error("Method {0} is not implemented.")] - NotImplemented(String), - #[error("Operation {0} failed.")] - OperationFailed(String), - #[error("Matrix has empty dimension {0:#?}.")] - MatrixIsEmpty((usize, usize)), - #[error("Dimension mismatch. Expected {expected:}. Actual {actual:}")] - SingleDimensionError { expected: usize, actual: usize }, - #[error("Index Layout error: {0}")] - IndexLayoutError(String), - #[error("MPI Rank does not exist. {0}")] - MpiRankError(i32), - #[error("Incompatible stride for Lapack.")] - IncompatibleStride, - #[error("Lapack error: {0}")] - LapackError(i32), - #[error("{0}")] - GeneralError(String), - #[error("I/O Error: {0}")] - IoError(String), - #[error("Umfpack Error Code: {0}")] - UmfpackError(i32), - #[error("Matrix is not square. Dimension: {0}x{1}")] - MatrixNotSquare(usize, usize), - #[error("Matrix is not Hermitian (complex conjugate symmetric).")] - MatrixNotHermitian, -} - -/// Alias for an Rlst Result type. -pub type RlstResult = std::result::Result; - -/// Data chunk of fixed size N. -/// The field `valid_entries` stores how many entries of the chunk -/// contain valid data. -pub struct DataChunk { - pub data: [Item; N], - pub start_index: usize, - pub valid_entries: usize, -} diff --git a/dense/Cargo.toml b/dense/Cargo.toml index 3e2152e4..dd970518 100644 --- a/dense/Cargo.toml +++ b/dense/Cargo.toml @@ -18,19 +18,20 @@ categories = ["mathematics", "science"] name = "rlst_dense" [dependencies] -num = "0.4" -cauchy = "0.4" +num = {version = "0.4", features = ["serde", "rand"] } rand = "0.8" itertools = "0.12" rand_distr = "0.4" rlst-blis = { path = "../blis" } approx = { version = "0.5", features = ["num-complex"] } -rlst-common = { path = "../common" } paste = "1" rand_chacha = "0.3" rlst-blis-src = { path = "../blis-src" } rlst-netlib-lapack-src = { path = "../netlib-lapack-src" } lapack = "0.19.*" +thiserror = "1.0.40" +serde = "1.*" + [dev-dependencies] criterion = { version = "0.3", features = ["html_reports"] } diff --git a/dense/src/array.rs b/dense/src/array.rs index a5fe6f1d..62379617 100644 --- a/dense/src/array.rs +++ b/dense/src/array.rs @@ -9,8 +9,8 @@ use crate::data_container::SliceContainer; use crate::data_container::SliceContainerMut; use crate::data_container::VectorContainer; use crate::traits::*; -use rlst_common::types::DataChunk; -use rlst_common::types::Scalar; +use crate::types::DataChunk; +use crate::types::RlstScalar; pub mod empty_axis; pub mod iterators; @@ -36,11 +36,11 @@ pub type SliceArrayMut<'a, Item, const NDIM: usize> = /// The basic tuple type defining an array. pub struct Array(ArrayImpl) where - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape; impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const NDIM: usize, > Array @@ -57,7 +57,7 @@ impl< } /// Create a new heap allocated array from a given shape. -impl DynamicArray { +impl DynamicArray { pub fn from_shape(shape: [usize; NDIM]) -> Self { let size = shape.iter().product(); Self::new(BaseArray::new(VectorContainer::new(size), shape)) @@ -65,7 +65,7 @@ impl DynamicArray { } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const NDIM: usize, > Shape for Array @@ -76,7 +76,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const NDIM: usize, > UnsafeRandomAccessByValue for Array @@ -89,7 +89,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + ChunkedAccess, const NDIM: usize, const N: usize, @@ -97,16 +97,13 @@ impl< { type Item = Item; #[inline] - fn get_chunk( - &self, - chunk_index: usize, - ) -> Option> { + fn get_chunk(&self, chunk_index: usize) -> Option> { self.0.get_chunk(chunk_index) } } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + UnsafeRandomAccessByRef, @@ -121,7 +118,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + UnsafeRandomAccessMut, @@ -136,7 +133,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + UnsafeRandomAccessByRef, @@ -151,7 +148,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + UnsafeRandomAccessByRef @@ -166,7 +163,7 @@ impl< } /// Create an empty chunk. -pub(crate) fn empty_chunk( +pub(crate) fn empty_chunk( chunk_index: usize, nelements: usize, ) -> Option> { @@ -188,7 +185,7 @@ pub(crate) fn empty_chunk( } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + RawAccess, const NDIM: usize, > RawAccess for Array @@ -201,7 +198,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + RawAccessMut, const NDIM: usize, > RawAccessMut for Array @@ -212,7 +209,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const NDIM: usize, > NumberOfElements for Array @@ -223,7 +220,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + Stride, const NDIM: usize, > Stride for Array @@ -234,7 +231,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + ResizeInPlace, const NDIM: usize, > ResizeInPlace for Array @@ -248,13 +245,13 @@ impl< /// /// Empty arrays serve as convenient containers for input into functions that /// resize an array before filling it with data. -pub fn empty_array() -> DynamicArray { +pub fn empty_array() -> DynamicArray { let shape = [0; NDIM]; let container = VectorContainer::new(0); Array::new(BaseArray::new(container, shape)) } -impl<'a, Item: Scalar, const NDIM: usize> SliceArray<'a, Item, NDIM> { +impl<'a, Item: RlstScalar, const NDIM: usize> SliceArray<'a, Item, NDIM> { /// Create a new array from a slice with a given `shape`. /// /// The `stride` is automatically assumed to be column major. @@ -276,7 +273,7 @@ impl<'a, Item: Scalar, const NDIM: usize> SliceArray<'a, Item, NDIM> { } } -impl<'a, Item: Scalar, const NDIM: usize> SliceArrayMut<'a, Item, NDIM> { +impl<'a, Item: RlstScalar, const NDIM: usize> SliceArrayMut<'a, Item, NDIM> { /// Create a new array from a slice with a given `shape`. /// /// The `stride` is automatically assumed to be column major. @@ -299,7 +296,7 @@ impl<'a, Item: Scalar, const NDIM: usize> SliceArrayMut<'a, Item, NDIM> { } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const NDIM: usize, > std::fmt::Debug for Array diff --git a/dense/src/array/empty_axis.rs b/dense/src/array/empty_axis.rs index fcf4c3e1..1737b37e 100644 --- a/dense/src/array/empty_axis.rs +++ b/dense/src/array/empty_axis.rs @@ -15,7 +15,7 @@ pub enum AxisPosition { /// Array implementation that provides an appended empty axis. pub struct ArrayAppendAxis< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const ADIM: usize, const NDIM: usize, @@ -27,7 +27,7 @@ pub struct ArrayAppendAxis< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const ADIM: usize, const NDIM: usize, @@ -41,7 +41,7 @@ where } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const ADIM: usize, const NDIM: usize, @@ -58,7 +58,7 @@ where } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + UnsafeRandomAccessByRef, @@ -77,7 +77,7 @@ where } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + UnsafeRandomAccessMut, @@ -96,7 +96,7 @@ where } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const ADIM: usize, const NDIM: usize, @@ -117,7 +117,7 @@ where } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + RawAccess @@ -135,7 +135,7 @@ where } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + RawAccessMut @@ -152,7 +152,7 @@ where } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + Stride, const ADIM: usize, const NDIM: usize, @@ -179,7 +179,7 @@ where } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const ADIM: usize, > Array @@ -196,7 +196,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + Stride diff --git a/dense/src/array/iterators.rs b/dense/src/array/iterators.rs index a881a0bb..60bf827a 100644 --- a/dense/src/array/iterators.rs +++ b/dense/src/array/iterators.rs @@ -2,7 +2,7 @@ use crate::array::*; use crate::layout::convert_1d_nd_from_shape; -use rlst_common::types::Scalar; +use crate::types::RlstScalar; use super::slice::ArraySlice; @@ -11,7 +11,7 @@ use super::slice::ArraySlice; /// This iterator returns elements of an array in standard column major order. pub struct ArrayDefaultIterator< 'a, - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const NDIM: usize, > { @@ -24,7 +24,7 @@ pub struct ArrayDefaultIterator< /// Mutable default iterator. Like [ArrayDefaultIterator] but with mutable access. pub struct ArrayDefaultIteratorMut< 'a, - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + UnsafeRandomAccessMut, @@ -70,7 +70,7 @@ impl, const NDIM: usize> AsMultiIndex + Shape, const NDIM: usize, > ArrayDefaultIterator<'a, Item, ArrayImpl, NDIM> @@ -87,7 +87,7 @@ impl< impl< 'a, - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + UnsafeRandomAccessMut, @@ -107,7 +107,7 @@ impl< impl< 'a, - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const NDIM: usize, > std::iter::Iterator for ArrayDefaultIterator<'a, Item, ArrayImpl, NDIM> @@ -132,7 +132,7 @@ impl< impl< 'a, - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + UnsafeRandomAccessMut + Shape, @@ -155,7 +155,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const NDIM: usize, > crate::traits::DefaultIterator for Array @@ -169,7 +169,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + UnsafeRandomAccessMut, @@ -186,7 +186,7 @@ impl< pub struct RowIterator< 'a, - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const NDIM: usize, > { @@ -197,7 +197,7 @@ pub struct RowIterator< pub struct RowIteratorMut< 'a, - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + UnsafeRandomAccessMut, @@ -208,7 +208,7 @@ pub struct RowIteratorMut< current_row: usize, } -impl<'a, Item: Scalar, ArrayImpl: UnsafeRandomAccessByValue<2, Item = Item> + Shape<2>> +impl<'a, Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue<2, Item = Item> + Shape<2>> std::iter::Iterator for RowIterator<'a, Item, ArrayImpl, 2> { type Item = Array, 2, 1>, 1>; @@ -224,7 +224,7 @@ impl<'a, Item: Scalar, ArrayImpl: UnsafeRandomAccessByValue<2, Item = Item> + Sh impl< 'a, - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue<2, Item = Item> + UnsafeRandomAccessMut<2, Item = Item> + Shape<2>, @@ -246,7 +246,7 @@ impl< } } -impl + Shape<2>> +impl + Shape<2>> Array { /// Return a row iterator for a two-dimensional array. @@ -260,7 +260,7 @@ impl + Shape< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue<2, Item = Item> + Shape<2> + UnsafeRandomAccessMut<2, Item = Item>, @@ -279,7 +279,7 @@ impl< pub struct ColIterator< 'a, - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const NDIM: usize, > { @@ -290,7 +290,7 @@ pub struct ColIterator< pub struct ColIteratorMut< 'a, - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + UnsafeRandomAccessMut, @@ -301,7 +301,7 @@ pub struct ColIteratorMut< current_col: usize, } -impl<'a, Item: Scalar, ArrayImpl: UnsafeRandomAccessByValue<2, Item = Item> + Shape<2>> +impl<'a, Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue<2, Item = Item> + Shape<2>> std::iter::Iterator for ColIterator<'a, Item, ArrayImpl, 2> { type Item = Array, 2, 1>, 1>; @@ -317,7 +317,7 @@ impl<'a, Item: Scalar, ArrayImpl: UnsafeRandomAccessByValue<2, Item = Item> + Sh impl< 'a, - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue<2, Item = Item> + UnsafeRandomAccessMut<2, Item = Item> + Shape<2>, @@ -339,7 +339,7 @@ impl< } } -impl + Shape<2>> +impl + Shape<2>> Array { /// Return a column iterator for a two-dimensional array. @@ -353,7 +353,7 @@ impl + Shape< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue<2, Item = Item> + Shape<2> + UnsafeRandomAccessMut<2, Item = Item>, diff --git a/dense/src/array/mult_into.rs b/dense/src/array/mult_into.rs index 10a619a7..b6169be3 100644 --- a/dense/src/array/mult_into.rs +++ b/dense/src/array/mult_into.rs @@ -1,14 +1,14 @@ //! Implementation of array multiplication. +use crate::gemm::Gemm; use crate::traits::MultInto; use crate::traits::MultIntoResize; -use rlst_blis::interface::gemm::Gemm; -pub use rlst_blis::interface::types::TransMode; +pub use crate::types::TransMode; use super::{empty_axis::AxisPosition, *}; impl< - Item: Scalar + Gemm, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue<2, Item = Item> + UnsafeRandomAccessMut<2, Item = Item> + Shape<2> @@ -38,7 +38,7 @@ impl< } impl< - Item: Scalar + Gemm, + Item: RlstScalar + Gemm, ArrayImpl: UnsafeRandomAccessByValue<1, Item = Item> + UnsafeRandomAccessMut<1, Item = Item> + Shape<1> @@ -77,7 +77,7 @@ impl< } impl< - Item: Scalar + Gemm, + Item: RlstScalar + Gemm, ArrayImpl: UnsafeRandomAccessByValue<1, Item = Item> + UnsafeRandomAccessMut<1, Item = Item> + Shape<1> @@ -118,7 +118,7 @@ impl< /// MultIntoResize impl< - Item: Scalar + Gemm, + Item: RlstScalar + Gemm, ArrayImpl: UnsafeRandomAccessByValue<2, Item = Item> + UnsafeRandomAccessMut<2, Item = Item> + Shape<2> @@ -157,7 +157,7 @@ impl< } impl< - Item: Scalar + Gemm, + Item: RlstScalar + Gemm, ArrayImpl: UnsafeRandomAccessByValue<1, Item = Item> + UnsafeRandomAccessMut<1, Item = Item> + Shape<1> @@ -205,7 +205,7 @@ impl< } impl< - Item: Scalar + Gemm, + Item: RlstScalar + Gemm, ArrayImpl: UnsafeRandomAccessByValue<1, Item = Item> + UnsafeRandomAccessMut<1, Item = Item> + Shape<1> diff --git a/dense/src/array/operations.rs b/dense/src/array/operations.rs index c7b04565..9e62ab1a 100644 --- a/dense/src/array/operations.rs +++ b/dense/src/array/operations.rs @@ -1,13 +1,13 @@ //! Operations on arrays. -use num::{Float, Zero}; -use rlst_common::types::RlstResult; +use crate::types::RlstResult; +use num::Zero; use crate::{layout::convert_1d_nd_from_shape, traits::MatrixSvd}; use super::*; impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + UnsafeRandomAccessMut, @@ -171,7 +171,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const NDIM: usize, > Array @@ -191,7 +191,7 @@ impl< } } -impl + Shape<1>> +impl + Shape<1>> Array where Item::Real: num::Float, @@ -217,29 +217,29 @@ where } /// Compute the maximum (or inf) norm of a vector. - pub fn norm_inf(self) -> ::Real { + pub fn norm_inf(self) -> ::Real { self.iter() - .map(|elem| ::abs(elem)) - .reduce(<::Real as Float>::max) + .map(|elem| ::abs(elem)) + .reduce(<::Real as num::Float>::max) .unwrap() } /// Compute the 1-norm of a vector. - pub fn norm_1(self) -> ::Real { + pub fn norm_1(self) -> ::Real { self.iter() - .map(|elem| ::abs(elem)) - .fold(<::Real as Zero>::zero(), |acc, elem| { + .map(|elem| ::abs(elem)) + .fold(<::Real as Zero>::zero(), |acc, elem| { acc + elem }) } /// Compute the 2-norm of a vector. - pub fn norm_2(self) -> ::Real { - Scalar::sqrt( + pub fn norm_2(self) -> ::Real { + RlstScalar::sqrt( self.iter() - .map(|elem| ::abs(elem)) + .map(|elem| ::abs(elem)) .map(|elem| elem * elem) - .fold(<::Real as Zero>::zero(), |acc, elem| { + .fold(<::Real as Zero>::zero(), |acc, elem| { acc + elem }), ) @@ -281,7 +281,7 @@ where } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue<2, Item = Item> + Shape<2> + Stride<2> @@ -293,10 +293,10 @@ where /// Compute the 2-norm of a matrix. /// /// This method allocates temporary memory during execution. - pub fn norm_2_alloc(self) -> RlstResult<::Real> { + pub fn norm_2_alloc(self) -> RlstResult<::Real> { let k = *self.shape().iter().min().unwrap(); - let mut singular_values = vec![<::Real as Zero>::zero(); k]; + let mut singular_values = vec![<::Real as Zero>::zero(); k]; self.into_singular_values_alloc(singular_values.as_mut_slice())?; @@ -304,16 +304,16 @@ where } } -impl + Shape<2>> +impl + Shape<2>> Array { /// Compute the Frobenius-norm of a matrix. pub fn norm_fro(self) -> Item::Real { - Scalar::sqrt( + RlstScalar::sqrt( self.iter() - .map(|elem| ::abs(elem)) + .map(|elem| ::abs(elem)) .map(|elem| elem * elem) - .fold(<::Real as Zero>::zero(), |acc, elem| { + .fold(<::Real as Zero>::zero(), |acc, elem| { acc + elem }), ) @@ -323,7 +323,7 @@ impl + Shape< pub fn norm_inf(self) -> Item::Real { self.row_iter() .map(|row| row.norm_1()) - .reduce(<::Real as Float>::max) + .reduce(<::Real as num::Float>::max) .unwrap() } @@ -331,7 +331,7 @@ impl + Shape< pub fn norm_1(self) -> Item::Real { self.col_iter() .map(|row| row.norm_1()) - .reduce(<::Real as Float>::max) + .reduce(<::Real as num::Float>::max) .unwrap() } } diff --git a/dense/src/array/operators/addition.rs b/dense/src/array/operators/addition.rs index ce4a2361..73821678 100644 --- a/dense/src/array/operators/addition.rs +++ b/dense/src/array/operators/addition.rs @@ -3,7 +3,7 @@ use crate::array::*; pub struct ArrayAddition< - Item: Scalar, + Item: RlstScalar, ArrayImpl1: UnsafeRandomAccessByValue + Shape, ArrayImpl2: UnsafeRandomAccessByValue + Shape, const NDIM: usize, @@ -13,7 +13,7 @@ pub struct ArrayAddition< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl1: UnsafeRandomAccessByValue + Shape, ArrayImpl2: UnsafeRandomAccessByValue + Shape, const NDIM: usize, @@ -38,7 +38,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl1: UnsafeRandomAccessByValue + Shape, ArrayImpl2: UnsafeRandomAccessByValue + Shape, const NDIM: usize, @@ -53,7 +53,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl1: UnsafeRandomAccessByValue + Shape + ChunkedAccess, ArrayImpl2: UnsafeRandomAccessByValue + Shape + ChunkedAccess, const NDIM: usize, @@ -62,10 +62,7 @@ impl< { type Item = Item; #[inline] - fn get_chunk( - &self, - chunk_index: usize, - ) -> Option> { + fn get_chunk(&self, chunk_index: usize) -> Option> { if let (Some(mut chunk1), Some(chunk2)) = ( self.operator1.get_chunk(chunk_index), self.operator2.get_chunk(chunk_index), @@ -81,7 +78,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl1: UnsafeRandomAccessByValue + Shape, ArrayImpl2: UnsafeRandomAccessByValue + Shape, const NDIM: usize, @@ -93,7 +90,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl1: UnsafeRandomAccessByValue + Shape, ArrayImpl2: UnsafeRandomAccessByValue + Shape, const NDIM: usize, diff --git a/dense/src/array/operators/cmp_wise_division.rs b/dense/src/array/operators/cmp_wise_division.rs index 44ff82d9..9fc52124 100644 --- a/dense/src/array/operators/cmp_wise_division.rs +++ b/dense/src/array/operators/cmp_wise_division.rs @@ -3,7 +3,7 @@ use crate::array::*; pub struct CmpWiseDivision< - Item: Scalar, + Item: RlstScalar, ArrayImpl1: UnsafeRandomAccessByValue + Shape, ArrayImpl2: UnsafeRandomAccessByValue + Shape, const NDIM: usize, @@ -13,7 +13,7 @@ pub struct CmpWiseDivision< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl1: UnsafeRandomAccessByValue + Shape, ArrayImpl2: UnsafeRandomAccessByValue + Shape, const NDIM: usize, @@ -38,7 +38,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl1: UnsafeRandomAccessByValue + Shape, ArrayImpl2: UnsafeRandomAccessByValue + Shape, const NDIM: usize, @@ -53,7 +53,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl1: UnsafeRandomAccessByValue + Shape + ChunkedAccess, ArrayImpl2: UnsafeRandomAccessByValue + Shape + ChunkedAccess, const NDIM: usize, @@ -62,10 +62,7 @@ impl< { type Item = Item; #[inline] - fn get_chunk( - &self, - chunk_index: usize, - ) -> Option> { + fn get_chunk(&self, chunk_index: usize) -> Option> { if let (Some(mut chunk1), Some(chunk2)) = ( self.operator1.get_chunk(chunk_index), self.operator2.get_chunk(chunk_index), @@ -81,7 +78,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl1: UnsafeRandomAccessByValue + Shape, ArrayImpl2: UnsafeRandomAccessByValue + Shape, const NDIM: usize, @@ -93,7 +90,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl1: UnsafeRandomAccessByValue + Shape, ArrayImpl2: UnsafeRandomAccessByValue + Shape, const NDIM: usize, diff --git a/dense/src/array/operators/cmp_wise_product.rs b/dense/src/array/operators/cmp_wise_product.rs index edd69086..17cc2f88 100644 --- a/dense/src/array/operators/cmp_wise_product.rs +++ b/dense/src/array/operators/cmp_wise_product.rs @@ -3,7 +3,7 @@ use crate::array::*; pub struct CmpWiseProduct< - Item: Scalar, + Item: RlstScalar, ArrayImpl1: UnsafeRandomAccessByValue + Shape, ArrayImpl2: UnsafeRandomAccessByValue + Shape, const NDIM: usize, @@ -13,7 +13,7 @@ pub struct CmpWiseProduct< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl1: UnsafeRandomAccessByValue + Shape, ArrayImpl2: UnsafeRandomAccessByValue + Shape, const NDIM: usize, @@ -38,7 +38,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl1: UnsafeRandomAccessByValue + Shape, ArrayImpl2: UnsafeRandomAccessByValue + Shape, const NDIM: usize, @@ -53,7 +53,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl1: UnsafeRandomAccessByValue + Shape + ChunkedAccess, ArrayImpl2: UnsafeRandomAccessByValue + Shape + ChunkedAccess, const NDIM: usize, @@ -62,10 +62,7 @@ impl< { type Item = Item; #[inline] - fn get_chunk( - &self, - chunk_index: usize, - ) -> Option> { + fn get_chunk(&self, chunk_index: usize) -> Option> { if let (Some(mut chunk1), Some(chunk2)) = ( self.operator1.get_chunk(chunk_index), self.operator2.get_chunk(chunk_index), @@ -81,7 +78,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl1: UnsafeRandomAccessByValue + Shape, ArrayImpl2: UnsafeRandomAccessByValue + Shape, const NDIM: usize, @@ -93,7 +90,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl1: UnsafeRandomAccessByValue + Shape, ArrayImpl2: UnsafeRandomAccessByValue + Shape, const NDIM: usize, diff --git a/dense/src/array/operators/conjugate.rs b/dense/src/array/operators/conjugate.rs index a95e9b6c..36d13b0b 100644 --- a/dense/src/array/operators/conjugate.rs +++ b/dense/src/array/operators/conjugate.rs @@ -3,7 +3,7 @@ use crate::array::*; pub struct ArrayConjugate< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const NDIM: usize, > { @@ -11,7 +11,7 @@ pub struct ArrayConjugate< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const NDIM: usize, > ArrayConjugate @@ -22,7 +22,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const NDIM: usize, > UnsafeRandomAccessByValue for ArrayConjugate @@ -35,7 +35,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const NDIM: usize, > Shape for ArrayConjugate @@ -46,7 +46,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const NDIM: usize, > Array @@ -57,7 +57,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + ChunkedAccess, const NDIM: usize, const N: usize, diff --git a/dense/src/array/operators/other.rs b/dense/src/array/operators/other.rs index 696dd30e..d62ccee8 100644 --- a/dense/src/array/operators/other.rs +++ b/dense/src/array/operators/other.rs @@ -7,7 +7,7 @@ use crate::array::*; use num::One; impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl1: UnsafeRandomAccessByValue + Shape, ArrayImpl2: UnsafeRandomAccessByValue + Shape, const NDIM: usize, diff --git a/dense/src/array/operators/scalar_mult.rs b/dense/src/array/operators/scalar_mult.rs index c1a5e7ae..1d365979 100644 --- a/dense/src/array/operators/scalar_mult.rs +++ b/dense/src/array/operators/scalar_mult.rs @@ -1,10 +1,10 @@ //! Container representing multiplication with a scalar use crate::array::*; -use rlst_common::types::*; +use crate::types::*; pub struct ArrayScalarMult< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const NDIM: usize, > { @@ -13,7 +13,7 @@ pub struct ArrayScalarMult< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const NDIM: usize, > ArrayScalarMult @@ -24,7 +24,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const NDIM: usize, > UnsafeRandomAccessByValue for ArrayScalarMult @@ -37,7 +37,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const NDIM: usize, > Shape for ArrayScalarMult @@ -48,7 +48,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + ChunkedAccess, const NDIM: usize, const N: usize, @@ -69,7 +69,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const NDIM: usize, > Array diff --git a/dense/src/array/operators/to_complex.rs b/dense/src/array/operators/to_complex.rs index 4e594094..e080479d 100644 --- a/dense/src/array/operators/to_complex.rs +++ b/dense/src/array/operators/to_complex.rs @@ -1,26 +1,26 @@ //! Container representing multiplication with a scalar -use rlst_common::types::{c32, c64}; +use crate::types::{c32, c64}; use std::marker::PhantomData; use crate::array::*; pub struct ArrayToComplex< - Item: Scalar, - ArrayImpl: UnsafeRandomAccessByValue::Real> + Shape, + Item: RlstScalar, + ArrayImpl: UnsafeRandomAccessByValue::Real> + Shape, const NDIM: usize, > { - operator: Array<::Real, ArrayImpl, NDIM>, + operator: Array<::Real, ArrayImpl, NDIM>, _marker: PhantomData, } impl< - Item: Scalar, - ArrayImpl: UnsafeRandomAccessByValue::Real> + Shape, + Item: RlstScalar, + ArrayImpl: UnsafeRandomAccessByValue::Real> + Shape, const NDIM: usize, > ArrayToComplex { - pub fn new(operator: Array<::Real, ArrayImpl, NDIM>) -> Self { + pub fn new(operator: Array<::Real, ArrayImpl, NDIM>) -> Self { Self { operator, _marker: PhantomData, @@ -29,21 +29,21 @@ impl< } impl< - Item: Scalar, - ArrayImpl: UnsafeRandomAccessByValue::Real> + Shape, + Item: RlstScalar, + ArrayImpl: UnsafeRandomAccessByValue::Real> + Shape, const NDIM: usize, > UnsafeRandomAccessByValue for ArrayToComplex { type Item = Item; #[inline] unsafe fn get_value_unchecked(&self, multi_index: [usize; NDIM]) -> Self::Item { - ::from_real(self.operator.get_value_unchecked(multi_index)) + ::from_real(self.operator.get_value_unchecked(multi_index)) } } impl< - Item: Scalar, - ArrayImpl: UnsafeRandomAccessByValue::Real> + Shape, + Item: RlstScalar, + ArrayImpl: UnsafeRandomAccessByValue::Real> + Shape, const NDIM: usize, > Shape for ArrayToComplex { @@ -53,10 +53,10 @@ impl< } impl< - Item: Scalar, - ArrayImpl: UnsafeRandomAccessByValue::Real> + Item: RlstScalar, + ArrayImpl: UnsafeRandomAccessByValue::Real> + Shape - + ChunkedAccess::Real>, + + ChunkedAccess::Real>, const NDIM: usize, const N: usize, > ChunkedAccess for ArrayToComplex @@ -68,7 +68,7 @@ impl< let mut data = [::zero(); N]; for (d, &c) in data.iter_mut().zip(chunk.data.iter()) { - *d = ::from_real(c); + *d = ::from_real(c); } Some(DataChunk:: { data, diff --git a/dense/src/array/operators/transpose.rs b/dense/src/array/operators/transpose.rs index a6ee7707..0cb1ef3c 100644 --- a/dense/src/array/operators/transpose.rs +++ b/dense/src/array/operators/transpose.rs @@ -3,7 +3,7 @@ use crate::{array::*, layout::convert_1d_nd_from_shape}; pub struct ArrayTranspose< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const NDIM: usize, > { @@ -13,7 +13,7 @@ pub struct ArrayTranspose< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const NDIM: usize, > ArrayTranspose @@ -34,7 +34,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const NDIM: usize, > UnsafeRandomAccessByValue for ArrayTranspose @@ -48,7 +48,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + UnsafeRandomAccessByRef, @@ -64,7 +64,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + UnsafeRandomAccessMut, @@ -80,7 +80,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const NDIM: usize, > Shape for ArrayTranspose @@ -91,7 +91,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const NDIM: usize, > Array @@ -115,7 +115,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + ChunkedAccess, const NDIM: usize, const N: usize, diff --git a/dense/src/array/random.rs b/dense/src/array/random.rs index fa3e19e3..2ae8f4c4 100644 --- a/dense/src/array/random.rs +++ b/dense/src/array/random.rs @@ -2,24 +2,24 @@ use crate::tools::*; use crate::traits::*; +use crate::types::RlstScalar; use rand::prelude::*; use rand_chacha::ChaCha8Rng; use rand_distr::Standard; use rand_distr::StandardNormal; -use rlst_common::types::Scalar; use super::Array; impl< - Item: Scalar + RandScalar, + Item: RlstScalar + RandScalar, ArrayImpl: UnsafeRandomAccessByValue + UnsafeRandomAccessMut + Shape, const NDIM: usize, > Array where - StandardNormal: Distribution<::Real>, - Standard: Distribution<::Real>, + StandardNormal: Distribution<::Real>, + Standard: Distribution<::Real>, { /// Fill an array with normally distributed random numbers. pub fn fill_from_standard_normal(&mut self, rng: &mut R) { diff --git a/dense/src/array/slice.rs b/dense/src/array/slice.rs index a870a75d..8429cc96 100644 --- a/dense/src/array/slice.rs +++ b/dense/src/array/slice.rs @@ -9,7 +9,7 @@ use super::*; /// Generic structure to store Array slices. pub struct ArraySlice< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const ADIM: usize, const NDIM: usize, @@ -26,7 +26,7 @@ pub struct ArraySlice< // Implementation of ArraySlice impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const ADIM: usize, const NDIM: usize, @@ -58,7 +58,7 @@ where } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const ADIM: usize, const NDIM: usize, @@ -77,7 +77,7 @@ where } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + UnsafeRandomAccessByRef, @@ -98,7 +98,7 @@ where } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const ADIM: usize, const NDIM: usize, @@ -120,7 +120,7 @@ where } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + RawAccess @@ -143,7 +143,7 @@ where } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + Stride, const ADIM: usize, const NDIM: usize, @@ -165,7 +165,7 @@ where } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const ADIM: usize, > Array @@ -195,7 +195,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + Stride, const ADIM: usize, const NDIM: usize, @@ -208,10 +208,7 @@ where type Item = Item; #[inline] - fn get_chunk( - &self, - chunk_index: usize, - ) -> Option> { + fn get_chunk(&self, chunk_index: usize) -> Option> { let nelements = self.shape().iter().product(); if let Some(mut chunk) = empty_chunk(chunk_index, nelements) { for count in 0..chunk.valid_entries { @@ -230,7 +227,7 @@ where } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + UnsafeRandomAccessByRef @@ -252,7 +249,7 @@ where } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + RawAccessMut diff --git a/dense/src/array/views.rs b/dense/src/array/views.rs index 1afa943f..f6d51216 100644 --- a/dense/src/array/views.rs +++ b/dense/src/array/views.rs @@ -7,12 +7,12 @@ use crate::layout::{check_multi_index_in_bounds, convert_1d_nd_from_shape}; use super::Array; use crate::traits::*; -use rlst_common::types::*; +use crate::types::*; /// Basic structure for a `View` pub struct ArrayView< 'a, - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const NDIM: usize, > { @@ -21,7 +21,7 @@ pub struct ArrayView< pub struct ArrayViewMut< 'a, - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + UnsafeRandomAccessMut, @@ -32,7 +32,7 @@ pub struct ArrayViewMut< impl< 'a, - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const NDIM: usize, > ArrayView<'a, Item, ArrayImpl, NDIM> @@ -44,7 +44,7 @@ impl< impl< 'a, - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + UnsafeRandomAccessMut, @@ -61,7 +61,7 @@ impl< impl< 'a, - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const NDIM: usize, > Shape for ArrayView<'a, Item, ArrayImpl, NDIM> @@ -73,7 +73,7 @@ impl< impl< 'a, - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + Stride, const NDIM: usize, > Stride for ArrayView<'a, Item, ArrayImpl, NDIM> @@ -85,7 +85,7 @@ impl< impl< 'a, - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + RawAccess @@ -102,7 +102,7 @@ impl< impl< 'a, - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const NDIM: usize, > UnsafeRandomAccessByValue for ArrayView<'a, Item, ArrayImpl, NDIM> @@ -116,7 +116,7 @@ impl< impl< 'a, - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + UnsafeRandomAccessByRef, @@ -132,17 +132,14 @@ impl< impl< 'a, - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + ChunkedAccess, const NDIM: usize, const N: usize, > ChunkedAccess for ArrayView<'a, Item, ArrayImpl, NDIM> { type Item = Item; - fn get_chunk( - &self, - chunk_index: usize, - ) -> Option> { + fn get_chunk(&self, chunk_index: usize) -> Option> { self.arr.get_chunk(chunk_index) } } @@ -152,7 +149,7 @@ impl< impl< 'a, - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + UnsafeRandomAccessMut, @@ -166,7 +163,7 @@ impl< impl< 'a, - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + Stride @@ -181,7 +178,7 @@ impl< impl< 'a, - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + RawAccess @@ -199,7 +196,7 @@ impl< impl< 'a, - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + RawAccess @@ -216,7 +213,7 @@ impl< impl< 'a, - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + UnsafeRandomAccessMut, @@ -232,7 +229,7 @@ impl< impl< 'a, - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + UnsafeRandomAccessByRef @@ -249,7 +246,7 @@ impl< impl< 'a, - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + ChunkedAccess @@ -259,17 +256,14 @@ impl< > ChunkedAccess for ArrayViewMut<'a, Item, ArrayImpl, NDIM> { type Item = Item; - fn get_chunk( - &self, - chunk_index: usize, - ) -> Option> { + fn get_chunk(&self, chunk_index: usize) -> Option> { self.arr.get_chunk(chunk_index) } } impl< 'a, - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + UnsafeRandomAccessByRef @@ -286,7 +280,7 @@ impl< impl< 'a, - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + UnsafeRandomAccessByRef @@ -303,7 +297,7 @@ impl< //////////////////////////////////////////// pub struct ArraySubView< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const NDIM: usize, > { @@ -313,7 +307,7 @@ pub struct ArraySubView< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const NDIM: usize, > ArraySubView @@ -340,7 +334,7 @@ impl< // Basic traits for ArraySubView impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const NDIM: usize, > Shape for ArraySubView @@ -351,7 +345,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + Stride, const NDIM: usize, > Stride for ArraySubView @@ -362,7 +356,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + RawAccess @@ -381,7 +375,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const NDIM: usize, > UnsafeRandomAccessByValue for ArraySubView @@ -396,7 +390,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + UnsafeRandomAccessByRef, @@ -413,7 +407,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + ChunkedAccess, const NDIM: usize, const N: usize, @@ -421,10 +415,7 @@ impl< { type Item = Item; #[inline] - fn get_chunk( - &self, - chunk_index: usize, - ) -> Option> { + fn get_chunk(&self, chunk_index: usize) -> Option> { if self.offset == [0; NDIM] && self.shape() == self.arr.shape() { // If the view is on the full array we can just pass on the chunk request self.arr.get_chunk(chunk_index) @@ -451,7 +442,7 @@ impl< // Basic traits for ArrayViewMut impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + UnsafeRandomAccessMut @@ -469,7 +460,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + UnsafeRandomAccessMut, @@ -486,7 +477,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const NDIM: usize, > Array @@ -502,7 +493,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape, const NDIM: usize, > Array @@ -514,7 +505,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue + Shape + UnsafeRandomAccessMut, diff --git a/dense/src/base_array.rs b/dense/src/base_array.rs index 1a393c29..50bddc6c 100644 --- a/dense/src/base_array.rs +++ b/dense/src/base_array.rs @@ -9,17 +9,17 @@ use crate::layout::{ check_multi_index_in_bounds, convert_1d_nd_from_shape, convert_nd_raw, stride_from_shape, }; use crate::traits::*; -use rlst_common::types::Scalar; +use crate::types::RlstScalar; /// Definition of a [BaseArray]. The `data` stores the actual array data, `shape` stores /// the shape of the array, and `stride` contains the `stride` of the underlying data. -pub struct BaseArray, const NDIM: usize> { +pub struct BaseArray, const NDIM: usize> { data: Data, shape: [usize; NDIM], stride: [usize; NDIM], } -impl, const NDIM: usize> +impl, const NDIM: usize> BaseArray { pub fn new(data: Data, shape: [usize; NDIM]) -> Self { @@ -58,7 +58,7 @@ impl, const NDIM: usize> } } -impl, const NDIM: usize> Shape +impl, const NDIM: usize> Shape for BaseArray { fn shape(&self) -> [usize; NDIM] { @@ -66,7 +66,7 @@ impl, const NDIM: usize> Shape, const NDIM: usize> +impl, const NDIM: usize> UnsafeRandomAccessByRef for BaseArray { type Item = Item; @@ -79,7 +79,7 @@ impl, const NDIM: usize> } } -impl, const NDIM: usize> +impl, const NDIM: usize> UnsafeRandomAccessByValue for BaseArray { type Item = Item; @@ -92,7 +92,7 @@ impl, const NDIM: usize> } } -impl, const NDIM: usize> +impl, const NDIM: usize> UnsafeRandomAccessMut for BaseArray { type Item = Item; @@ -105,16 +105,13 @@ impl, const NDIM: usize> } } -impl, const N: usize, const NDIM: usize> +impl, const N: usize, const NDIM: usize> ChunkedAccess for BaseArray { type Item = Item; #[inline] - fn get_chunk( - &self, - chunk_index: usize, - ) -> Option> { + fn get_chunk(&self, chunk_index: usize) -> Option> { let nelements = self.shape().iter().product(); if let Some(mut chunk) = empty_chunk(chunk_index, nelements) { for count in 0..chunk.valid_entries { @@ -132,7 +129,7 @@ impl, const N: usize, const ND } } -impl, const NDIM: usize> RawAccess +impl, const NDIM: usize> RawAccess for BaseArray { type Item = Item; @@ -142,7 +139,7 @@ impl, const NDIM: usize> RawAcces } } -impl, const NDIM: usize> RawAccessMut +impl, const NDIM: usize> RawAccessMut for BaseArray { fn data_mut(&mut self) -> &mut [Self::Item] { @@ -150,7 +147,7 @@ impl, const NDIM: usize> RawAc } } -impl, const NDIM: usize> Stride +impl, const NDIM: usize> Stride for BaseArray { fn stride(&self) -> [usize; NDIM] { @@ -158,7 +155,7 @@ impl, const NDIM: usize> Stride, const NDIM: usize> +impl, const NDIM: usize> ResizeInPlace for BaseArray { fn resize_in_place(&mut self, shape: [usize; NDIM]) { diff --git a/dense/src/data_container.rs b/dense/src/data_container.rs index 979c525a..3a25c487 100644 --- a/dense/src/data_container.rs +++ b/dense/src/data_container.rs @@ -10,12 +10,12 @@ //! - [SliceContainerMut] - Like [SliceContainer] but provides also mutable access. //! +use crate::types::RlstScalar; use num::Zero; -use rlst_common::types::Scalar; /// Defines the basic behaviour of a data container. pub trait DataContainer { - type Item: Scalar; + type Item: RlstScalar; /// Access the container unchecked. /// @@ -57,7 +57,7 @@ pub trait ResizeableDataContainerMut: DataContainerMut { /// A container that uses dynamic vectors. #[derive(Clone)] -pub struct VectorContainer { +pub struct VectorContainer { data: Vec, } @@ -66,21 +66,21 @@ pub struct VectorContainer { /// The size of this container needs to be known at compile time. /// It is useful for data structures that should be stack allocated. #[derive(Clone)] -pub struct ArrayContainer { +pub struct ArrayContainer { data: [Item; N], } /// A container that takes a reference to a slice. -pub struct SliceContainer<'a, Item: Scalar> { +pub struct SliceContainer<'a, Item: RlstScalar> { data: &'a [Item], } /// A container that takes a mutable reference to a slice. -pub struct SliceContainerMut<'a, Item: Scalar> { +pub struct SliceContainerMut<'a, Item: RlstScalar> { data: &'a mut [Item], } -impl VectorContainer { +impl VectorContainer { /// New vector container by specifying the number of elements. /// /// The container is initialized with zeros by default. @@ -91,7 +91,7 @@ impl VectorContainer { } } -impl ArrayContainer { +impl ArrayContainer { pub fn new() -> ArrayContainer { ArrayContainer:: { data: [num::cast::(0.0).unwrap(); N], @@ -99,27 +99,27 @@ impl ArrayContainer { } } -impl Default for ArrayContainer { +impl Default for ArrayContainer { fn default() -> Self { Self::new() } } -impl<'a, Item: Scalar> SliceContainer<'a, Item> { +impl<'a, Item: RlstScalar> SliceContainer<'a, Item> { /// New slice container from a reference. pub fn new(slice: &'a [Item]) -> SliceContainer { SliceContainer:: { data: slice } } } -impl<'a, Item: Scalar> SliceContainerMut<'a, Item> { +impl<'a, Item: RlstScalar> SliceContainerMut<'a, Item> { /// New mutable slice container from mutable reference. pub fn new(slice: &'a mut [Item]) -> SliceContainerMut { SliceContainerMut:: { data: slice } } } -impl DataContainer for VectorContainer { +impl DataContainer for VectorContainer { type Item = Item; unsafe fn get_unchecked_value(&self, index: usize) -> Self::Item { @@ -141,7 +141,7 @@ impl DataContainer for VectorContainer { } } -impl DataContainerMut for VectorContainer { +impl DataContainerMut for VectorContainer { unsafe fn get_unchecked_mut(&mut self, index: usize) -> &mut Self::Item { debug_assert!(index < self.number_of_elements()); self.data.get_unchecked_mut(index) @@ -152,13 +152,13 @@ impl DataContainerMut for VectorContainer { } } -impl ResizeableDataContainerMut for VectorContainer { +impl ResizeableDataContainerMut for VectorContainer { fn resize(&mut self, new_len: usize) { self.data.resize(new_len, ::zero()); } } -impl DataContainer for ArrayContainer { +impl DataContainer for ArrayContainer { type Item = Item; unsafe fn get_unchecked_value(&self, index: usize) -> Self::Item { @@ -180,7 +180,7 @@ impl DataContainer for ArrayContainer { } } -impl DataContainerMut for ArrayContainer { +impl DataContainerMut for ArrayContainer { unsafe fn get_unchecked_mut(&mut self, index: usize) -> &mut Self::Item { debug_assert!(index < self.number_of_elements()); self.data.get_unchecked_mut(index) @@ -191,7 +191,7 @@ impl DataContainerMut for ArrayContainer } } -impl<'a, Item: Scalar> DataContainer for SliceContainer<'a, Item> { +impl<'a, Item: RlstScalar> DataContainer for SliceContainer<'a, Item> { type Item = Item; unsafe fn get_unchecked_value(&self, index: usize) -> Self::Item { @@ -213,7 +213,7 @@ impl<'a, Item: Scalar> DataContainer for SliceContainer<'a, Item> { } } -impl<'a, Item: Scalar> DataContainer for SliceContainerMut<'a, Item> { +impl<'a, Item: RlstScalar> DataContainer for SliceContainerMut<'a, Item> { type Item = Item; unsafe fn get_unchecked_value(&self, index: usize) -> Self::Item { @@ -235,7 +235,7 @@ impl<'a, Item: Scalar> DataContainer for SliceContainerMut<'a, Item> { } } -impl<'a, Item: Scalar> DataContainerMut for SliceContainerMut<'a, Item> { +impl<'a, Item: RlstScalar> DataContainerMut for SliceContainerMut<'a, Item> { unsafe fn get_unchecked_mut(&mut self, index: usize) -> &mut Self::Item { debug_assert!(index < self.number_of_elements()); self.data.get_unchecked_mut(index) diff --git a/dense/src/gemm.rs b/dense/src/gemm.rs new file mode 100644 index 00000000..418ac28b --- /dev/null +++ b/dense/src/gemm.rs @@ -0,0 +1,26 @@ +//! Gemm trait for matrix multiplication +use crate::types::TransMode; + +pub trait Gemm: Sized { + #[allow(clippy::too_many_arguments)] + fn gemm( + transa: TransMode, + transb: TransMode, + m: usize, + n: usize, + k: usize, + alpha: Self, + a: &[Self], + rsa: usize, + csa: usize, + b: &[Self], + rsb: usize, + csb: usize, + beta: Self, + c: &mut [Self], + rsc: usize, + csc: usize, + ); +} + +pub mod blis; diff --git a/dense/src/gemm/blis.rs b/dense/src/gemm/blis.rs new file mode 100644 index 00000000..9a9aaa6e --- /dev/null +++ b/dense/src/gemm/blis.rs @@ -0,0 +1,102 @@ +//! Blis implementation of Gemm + +extern crate rlst_blis_src; +use crate::gemm::Gemm; +use crate::types::{c32, c64, TransMode}; +use rlst_blis::raw; + +/// Compute expected size of a data slice from stride and shape. +fn get_expected_data_size(stride: [usize; 2], shape: [usize; 2]) -> usize { + if shape[0] == 0 || shape[1] == 0 { + return 0; + } + + 1 + (shape[0] - 1) * stride[0] + (shape[1] - 1) * stride[1] +} + +/// Panic if expected data size is not identical to actual data size. +fn assert_data_size(nelems: usize, stride: [usize; 2], shape: [usize; 2]) { + let expected = get_expected_data_size(stride, shape); + + assert_eq!( + expected, nelems, + "Wrong size for data slice. Actual size {}. Expected size {}.", + nelems, expected + ); +} + +fn convert_trans_mode(trans: TransMode) -> u32 { + match trans { + TransMode::ConjNoTrans => raw::trans_t_BLIS_CONJ_NO_TRANSPOSE, + TransMode::NoTrans => raw::trans_t_BLIS_NO_TRANSPOSE, + TransMode::Trans => raw::trans_t_BLIS_TRANSPOSE, + TransMode::ConjTrans => raw::trans_t_BLIS_CONJ_TRANSPOSE, + } +} + +macro_rules! impl_gemm { + ($scalar:ty, $blis_gemm:ident, $blis_type:ty) => { + impl Gemm for $scalar { + fn gemm( + transa: TransMode, + transb: TransMode, + m: usize, + n: usize, + k: usize, + alpha: Self, + a: &[Self], + rsa: usize, + csa: usize, + b: &[Self], + rsb: usize, + csb: usize, + beta: Self, + c: &mut [Self], + rsc: usize, + csc: usize, + ) { + match transa { + TransMode::NoTrans => assert_data_size(a.len(), [rsa, csa], [m, k]), + TransMode::ConjNoTrans => assert_data_size(a.len(), [rsa, csa], [m, k]), + TransMode::Trans => assert_data_size(a.len(), [rsa, csa], [k, m]), + TransMode::ConjTrans => assert_data_size(a.len(), [rsa, csa], [k, m]), + } + + match transb { + TransMode::NoTrans => assert_data_size(b.len(), [rsb, csb], [k, n]), + TransMode::ConjNoTrans => assert_data_size(b.len(), [rsb, csb], [k, n]), + TransMode::Trans => assert_data_size(b.len(), [rsb, csb], [n, k]), + TransMode::ConjTrans => assert_data_size(b.len(), [rsb, csb], [n, k]), + } + + assert_data_size(c.len(), [rsc, csc], [m, n]); + + unsafe { + raw::$blis_gemm( + convert_trans_mode(transa), + convert_trans_mode(transb), + m as i64, + n as i64, + k as i64, + &alpha as *const _ as *const $blis_type, + a.as_ptr() as *const _ as *const $blis_type, + rsa as i64, + csa as i64, + b.as_ptr() as *const _ as *const $blis_type, + rsb as i64, + csb as i64, + &beta as *const _ as *const $blis_type, + c.as_mut_ptr() as *mut _ as *mut $blis_type, + rsc as i64, + csc as i64, + ) + }; + } + } + }; +} + +impl_gemm!(f32, bli_sgemm, f32); +impl_gemm!(f64, bli_dgemm, f64); +impl_gemm!(c32, bli_cgemm, raw::scomplex); +impl_gemm!(c64, bli_zgemm, raw::dcomplex); diff --git a/dense/src/lib.rs b/dense/src/lib.rs index 84ee7444..e7c67e8e 100644 --- a/dense/src/lib.rs +++ b/dense/src/lib.rs @@ -13,6 +13,8 @@ pub mod tools; pub mod traits; pub mod array; +pub mod gemm; pub mod layout; pub mod macros; pub mod matrix_multiply; +pub mod types; diff --git a/dense/src/linalg/inverse.rs b/dense/src/linalg/inverse.rs index 0efafd77..32d6d32b 100644 --- a/dense/src/linalg/inverse.rs +++ b/dense/src/linalg/inverse.rs @@ -3,9 +3,9 @@ //! use crate::array::Array; use crate::traits::*; +use crate::types::{c32, c64, RlstError, RlstResult, RlstScalar}; use lapack::{cgetrf, cgetri, dgetrf, dgetri, sgetrf, sgetri, zgetrf, zgetri}; use num::traits::Zero; -use rlst_common::types::{c32, c64, RlstError, RlstResult, Scalar}; use super::assert_lapack_stride; diff --git a/dense/src/linalg/lu.rs b/dense/src/linalg/lu.rs index b951d945..591838ff 100644 --- a/dense/src/linalg/lu.rs +++ b/dense/src/linalg/lu.rs @@ -2,9 +2,9 @@ use super::assert_lapack_stride; use crate::array::Array; use crate::traits::*; +use crate::types::*; use lapack::{cgetrf, cgetrs, dgetrf, dgetrs, sgetrf, sgetrs, zgetrf, zgetrs}; use num::One; -use rlst_common::types::*; /// Compute the LU decomposition of a matrix. /// @@ -13,7 +13,7 @@ use rlst_common::types::*; /// `L` is a `(m, k)` unit lower triangular matrix, and `U` is /// an `(k, n)` upper triangular matrix. pub trait MatrixLuDecomposition: Sized { - type Item: Scalar; + type Item: RlstScalar; type ArrayImpl: UnsafeRandomAccessByValue<2, Item = Self::Item> + Stride<2> + RawAccessMut @@ -155,7 +155,7 @@ pub enum LuTrans { /// Container for the LU Decomposition of a matrix. pub struct LuDecomposition< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue<2, Item = Item> + Stride<2> + Shape<2> + RawAccessMut, > { arr: Array, diff --git a/dense/src/linalg/pseudo_inverse.rs b/dense/src/linalg/pseudo_inverse.rs index de34ef7d..5ba1f687 100644 --- a/dense/src/linalg/pseudo_inverse.rs +++ b/dense/src/linalg/pseudo_inverse.rs @@ -3,14 +3,14 @@ use crate::array::Array; use crate::rlst_dynamic_array2; use crate::traits::*; +use crate::types::{c32, c64, RlstResult, RlstScalar}; use itertools::Itertools; use num::traits::{One, Zero}; -use rlst_common::types::{c32, c64, RlstResult, Scalar}; use crate::linalg::svd::*; pub trait MatrixPseudoInverse { - type Item: Scalar; + type Item: RlstScalar; /// Compute the pseudo inverse into the array `pinv`. /// @@ -29,7 +29,7 @@ pub trait MatrixPseudoInverse { >( self, pinv: Array, - tol: ::Real, + tol: ::Real, ) -> RlstResult<()>; /// Compute the pseudo inverse into the array `pinv`. @@ -50,7 +50,7 @@ pub trait MatrixPseudoInverse { >( self, pinv: Array, - tol: ::Real, + tol: ::Real, ) -> RlstResult<()>; } @@ -75,7 +75,7 @@ macro_rules! impl_pinv { >( self, mut pinv: Array<$scalar, ArrayImplPInv, 2>, - tol: <$scalar as Scalar>::Real, + tol: <$scalar as RlstScalar>::Real, ) -> RlstResult<()> { let shape = pinv.shape(); @@ -94,12 +94,12 @@ macro_rules! impl_pinv { >( self, mut pinv: Array<$scalar, ArrayImplPInv, 2>, - tol: <$scalar as Scalar>::Real, + tol: <$scalar as RlstScalar>::Real, ) -> RlstResult<()> { let shape = self.shape(); let k = std::cmp::min(self.shape()[0], self.shape()[1]); let mode = crate::linalg::svd::SvdMode::Reduced; - let mut singvals = vec![<<$scalar as Scalar>::Real as Zero>::zero(); k]; + let mut singvals = vec![<<$scalar as RlstScalar>::Real as Zero>::zero(); k]; let mut u = rlst_dynamic_array2!($scalar, [self.shape()[0], k]); let mut vt = rlst_dynamic_array2!($scalar, [k, self.shape()[1]]); @@ -136,7 +136,7 @@ macro_rules! impl_pinv { for (col_index, &singval) in singvals.iter().take(index).enumerate() { v.view_mut().slice(1, col_index).scale_in_place( - (<<$scalar as Scalar>::Real as One>::one() / singval).into(), + (<<$scalar as RlstScalar>::Real as One>::one() / singval).into(), ); } diff --git a/dense/src/linalg/qr.rs b/dense/src/linalg/qr.rs index 6a80fdfb..3e39df43 100644 --- a/dense/src/linalg/qr.rs +++ b/dense/src/linalg/qr.rs @@ -6,8 +6,8 @@ use crate::traits::*; use itertools::Itertools; use lapack::{cgeqp3, cunmqr, dgeqp3, dormqr, sgeqp3, sormqr, zgeqp3, zunmqr}; +use crate::types::*; use num::Zero; -use rlst_common::types::*; #[derive(Clone, Copy)] #[repr(u8)] @@ -25,7 +25,7 @@ pub enum ApplyQSide { /// `k=min(m, n)`. `P` is of dimension `(n, n)`. The user can choose in the /// method [MatrixQrDecomposition::get_q_alloc] how many columns of `Q` to return. pub trait MatrixQrDecomposition: Sized { - type Item: Scalar; + type Item: RlstScalar; type ArrayImpl: UnsafeRandomAccessByValue<2, Item = Self::Item> + Stride<2> + RawAccessMut @@ -110,7 +110,7 @@ pub enum ApplyQTrans { } pub struct QrDecomposition< - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue<2, Item = Item> + Stride<2> + Shape<2> + RawAccessMut, > { arr: Array, @@ -412,7 +412,8 @@ macro_rules! implement_qr_complex { let mut jpvt = vec![0 as i32; n as usize]; let mut tau = vec![<$scalar as Zero>::zero(); k]; - let mut rwork = vec![<<$scalar as Scalar>::Real as Zero>::zero(); 2 * n as usize]; + let mut rwork = + vec![<<$scalar as RlstScalar>::Real as Zero>::zero(); 2 * n as usize]; let mut work_query = [<$scalar as Zero>::zero()]; let lwork = -1; @@ -660,7 +661,6 @@ implement_qr_complex!(c32, cgeqp3, cunmqr); mod test { use crate::{assert_array_abs_diff_eq, assert_array_relative_eq, traits::*}; - use rlst_common::types::*; use crate::array::empty_array; use crate::rlst_dynamic_array2; @@ -703,8 +703,8 @@ mod test { assert_array_relative_eq!(actual, mat2, $tol); let qtq = empty_array::<$scalar, 2>().mult_into_resize( - rlst_blis::interface::types::TransMode::ConjTrans, - rlst_blis::interface::types::TransMode::NoTrans, + crate::types::TransMode::ConjTrans, + crate::types::TransMode::NoTrans, 1.0.into(), q_mat.view(), q_mat.view(), @@ -745,8 +745,8 @@ mod test { assert_array_relative_eq!(actual, mat2, $tol); let qtq = empty_array::<$scalar, 2>().mult_into_resize( - rlst_blis::interface::types::TransMode::ConjTrans, - rlst_blis::interface::types::TransMode::NoTrans, + crate::types::TransMode::ConjTrans, + crate::types::TransMode::NoTrans, 1.0.into(), q_mat.view(), q_mat.view(), diff --git a/dense/src/linalg/svd.rs b/dense/src/linalg/svd.rs index 825b5e5a..c0f5d597 100644 --- a/dense/src/linalg/svd.rs +++ b/dense/src/linalg/svd.rs @@ -1,14 +1,14 @@ //! Singular Value Decomposition. use crate::array::Array; use crate::traits::*; +use crate::types::{c32, c64, RlstError, RlstResult, RlstScalar}; use lapack::{cgesvd, dgesvd, sgesvd, zgesvd}; use num::traits::Zero; -use rlst_common::types::{c32, c64, RlstError, RlstResult, Scalar}; use super::assert_lapack_stride; pub trait MatrixSvd { - type Item: Scalar; + type Item: RlstScalar; /// Compute the singular values of the matrix. /// @@ -18,7 +18,7 @@ pub trait MatrixSvd { /// This method allocates temporary memory during execution. fn into_singular_values_alloc( self, - singular_values: &mut [::Real], + singular_values: &mut [::Real], ) -> RlstResult<()>; /// Compute the singular value decomposition. @@ -52,7 +52,7 @@ pub trait MatrixSvd { self, u: Array, vt: Array, - singular_values: &mut [::Real], + singular_values: &mut [::Real], mode: SvdMode, ) -> RlstResult<()>; } @@ -75,7 +75,7 @@ macro_rules! impl_svd_real { fn into_singular_values_alloc( mut self, - singular_values: &mut [<$scalar as Scalar>::Real], + singular_values: &mut [<$scalar as RlstScalar>::Real], ) -> RlstResult<()> { let lwork: i32 = -1; let mut work = [<$scalar as Zero>::zero(); 1]; @@ -158,7 +158,7 @@ macro_rules! impl_svd_real { mut self, mut u: Array<$scalar, ArrayImplU, 2>, mut vt: Array<$scalar, ArrayImplVt, 2>, - singular_values: &mut [<$scalar as Scalar>::Real], + singular_values: &mut [<$scalar as RlstScalar>::Real], mode: SvdMode, ) -> RlstResult<()> { let lwork: i32 = -1; @@ -267,7 +267,7 @@ macro_rules! impl_svd_complex { fn into_singular_values_alloc( mut self, - singular_values: &mut [<$scalar as Scalar>::Real], + singular_values: &mut [<$scalar as RlstScalar>::Real], ) -> RlstResult<()> { let lwork: i32 = -1; let mut work = [<$scalar as Zero>::zero(); 1]; @@ -285,7 +285,8 @@ macro_rules! impl_svd_complex { let ldvt = 1; let mut info = 0; - let mut rwork = vec![<<$scalar as Scalar>::Real as Zero>::zero(); 5 * k as usize]; + let mut rwork = + vec![<<$scalar as RlstScalar>::Real as Zero>::zero(); 5 * k as usize]; unsafe { $gesvd( @@ -354,7 +355,7 @@ macro_rules! impl_svd_complex { mut self, mut u: Array<$scalar, ArrayImplU, 2>, mut vt: Array<$scalar, ArrayImplVt, 2>, - singular_values: &mut [<$scalar as Scalar>::Real], + singular_values: &mut [<$scalar as RlstScalar>::Real], mode: SvdMode, ) -> RlstResult<()> { let lwork: i32 = -1; @@ -395,7 +396,8 @@ macro_rules! impl_svd_complex { let ldvt = vt.stride()[1] as i32; let mut info = 0; - let mut rwork = vec![<<$scalar as Scalar>::Real as Zero>::zero(); 5 * k as usize]; + let mut rwork = + vec![<<$scalar as RlstScalar>::Real as Zero>::zero(); 5 * k as usize]; unsafe { $gesvd( @@ -491,7 +493,7 @@ mod test { [](5, 10, SvdMode::Full, $tol); } - fn [](m: usize, n: usize, tol: <$scalar as Scalar>::Real) { + fn [](m: usize, n: usize, tol: <$scalar as RlstScalar>::Real) { let k = std::cmp::min(m, n); let mut mat = rlst_dynamic_array2!($scalar, [m, m]); let mut q = rlst_dynamic_array2!($scalar, [m, m]); @@ -502,12 +504,12 @@ mod test { qr.get_q_alloc(q.view_mut()).unwrap(); for index in 0..k { - sigma[[index, index]] = ((k - index) as <$scalar as Scalar>::Real).into(); + sigma[[index, index]] = ((k - index) as <$scalar as RlstScalar>::Real).into(); } let a = empty_array::<$scalar, 2>().simple_mult_into_resize(q.view(), sigma.view()); - let mut singvals = rlst_dynamic_array1!(<$scalar as Scalar>::Real, [k]); + let mut singvals = rlst_dynamic_array1!(<$scalar as RlstScalar>::Real, [k]); a.into_singular_values_alloc(singvals.data_mut()).unwrap(); @@ -516,7 +518,7 @@ mod test { } } - fn [](m: usize, n: usize, mode: SvdMode, tol: <$scalar as Scalar>::Real) { + fn [](m: usize, n: usize, mode: SvdMode, tol: <$scalar as RlstScalar>::Real) { let k = std::cmp::min(m, n); let mut mat_u; @@ -552,7 +554,7 @@ mod test { qr.get_q_alloc(vt.view_mut()).unwrap(); for index in 0..k { - sigma[[index, index]] = ((k - index) as <$scalar as Scalar>::Real).into(); + sigma[[index, index]] = ((k - index) as <$scalar as RlstScalar>::Real).into(); } let a = empty_array::<$scalar, 2>().simple_mult_into_resize( @@ -567,7 +569,7 @@ mod test { vt.set_zero(); sigma.set_zero(); - let mut singvals = rlst_dynamic_array1!(<$scalar as Scalar>::Real, [k]); + let mut singvals = rlst_dynamic_array1!(<$scalar as RlstScalar>::Real, [k]); a.into_svd_alloc(u.view_mut(), vt.view_mut(), singvals.data_mut(), mode) .unwrap(); diff --git a/dense/src/matrix_multiply.rs b/dense/src/matrix_multiply.rs index 7d759f06..bef0d127 100644 --- a/dense/src/matrix_multiply.rs +++ b/dense/src/matrix_multiply.rs @@ -3,13 +3,13 @@ //! This module implements the matrix multiplication. The current implementation //! uses the [rlst-blis] crate. +use crate::gemm::Gemm; use crate::traits::*; -use rlst_blis::interface::gemm::Gemm; -use rlst_blis::interface::types::TransMode; -use rlst_common::types::Scalar; +use crate::types::RlstScalar; +use crate::types::TransMode; pub fn matrix_multiply< - Item: Scalar + Gemm, + Item: RlstScalar + Gemm, MatA: RawAccess + Shape<2> + Stride<2>, MatB: RawAccess + Shape<2> + Stride<2>, MatC: RawAccessMut + Shape<2> + Stride<2>, @@ -75,7 +75,7 @@ pub fn matrix_multiply< mod test { use crate::assert_array_relative_eq; - use rlst_common::types::{c32, c64}; + use crate::types::{c32, c64}; use super::*; use crate::array::Array; @@ -103,19 +103,19 @@ mod test { matrix_multiply( transa, transb, - <$ScalarType as Scalar>::from_real(1.), + <$ScalarType as RlstScalar>::from_real(1.), &mat_a, &mat_b, - <$ScalarType as Scalar>::from_real(1.), + <$ScalarType as RlstScalar>::from_real(1.), &mut mat_c, ); matrix_multiply_compare( transa, transb, - <$ScalarType as Scalar>::from_real(1.), + <$ScalarType as RlstScalar>::from_real(1.), &mat_a, &mat_b, - <$ScalarType as Scalar>::from_real(1.), + <$ScalarType as RlstScalar>::from_real(1.), &mut expected, ); @@ -145,7 +145,7 @@ mod test { }; } - fn matrix_multiply_compare( + fn matrix_multiply_compare( transa: TransMode, transb: TransMode, alpha: Item, diff --git a/dense/src/tools.rs b/dense/src/tools.rs index d67b0226..1f2e0cc7 100644 --- a/dense/src/tools.rs +++ b/dense/src/tools.rs @@ -1,25 +1,25 @@ //! Useful library tools. use crate::traits::*; +use crate::types::*; use rand::prelude::*; use rand_distr::Distribution; -use rlst_common::types::*; /// This trait implements a simple convenient function to return random scalars /// from a given random number generator and distribution. For complex types the /// generator and distribution are separately applied to obtain the real and imaginary /// part of the random number. -pub trait RandScalar: Scalar { +pub trait RandScalar: RlstScalar { /// Returns a random number from a given random number generator `rng` and associated /// distribution `dist`. - fn random_scalar::Real>>( + fn random_scalar::Real>>( rng: &mut R, dist: &D, ) -> Self; } impl RandScalar for f32 { - fn random_scalar::Real>>( + fn random_scalar::Real>>( rng: &mut R, dist: &D, ) -> Self { @@ -28,7 +28,7 @@ impl RandScalar for f32 { } impl RandScalar for f64 { - fn random_scalar::Real>>( + fn random_scalar::Real>>( rng: &mut R, dist: &D, ) -> Self { @@ -37,7 +37,7 @@ impl RandScalar for f64 { } impl RandScalar for c32 { - fn random_scalar::Real>>( + fn random_scalar::Real>>( rng: &mut R, dist: &D, ) -> Self { @@ -46,7 +46,7 @@ impl RandScalar for c32 { } impl RandScalar for c64 { - fn random_scalar::Real>>( + fn random_scalar::Real>>( rng: &mut R, dist: &D, ) -> Self { @@ -127,7 +127,7 @@ macro_rules! assert_matrix_ulps_eq { }}; } -pub trait PrettyPrint { +pub trait PrettyPrint { fn pretty_print(&self); fn pretty_print_with_dimension(&self, rows: usize, cols: usize); fn pretty_print_advanced( @@ -215,7 +215,7 @@ pretty_print_impl!(c32, fmt_complex); pretty_print_impl!(c64, fmt_complex); // https://stackoverflow.com/a/65266882 -fn fmt_real(num: T::Real, width: usize, precision: usize, exp_pad: usize) -> String { +fn fmt_real(num: T::Real, width: usize, precision: usize, exp_pad: usize) -> String { let mut num = format!("{:.precision$e}", num, precision = precision); // Safe to `unwrap` as `num` is guaranteed to contain `'e'` let exp = num.split_off(num.find('e').unwrap()); @@ -229,7 +229,7 @@ fn fmt_real(num: T::Real, width: usize, precision: usize, exp_pad: us format!("{:>width$}", num, width = width) } -fn fmt_complex(num: T, width: usize, precision: usize, exp_pad: usize) -> String { +fn fmt_complex(num: T, width: usize, precision: usize, exp_pad: usize) -> String { let sign = if num.im() < ::zero() { "-" } else { diff --git a/dense/src/traits.rs b/dense/src/traits.rs index 67d4b9d5..848224dd 100644 --- a/dense/src/traits.rs +++ b/dense/src/traits.rs @@ -8,8 +8,8 @@ pub use crate::linalg::{ }; pub use accessors::*; -use rlst_blis::interface::types::TransMode; -use rlst_common::types::*; +use crate::types::TransMode; +use crate::types::*; /// Return the shape of the object. pub trait Shape { @@ -44,7 +44,7 @@ pub trait ResizeInPlace { /// Multiply First * Second and sum into Self pub trait MultInto { - type Item: Scalar; + type Item: RlstScalar; fn simple_mult_into(self, arr_a: First, arr_b: Second) -> Self where Self: Sized, @@ -71,7 +71,7 @@ pub trait MultInto { /// Multiply First * Second and sum into Self. Allow to resize Self if necessary pub trait MultIntoResize { - type Item: Scalar; + type Item: RlstScalar; fn simple_mult_into_resize(self, arr_a: First, arr_b: Second) -> Self where Self: Sized, @@ -98,7 +98,7 @@ pub trait MultIntoResize { /// Default iterator. pub trait DefaultIterator { - type Item: Scalar; + type Item: RlstScalar; type Iter<'a>: std::iter::Iterator where Self: 'a; @@ -108,7 +108,7 @@ pub trait DefaultIterator { /// Mutable default iterator. pub trait DefaultIteratorMut { - type Item: Scalar; + type Item: RlstScalar; type IterMut<'a>: std::iter::Iterator where Self: 'a; @@ -120,7 +120,7 @@ pub trait DefaultIteratorMut { /// `i` is row, `j` is column, and `data` is the corresponding /// element. pub trait AijIterator { - type Item: Scalar; + type Item: RlstScalar; type Iter<'a>: std::iter::Iterator where Self: 'a; diff --git a/dense/src/traits/accessors.rs b/dense/src/traits/accessors.rs index eaa27689..3f95e33c 100644 --- a/dense/src/traits/accessors.rs +++ b/dense/src/traits/accessors.rs @@ -17,11 +17,11 @@ //! To get raw access to the underlying data use the [`RawAccess`] and [`RawAccessMut`] traits. use crate::traits::Shape; -use rlst_common::types::{DataChunk, Scalar}; +use crate::types::{DataChunk, RlstScalar}; /// This trait provides unsafe access by value to the underlying data. pub trait UnsafeRandomAccessByValue { - type Item: Scalar; + type Item: RlstScalar; /// Return the element at position determined by `multi_index`. /// @@ -32,7 +32,7 @@ pub trait UnsafeRandomAccessByValue { /// This trait provides unsafe access by reference to the underlying data. pub trait UnsafeRandomAccessByRef { - type Item: Scalar; + type Item: RlstScalar; /// Return a mutable reference to the element at position determined by `multi_index`. /// @@ -43,7 +43,7 @@ pub trait UnsafeRandomAccessByRef { /// This trait provides unsafe mutable access to the underlying data. pub trait UnsafeRandomAccessMut { - type Item: Scalar; + type Item: RlstScalar; /// Return a mutable reference to the element at position determined by `multi_index`. /// @@ -72,13 +72,13 @@ pub trait RandomAccessMut: UnsafeRandomAccessMut { /// Return chunks of data of size N; pub trait ChunkedAccess { - type Item: Scalar; + type Item: RlstScalar; fn get_chunk(&self, chunk_index: usize) -> Option>; } /// Get raw access to the underlying data. pub trait RawAccess { - type Item: Scalar; + type Item: RlstScalar; /// Get a slice of the whole data. fn data(&self) -> &[Self::Item]; @@ -100,7 +100,7 @@ fn check_dimension(multi_index: [usize; NDIM], shape: [usize; } impl< - Item: Scalar, + Item: RlstScalar, Mat: UnsafeRandomAccessByValue + Shape, const NDIM: usize, > RandomAccessByValue for Mat @@ -115,7 +115,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, Mat: UnsafeRandomAccessMut + Shape, const NDIM: usize, > RandomAccessMut for Mat @@ -130,7 +130,7 @@ impl< } impl< - Item: Scalar, + Item: RlstScalar, Mat: UnsafeRandomAccessByRef + Shape, const NDIM: usize, > RandomAccessByRef for Mat diff --git a/dense/src/types.rs b/dense/src/types.rs new file mode 100644 index 00000000..9f09aea7 --- /dev/null +++ b/dense/src/types.rs @@ -0,0 +1,379 @@ +//! Basic types. + +use thiserror::Error; + +/// The Rlst error type. +#[derive(Error, Debug)] +pub enum RlstError { + #[error("Method {0} is not implemented.")] + NotImplemented(String), + #[error("Operation {0} failed.")] + OperationFailed(String), + #[error("Matrix has empty dimension {0:#?}.")] + MatrixIsEmpty((usize, usize)), + #[error("Dimension mismatch. Expected {expected:}. Actual {actual:}")] + SingleDimensionError { expected: usize, actual: usize }, + #[error("Index Layout error: {0}")] + IndexLayoutError(String), + #[error("MPI Rank does not exist. {0}")] + MpiRankError(i32), + #[error("Incompatible stride for Lapack.")] + IncompatibleStride, + #[error("Lapack error: {0}")] + LapackError(i32), + #[error("{0}")] + GeneralError(String), + #[error("I/O Error: {0}")] + IoError(String), + #[error("Umfpack Error Code: {0}")] + UmfpackError(i32), + #[error("Matrix is not square. Dimension: {0}x{1}")] + MatrixNotSquare(usize, usize), + #[error("Matrix is not Hermitian (complex conjugate symmetric).")] + MatrixNotHermitian, +} + +/// Alias for an Rlst Result type. +pub type RlstResult = std::result::Result; + +/// Data chunk of fixed size N. +/// The field `valid_entries` stores how many entries of the chunk +/// contain valid data. +pub struct DataChunk { + pub data: [Item; N], + pub start_index: usize, + pub valid_entries: usize, +} + +use num::complex::Complex; +use num::traits::{Float, FromPrimitive, NumAssign, NumCast, NumOps, ToPrimitive, Zero}; +use rand::{distributions::Standard, prelude::*}; +use serde::{Deserialize, Serialize}; +use std::fmt::{Debug, Display, LowerExp, UpperExp}; +use std::iter::{Product, Sum}; +use std::ops::Neg; + +pub use num::complex::Complex32 as c32; +pub use num::complex::Complex64 as c64; + +use crate::gemm::Gemm; + +// The following RlstScalar trait and is implementation for f32, f64, c32, c64 +// is a modifed version of the Scalar trait from the Rust Cauchy package +// (https://github.com/rust-math/cauchy), which is MIT licensed. For the full license text see +// https://github.com/rust-math/cauchy/blob/master/LICENSE. +pub trait RlstScalar: + NumAssign + + FromPrimitive + + NumCast + + Neg + + Copy + + Clone + + Display + + Debug + + LowerExp + + UpperExp + + Sum + + Product + + Serialize + + Send + + Sync + + Gemm + + for<'de> Deserialize<'de> + + 'static +{ + type Real: RlstScalar + + NumOps + + Float; + type Complex: RlstScalar + + NumOps + + NumOps; + + /// Create a new real number + fn real(re: T) -> Self::Real; + /// Create a new complex number + fn complex(re: T, im: T) -> Self::Complex; + + fn from_real(re: Self::Real) -> Self; + + fn add_real(self, re: Self::Real) -> Self; + fn sub_real(self, re: Self::Real) -> Self; + fn mul_real(self, re: Self::Real) -> Self; + fn div_real(self, re: Self::Real) -> Self; + + fn add_complex(self, im: Self::Complex) -> Self::Complex; + fn sub_complex(self, im: Self::Complex) -> Self::Complex; + fn mul_complex(self, im: Self::Complex) -> Self::Complex; + fn div_complex(self, im: Self::Complex) -> Self::Complex; + + fn pow(self, n: Self) -> Self; + fn powi(self, n: i32) -> Self; + fn powf(self, n: Self::Real) -> Self; + fn powc(self, n: Self::Complex) -> Self::Complex; + + /// Real part + fn re(&self) -> Self::Real; + /// Imaginary part + fn im(&self) -> Self::Real; + /// As a complex number + fn as_c(&self) -> Self::Complex; + /// Complex conjugate + fn conj(&self) -> Self; + + /// Absolute value + fn abs(self) -> Self::Real; + /// Sqaure of absolute value + fn square(self) -> Self::Real; + + fn sqrt(self) -> Self; + fn exp(self) -> Self; + fn ln(self) -> Self; + fn sin(self) -> Self; + fn cos(self) -> Self; + fn tan(self) -> Self; + fn asin(self) -> Self; + fn acos(self) -> Self; + fn atan(self) -> Self; + fn sinh(self) -> Self; + fn cosh(self) -> Self; + fn tanh(self) -> Self; + fn asinh(self) -> Self; + fn acosh(self) -> Self; + fn atanh(self) -> Self; + + /// Generate an random number from + /// [rand::distributions::Standard](https://docs.rs/rand/0.7.2/rand/distributions/struct.Standard.html) + fn rand(rng: &mut impl Rng) -> Self; +} + +macro_rules! impl_float { + ($name:ident) => { + #[inline] + fn $name(self) -> Self { + Float::$name(self) + } + }; +} + +macro_rules! impl_complex { + ($name:ident) => { + #[inline] + fn $name(self) -> Self { + Complex::$name(self) + } + }; +} + +macro_rules! impl_with_real { + ($name:ident, $op:tt) => { + #[inline] + fn $name(self, re: Self::Real) -> Self { + self $op re + } + } +} + +macro_rules! impl_with_complex { + ($name:ident, $op:tt) => { + #[inline] + fn $name(self, im: Self::Complex) -> Self::Complex { + self $op im + } + } +} + +macro_rules! impl_scalar { + ($real:ty, $complex:ty) => { + impl RlstScalar for $real { + type Real = $real; + type Complex = $complex; + + #[inline] + fn re(&self) -> Self::Real { + *self + } + #[inline] + fn im(&self) -> Self::Real { + 0.0 + } + + #[inline] + fn from_real(re: Self::Real) -> Self { + re + } + + fn pow(self, n: Self) -> Self { + self.powf(n) + } + fn powi(self, n: i32) -> Self { + Float::powi(self, n) + } + fn powf(self, n: Self::Real) -> Self { + Float::powf(self, n) + } + fn powc(self, n: Self::Complex) -> Self::Complex { + self.as_c().powc(n) + } + + #[inline] + fn real(re: T) -> Self::Real { + NumCast::from(re).unwrap() + } + #[inline] + fn complex(re: T, im: T) -> Self::Complex { + Complex { + re: NumCast::from(re).unwrap(), + im: NumCast::from(im).unwrap(), + } + } + #[inline] + fn as_c(&self) -> Self::Complex { + Complex::new(*self, 0.0) + } + #[inline] + fn conj(&self) -> Self { + *self + } + #[inline] + fn square(self) -> Self::Real { + self * self + } + + fn rand(rng: &mut impl Rng) -> Self { + rng.sample(Standard) + } + + impl_with_real!(add_real, +); + impl_with_real!(sub_real, -); + impl_with_real!(mul_real, *); + impl_with_real!(div_real, /); + impl_with_complex!(add_complex, +); + impl_with_complex!(sub_complex, -); + impl_with_complex!(mul_complex, *); + impl_with_complex!(div_complex, /); + + impl_float!(sqrt); + impl_float!(abs); + impl_float!(exp); + impl_float!(ln); + impl_float!(sin); + impl_float!(cos); + impl_float!(tan); + impl_float!(sinh); + impl_float!(cosh); + impl_float!(tanh); + impl_float!(asin); + impl_float!(acos); + impl_float!(atan); + impl_float!(asinh); + impl_float!(acosh); + impl_float!(atanh); + } + + impl RlstScalar for $complex { + type Real = $real; + type Complex = $complex; + + #[inline] + fn re(&self) -> Self::Real { + self.re + } + #[inline] + fn im(&self) -> Self::Real { + self.im + } + + #[inline] + fn from_real(re: Self::Real) -> Self { + Self::new(re, Zero::zero()) + } + + fn pow(self, n: Self) -> Self { + self.powc(n) + } + fn powi(self, n: i32) -> Self { + self.powf(n as Self::Real) + } + fn powf(self, n: Self::Real) -> Self { + self.powf(n) + } + fn powc(self, n: Self::Complex) -> Self::Complex { + self.powc(n) + } + + #[inline] + fn real(re: T) -> Self::Real { + NumCast::from(re).unwrap() + } + #[inline] + fn complex(re: T, im: T) -> Self::Complex { + Complex { + re: NumCast::from(re).unwrap(), + im: NumCast::from(im).unwrap(), + } + } + #[inline] + fn as_c(&self) -> Self::Complex { + *self + } + #[inline] + fn conj(&self) -> Self { + Complex::conj(self) + } + #[inline] + fn square(self) -> Self::Real { + Complex::norm_sqr(&self) + } + #[inline] + fn abs(self) -> Self::Real { + Complex::norm(self) + } + + fn rand(rng: &mut impl Rng) -> Self { + rng.sample(Standard) + } + + impl_with_real!(add_real, +); + impl_with_real!(sub_real, -); + impl_with_real!(mul_real, *); + impl_with_real!(div_real, /); + impl_with_complex!(add_complex, +); + impl_with_complex!(sub_complex, -); + impl_with_complex!(mul_complex, *); + impl_with_complex!(div_complex, /); + + impl_complex!(sqrt); + impl_complex!(exp); + impl_complex!(ln); + impl_complex!(sin); + impl_complex!(cos); + impl_complex!(tan); + impl_complex!(sinh); + impl_complex!(cosh); + impl_complex!(tanh); + impl_complex!(asin); + impl_complex!(acos); + impl_complex!(atan); + impl_complex!(asinh); + impl_complex!(acosh); + impl_complex!(atanh); + } + } +} + +impl_scalar!(f32, c32); +impl_scalar!(f64, c64); + +/// Transposition Mode. +#[derive(Clone, Copy, PartialEq)] +#[repr(u32)] +pub enum TransMode { + /// Complex conjugate of matrix. + ConjNoTrans, + /// No modification of matrix. + NoTrans, + /// Transposition of matrix. + Trans, + /// Conjugate transpose of matrix. + ConjTrans, +} diff --git a/io/Cargo.toml b/io/Cargo.toml index 593f16f3..8c735e91 100644 --- a/io/Cargo.toml +++ b/io/Cargo.toml @@ -19,14 +19,12 @@ name = "rlst_io" [dependencies] num = "0.4" -cauchy = "0.4" ndarray = "0.15" rand = "0.8" itertools = "0.12" rand_distr = "0.4" matrixmultiply = "0.3" approx = { version = "0.5", features=["num-complex"] } -rlst-common = {path = "../common"} rlst-sparse = {path = "../sparse"} rlst-dense = {path = "../dense"} diff --git a/io/src/matrix_market.rs b/io/src/matrix_market.rs index 04a0571e..47a60d47 100644 --- a/io/src/matrix_market.rs +++ b/io/src/matrix_market.rs @@ -5,10 +5,10 @@ //! files must use the `general` property. To write out matrices the functions [write_array_mm] and //! [write_coordinate_mm] are provided. -use rlst_common::types::Scalar; -use rlst_common::types::{RlstError, RlstResult}; use rlst_dense::array::DynamicArray; use rlst_dense::traits::RawAccessMut; +use rlst_dense::types::RlstScalar; +use rlst_dense::types::{RlstError, RlstResult}; use rlst_sparse::sparse::csr_mat::CsrMatrix; use std::fs::File; use std::io::Write; @@ -52,11 +52,11 @@ impl MmIdentifier for f64 { const MMTYPE: &'static str = "real"; } -impl MmIdentifier for rlst_common::types::c32 { +impl MmIdentifier for rlst_dense::types::c32 { const MMTYPE: &'static str = "complex"; } -impl MmIdentifier for rlst_common::types::c64 { +impl MmIdentifier for rlst_dense::types::c64 { const MMTYPE: &'static str = "complex"; } @@ -86,7 +86,7 @@ pub struct MatrixMarketInfo { /// [rlst_dense::traits::Shape] traits. Any object satisfying these traits can be written /// out with this function. pub fn write_coordinate_mm< - T: Scalar + MmIdentifier, + T: RlstScalar + MmIdentifier, Mat: rlst_dense::traits::AijIterator + rlst_dense::traits::Shape<2>, >( mat: &Mat, @@ -121,7 +121,7 @@ pub fn write_coordinate_mm< /// [rlst_dense::traits::Shape] traits. Any object satisfying these traits can be written /// out with this function. pub fn write_array_mm< - T: Scalar + MmIdentifier, + T: RlstScalar + MmIdentifier, Mat: rlst_dense::traits::DefaultIterator + rlst_dense::traits::Shape<2>, >( mat: &Mat, @@ -147,7 +147,7 @@ pub fn write_array_mm< /// /// The function returns a [DynamicArray] object representing the data in the file. /// Currently only `general` matrices are supported without special symmetry. -pub fn read_array_mm(fname: &str) -> RlstResult> { +pub fn read_array_mm(fname: &str) -> RlstResult> { let mut reader = open_file(fname).unwrap(); let mm_info = parse_header(&mut reader).unwrap(); @@ -209,7 +209,7 @@ pub fn read_array_mm(fname: &str) -> RlstResult> { /// Returns a [rlst_sparse::sparse::csr_mat::CsrMatrix] sparse matrix object representing /// the data in the file. /// Currently only `general` matrices are supported without special symmetry. -pub fn read_coordinate_mm(fname: &str) -> RlstResult> { +pub fn read_coordinate_mm(fname: &str) -> RlstResult> { let mut reader = open_file(fname).unwrap(); let mm_info = parse_header(&mut reader).unwrap(); @@ -344,7 +344,7 @@ fn parse_header(reader: &mut io::Lines>) -> RlstResult( +fn parse_array( reader: &mut io::Lines>, buf: &mut [T], nelems: usize, @@ -395,7 +395,7 @@ fn parse_array( } /// Parse coordinate information. -fn parse_coordinate( +fn parse_coordinate( reader: &mut io::Lines>, rows: &mut [usize], cols: &mut [usize], diff --git a/lapack/Cargo.toml b/lapack/Cargo.toml index 8500c47e..4e1af9f0 100644 --- a/lapack/Cargo.toml +++ b/lapack/Cargo.toml @@ -19,7 +19,6 @@ name = "rlst_lapack" [dependencies] num = "0.4" -cauchy = "0.4" rand = "0.8" itertools = "0.12" rand_distr = "0.4" diff --git a/operator/Cargo.toml b/operator/Cargo.toml index bed9af94..8981d11e 100644 --- a/operator/Cargo.toml +++ b/operator/Cargo.toml @@ -9,12 +9,10 @@ edition = "2021" # See more keys and their definitions at https://doc.rust-lang.org/cargo/reference/manifest.html [dependencies] -rlst-common = {path = "../common"} rlst-sparse = {path = "../sparse"} rlst-dense = {path = "../dense"} rlst-blis = {path = "../blis"} num = "0.4" -cauchy = "0.4" thiserror = "1.0.40" [dev-dependencies] diff --git a/operator/src/interface/array_vector_space.rs b/operator/src/interface/array_vector_space.rs index 65fe32ed..07fc8414 100644 --- a/operator/src/interface/array_vector_space.rs +++ b/operator/src/interface/array_vector_space.rs @@ -3,7 +3,7 @@ use std::marker::PhantomData; use crate::space::*; -use rlst_common::types::Scalar; +use rlst_dense::types::RlstScalar; use rlst_dense::{ array::{ views::{ArrayView, ArrayViewMut}, @@ -14,16 +14,16 @@ use rlst_dense::{ rlst_dynamic_array1, }; -pub struct ArrayVectorSpace { +pub struct ArrayVectorSpace { dimension: usize, _marker: PhantomData, } -pub struct ArrayVectorSpaceElement { +pub struct ArrayVectorSpaceElement { elem: DynamicArray, } -impl ArrayVectorSpaceElement { +impl ArrayVectorSpaceElement { pub fn new(space: &ArrayVectorSpace) -> Self { Self { elem: rlst_dynamic_array1!(Item, [space.dimension()]), @@ -31,7 +31,7 @@ impl ArrayVectorSpaceElement { } } -impl ArrayVectorSpace { +impl ArrayVectorSpace { pub fn new(dimension: usize) -> Self { Self { dimension, @@ -40,13 +40,13 @@ impl ArrayVectorSpace { } } -impl IndexableSpace for ArrayVectorSpace { +impl IndexableSpace for ArrayVectorSpace { fn dimension(&self) -> usize { self.dimension } } -impl LinearSpace for ArrayVectorSpace { +impl LinearSpace for ArrayVectorSpace { type E = ArrayVectorSpaceElement; type F = Item; @@ -56,13 +56,13 @@ impl LinearSpace for ArrayVectorSpace { } } -impl InnerProductSpace for ArrayVectorSpace { +impl InnerProductSpace for ArrayVectorSpace { fn inner(&self, x: &Self::E, other: &Self::E) -> Self::F { x.view().inner(other.view()) } } -impl Element for ArrayVectorSpaceElement { +impl Element for ArrayVectorSpaceElement { type F = Item; type Space = ArrayVectorSpace; diff --git a/operator/src/interface/dense_matrix_operator.rs b/operator/src/interface/dense_matrix_operator.rs index 75153eb4..a05692ba 100644 --- a/operator/src/interface/dense_matrix_operator.rs +++ b/operator/src/interface/dense_matrix_operator.rs @@ -1,6 +1,5 @@ use crate::{space::*, AsApply, OperatorBase}; -use rlst_blis::interface::gemm::Gemm; -use rlst_common::types::Scalar; +use rlst_dense::types::RlstScalar; use rlst_dense::{ array::Array, traits::{MultInto, RawAccess, Shape, Stride, UnsafeRandomAccessByValue}, @@ -10,7 +9,7 @@ use super::array_vector_space::ArrayVectorSpace; pub struct DenseMatrixOperator< 'a, - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue<2, Item = Item> + Shape<2> + Stride<2> + RawAccess, > { arr: Array, @@ -20,7 +19,7 @@ pub struct DenseMatrixOperator< impl< 'a, - Item: Scalar, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue<2, Item = Item> + Shape<2> + Stride<2> + RawAccess, > std::fmt::Debug for DenseMatrixOperator<'a, Item, ArrayImpl> { @@ -33,7 +32,7 @@ impl< impl< 'a, - Item: Scalar + Gemm, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue<2, Item = Item> + Shape<2> + Stride<2> @@ -55,7 +54,7 @@ impl< impl< 'a, - Item: Scalar + Gemm, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue<2, Item = Item> + Shape<2> + Stride<2> + RawAccess, > OperatorBase for DenseMatrixOperator<'a, Item, ArrayImpl> { @@ -73,7 +72,7 @@ impl< impl< 'a, - Item: Scalar + Gemm, + Item: RlstScalar, ArrayImpl: UnsafeRandomAccessByValue<2, Item = Item> + Shape<2> + Stride<2> + RawAccess, > AsApply for DenseMatrixOperator<'a, Item, ArrayImpl> { @@ -83,7 +82,7 @@ impl< x: &::E, beta: ::F, y: &mut ::E, - ) -> rlst_common::types::RlstResult<()> { + ) -> rlst_dense::types::RlstResult<()> { use rlst_dense::array::mult_into::TransMode; y.view_mut().mult_into( TransMode::NoTrans, diff --git a/operator/src/interface/sparse_operator.rs b/operator/src/interface/sparse_operator.rs index c0065306..16b2cc63 100644 --- a/operator/src/interface/sparse_operator.rs +++ b/operator/src/interface/sparse_operator.rs @@ -1,18 +1,18 @@ use crate::{space::*, AsApply, OperatorBase}; -use rlst_common::types::Scalar; use rlst_dense::traits::{RawAccess, RawAccessMut, Shape}; +use rlst_dense::types::RlstScalar; use rlst_sparse::sparse::csc_mat::CscMatrix; use rlst_sparse::sparse::csr_mat::CsrMatrix; use super::array_vector_space::ArrayVectorSpace; -pub struct CsrMatrixOperator<'a, Item: Scalar> { +pub struct CsrMatrixOperator<'a, Item: RlstScalar> { csr_mat: &'a CsrMatrix, domain: &'a ArrayVectorSpace, range: &'a ArrayVectorSpace, } -impl std::fmt::Debug for CsrMatrixOperator<'_, Item> { +impl std::fmt::Debug for CsrMatrixOperator<'_, Item> { fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result { f.debug_struct("CsrMatrixOperator") .field("Dimension", &self.csr_mat.shape()) @@ -21,7 +21,7 @@ impl std::fmt::Debug for CsrMatrixOperator<'_, Item> { } } -impl<'a, Item: Scalar> CsrMatrixOperator<'a, Item> { +impl<'a, Item: RlstScalar> CsrMatrixOperator<'a, Item> { pub fn new( csr_mat: &'a CsrMatrix, domain: &'a ArrayVectorSpace, @@ -38,7 +38,7 @@ impl<'a, Item: Scalar> CsrMatrixOperator<'a, Item> { } } -impl OperatorBase for CsrMatrixOperator<'_, Item> { +impl OperatorBase for CsrMatrixOperator<'_, Item> { type Domain = ArrayVectorSpace; type Range = ArrayVectorSpace; @@ -52,27 +52,27 @@ impl OperatorBase for CsrMatrixOperator<'_, Item> { } } -impl AsApply for CsrMatrixOperator<'_, Item> { +impl AsApply for CsrMatrixOperator<'_, Item> { fn apply_extended( &self, alpha: ::F, x: &::E, beta: ::F, y: &mut ::E, - ) -> rlst_common::types::RlstResult<()> { + ) -> rlst_dense::types::RlstResult<()> { self.csr_mat .matmul(alpha, x.view().data(), beta, y.view_mut().data_mut()); Ok(()) } } -pub struct CscMatrixOperator<'a, Item: Scalar> { +pub struct CscMatrixOperator<'a, Item: RlstScalar> { csc_mat: &'a CscMatrix, domain: &'a ArrayVectorSpace, range: &'a ArrayVectorSpace, } -impl std::fmt::Debug for CscMatrixOperator<'_, Item> { +impl std::fmt::Debug for CscMatrixOperator<'_, Item> { fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result { f.debug_struct("CscMatrixOperator") .field("Dimension", &self.csc_mat.shape()) @@ -81,7 +81,7 @@ impl std::fmt::Debug for CscMatrixOperator<'_, Item> { } } -impl<'a, Item: Scalar> CscMatrixOperator<'a, Item> { +impl<'a, Item: RlstScalar> CscMatrixOperator<'a, Item> { pub fn new( csc_mat: &'a CscMatrix, domain: &'a ArrayVectorSpace, @@ -98,7 +98,7 @@ impl<'a, Item: Scalar> CscMatrixOperator<'a, Item> { } } -impl OperatorBase for CscMatrixOperator<'_, Item> { +impl OperatorBase for CscMatrixOperator<'_, Item> { type Domain = ArrayVectorSpace; type Range = ArrayVectorSpace; @@ -112,14 +112,14 @@ impl OperatorBase for CscMatrixOperator<'_, Item> { } } -impl AsApply for CscMatrixOperator<'_, Item> { +impl AsApply for CscMatrixOperator<'_, Item> { fn apply_extended( &self, alpha: ::F, x: &::E, beta: ::F, y: &mut ::E, - ) -> rlst_common::types::RlstResult<()> { + ) -> rlst_dense::types::RlstResult<()> { self.csc_mat .matmul(alpha, x.view().data(), beta, y.view_mut().data_mut()); Ok(()) diff --git a/operator/src/operations/conjugate_gradients.rs b/operator/src/operations/conjugate_gradients.rs index b43b2d16..3f3c9699 100644 --- a/operator/src/operations/conjugate_gradients.rs +++ b/operator/src/operations/conjugate_gradients.rs @@ -1,7 +1,7 @@ //! Arnoldi Iteration use crate::{AsApply, Element, InnerProductSpace, LinearSpace, NormedSpace}; use num::One; -use rlst_common::types::Scalar; +use rlst_dense::types::RlstScalar; pub struct CgIteration<'a, Space: InnerProductSpace, Op: AsApply> { operator: &'a Op, @@ -9,7 +9,7 @@ pub struct CgIteration<'a, Space: InnerProductSpace, Op: AsApply::Real, + tol: ::Real, #[allow(clippy::type_complexity)] callable: Option::E, &::E) + 'a>>, print_debug: bool, @@ -25,7 +25,7 @@ impl<'a, Space: InnerProductSpace, Op: AsApply> rhs, x: op.domain().zero(), max_iter: 1000, - tol: num::cast::::Real>(1E-6).unwrap(), + tol: num::cast::::Real>(1E-6).unwrap(), callable: None, print_debug: false, } @@ -36,7 +36,7 @@ impl<'a, Space: InnerProductSpace, Op: AsApply> self } - pub fn set_tol(mut self, tol: ::Real) -> Self { + pub fn set_tol(mut self, tol: ::Real) -> Self { self.tol = tol; self } @@ -56,15 +56,15 @@ impl<'a, Space: InnerProductSpace, Op: AsApply> self } - pub fn run(mut self) -> (Space::E, ::Real) { - fn print_success(it_count: usize, rel_res: T) { + pub fn run(mut self) -> (Space::E, ::Real) { + fn print_success(it_count: usize, rel_res: T) { println!( "CG converged in {} iterations with relative residual {:+E}.", it_count, rel_res ); } - fn print_fail(it_count: usize, rel_res: T) { + fn print_fail(it_count: usize, rel_res: T) { println!( "CG did not converge in {} iterations. Final relative residual is {:+E}.", it_count, rel_res diff --git a/operator/src/operations/modified_gram_schmidt.rs b/operator/src/operations/modified_gram_schmidt.rs index de5a3283..eceb6ad2 100644 --- a/operator/src/operations/modified_gram_schmidt.rs +++ b/operator/src/operations/modified_gram_schmidt.rs @@ -1,7 +1,7 @@ //! Gram Schmidt orthogonalization use crate::{frame::Frame, Element, InnerProductSpace, NormedSpace}; use num::One; -use rlst_common::types::Scalar; +use rlst_dense::types::RlstScalar; use rlst_dense::{ array::DynamicArray, traits::{RandomAccessMut, Shape}, @@ -10,7 +10,7 @@ pub struct ModifiedGramSchmidt; impl ModifiedGramSchmidt { pub fn orthogonalize< - Item: Scalar, + Item: RlstScalar, Elem: Element, Space: InnerProductSpace, FrameType: Frame, @@ -32,7 +32,7 @@ impl ModifiedGramSchmidt { *r_mat.get_mut([other_index, elem_index]).unwrap() = inner; elem.axpy_inplace(-inner, other_elem); } - let norm = ::from_real(space.norm(&elem)); + let norm = ::from_real(space.norm(&elem)); *r_mat.get_mut([elem_index, elem_index]).unwrap() = norm; elem.scale_inplace(::one() / norm); frame.get_mut(elem_index).unwrap().fill_inplace(&elem); @@ -45,8 +45,8 @@ mod test { use approx::{assert_abs_diff_eq, assert_relative_eq}; use num::Zero; - use rlst_common::types::c64; use rlst_dense::rlst_dynamic_array2; + use rlst_dense::types::c64; use rlst_dense::{rlst_dynamic_array1, traits::RawAccess}; use crate::space::frame::VectorFrame; diff --git a/operator/src/operator.rs b/operator/src/operator.rs index 4fd73b0e..530591d6 100644 --- a/operator/src/operator.rs +++ b/operator/src/operator.rs @@ -2,7 +2,7 @@ use crate::{FieldType, LinearSpace}; use num::{One, Zero}; -use rlst_common::types::*; +use rlst_dense::types::*; use std::fmt::Debug; // A base operator trait. diff --git a/operator/src/space/dual_space.rs b/operator/src/space/dual_space.rs index b9654217..f9460490 100644 --- a/operator/src/space/dual_space.rs +++ b/operator/src/space/dual_space.rs @@ -1,5 +1,5 @@ use super::{ElementView, LinearSpace}; -use rlst_common::types::RlstResult; +use rlst_dense::types::RlstResult; pub trait DualSpace: LinearSpace { type Space: LinearSpace; diff --git a/operator/src/space/element.rs b/operator/src/space/element.rs index 8dc0534f..7c30c47a 100644 --- a/operator/src/space/element.rs +++ b/operator/src/space/element.rs @@ -1,5 +1,5 @@ use num::One; -use rlst_common::types::Scalar; +use rlst_dense::types::RlstScalar; use super::LinearSpace; @@ -9,7 +9,7 @@ pub trait Element { type Space: LinearSpace; - type F: Scalar; + type F: RlstScalar; type View<'b> where diff --git a/operator/src/space/linear_space.rs b/operator/src/space/linear_space.rs index 48828911..62b26f53 100644 --- a/operator/src/space/linear_space.rs +++ b/operator/src/space/linear_space.rs @@ -2,7 +2,7 @@ use super::Element; use num::One; -use rlst_common::types::Scalar; +use rlst_dense::types::RlstScalar; /// Definition of a linear space /// @@ -10,7 +10,7 @@ use rlst_common::types::Scalar; /// elements of the space. pub trait LinearSpace { /// Field Type. - type F: Scalar; + type F: RlstScalar; /// Type associated with elements of the space. type E: Element; diff --git a/operator/src/space/normed_space.rs b/operator/src/space/normed_space.rs index bffb61ab..3cb469a8 100644 --- a/operator/src/space/normed_space.rs +++ b/operator/src/space/normed_space.rs @@ -2,15 +2,15 @@ use crate::ElementType; use crate::InnerProductSpace; use super::LinearSpace; -use rlst_common::types::Scalar; +use rlst_dense::types::RlstScalar; pub trait NormedSpace: LinearSpace { /// Norm of a vector. - fn norm(&self, x: &ElementType) -> ::Real; + fn norm(&self, x: &ElementType) -> ::Real; } impl NormedSpace for S { - fn norm(&self, x: &ElementType) -> ::Real { + fn norm(&self, x: &ElementType) -> ::Real { let abs_square = self.inner(x, x).abs(); abs_square.sqrt() } diff --git a/rlst/Cargo.toml b/rlst/Cargo.toml index 56d03890..85c92f66 100644 --- a/rlst/Cargo.toml +++ b/rlst/Cargo.toml @@ -19,7 +19,6 @@ name = "rlst" [dependencies] rlst-operator = {path = "../operator"} -rlst-common = {path = "../common"} rlst-sparse = {path = "../sparse"} rlst-dense = { path = "../dense"} rlst-blis = { path = "../blis"} diff --git a/rlst/examples/fixed_size_array.rs b/rlst/examples/fixed_size_array.rs index 01e6bdbc..3b6f6a00 100644 --- a/rlst/examples/fixed_size_array.rs +++ b/rlst/examples/fixed_size_array.rs @@ -1,7 +1,7 @@ //! Fixed size matrix example //! -use rlst_common::types::c64; +use rlst_dense::types::c64; use rlst_proc_macro::rlst_static_array; use rlst_proc_macro::rlst_static_type; diff --git a/rlst/tests/array_operations.rs b/rlst/tests/array_operations.rs index 23288330..3674831c 100644 --- a/rlst/tests/array_operations.rs +++ b/rlst/tests/array_operations.rs @@ -2,7 +2,7 @@ use approx::assert_relative_eq; use rlst::rlst_dynamic_array3; -use rlst_common::types::*; +use rlst_dense::types::*; use rlst_dense::{array::iterators::AsMultiIndex, layout::convert_1d_nd_from_shape}; use rlst_dense::{assert_array_relative_eq, traits::*}; diff --git a/sparse/Cargo.toml b/sparse/Cargo.toml index 02ec667e..7a53da12 100644 --- a/sparse/Cargo.toml +++ b/sparse/Cargo.toml @@ -6,7 +6,6 @@ edition = "2021" # See more keys and their definitions at https://doc.rust-lang.org/cargo/reference/manifest.html [dependencies] -rlst-common = {path = "../common"} rlst-dense = {path = "../dense"} mpi = { version="0.7.0", optional=true } mpi-sys = { version="0.2.1", optional=true } @@ -22,7 +21,6 @@ rlst-umfpack = { path = "../umfpack"} [dev-dependencies] cauchy = "0.4" float_eq = { version = "1", features = ["num"] } -rlst-io = { path = "../io"} rand_chacha = "0.3" [lib] diff --git a/sparse/src/distributed_vector.rs b/sparse/src/distributed_vector.rs index e357f9f0..f77ecf78 100644 --- a/sparse/src/distributed_vector.rs +++ b/sparse/src/distributed_vector.rs @@ -3,18 +3,18 @@ use crate::traits::index_layout::IndexLayout; use mpi::datatype::{Partition, PartitionMut}; use mpi::traits::*; -use rlst_common::types::{RlstResult, Scalar}; use rlst_dense::array::DynamicArray; +use rlst_dense::types::{RlstResult, RlstScalar}; use rlst_dense::{rlst_dynamic_array1, traits::*}; use crate::index_layout::DefaultMpiIndexLayout; -pub struct DistributedVector<'a, Item: Scalar + Equivalence, C: Communicator> { +pub struct DistributedVector<'a, Item: RlstScalar + Equivalence, C: Communicator> { index_layout: &'a DefaultMpiIndexLayout<'a, C>, local: DynamicArray, } -impl<'a, Item: Scalar + Equivalence, C: Communicator> DistributedVector<'a, Item, C> { +impl<'a, Item: RlstScalar + Equivalence, C: Communicator> DistributedVector<'a, Item, C> { pub fn new(index_layout: &'a DefaultMpiIndexLayout<'a, C>) -> Self { DistributedVector { index_layout, diff --git a/sparse/src/ghost_communicator.rs b/sparse/src/ghost_communicator.rs index 991d4baf..73f71f60 100644 --- a/sparse/src/ghost_communicator.rs +++ b/sparse/src/ghost_communicator.rs @@ -5,7 +5,7 @@ use std::os::raw::c_void; use mpi::topology::SimpleCommunicator; use mpi::traits::{AsRaw, Communicator, CommunicatorCollectives, Equivalence, FromRaw}; use mpi_sys; -use rlst_common::types::Scalar; +use rlst_dense::types::RlstScalar; use crate::traits::index_layout::IndexLayout; @@ -151,7 +151,7 @@ impl GhostCommunicator { } } - pub fn forward_send_ghosts( + pub fn forward_send_ghosts( &self, local_values: &[T], ghost_values: &mut [T], @@ -176,7 +176,7 @@ impl GhostCommunicator { ); } } - pub fn backward_send_ghosts( + pub fn backward_send_ghosts( &self, local_values: &mut [T], ghost_values: &[T], diff --git a/sparse/src/index_layout/default_mpi_index_layout.rs b/sparse/src/index_layout/default_mpi_index_layout.rs index d439ea92..99e4915c 100644 --- a/sparse/src/index_layout/default_mpi_index_layout.rs +++ b/sparse/src/index_layout/default_mpi_index_layout.rs @@ -1,6 +1,6 @@ use crate::traits::index_layout::IndexLayout; use mpi::traits::Communicator; -use rlst_common::types::RlstResult; +use rlst_dense::types::RlstResult; pub struct DefaultMpiIndexLayout<'a, C: Communicator> { size: usize, @@ -89,7 +89,7 @@ impl<'a, C: Communicator> IndexLayout for DefaultMpiIndexLayout<'a, C> { if rank < self.comm.size() as usize { Ok((self.counts[rank], self.counts[1 + rank])) } else { - Err(rlst_common::types::RlstError::MpiRankError(rank as i32)) + Err(rlst_dense::types::RlstError::MpiRankError(rank as i32)) } } diff --git a/sparse/src/index_layout/default_serial_index_layout.rs b/sparse/src/index_layout/default_serial_index_layout.rs index 46a781d7..430a275e 100644 --- a/sparse/src/index_layout/default_serial_index_layout.rs +++ b/sparse/src/index_layout/default_serial_index_layout.rs @@ -1,5 +1,5 @@ use crate::traits::index_layout::IndexLayout; -use rlst_common::types::{RlstError, RlstResult}; +use rlst_dense::types::{RlstError, RlstResult}; pub struct DefaultSerialIndexLayout { size: usize, diff --git a/sparse/src/sparse/csc_mat.rs b/sparse/src/sparse/csc_mat.rs index 601bb843..b08ec2d1 100644 --- a/sparse/src/sparse/csc_mat.rs +++ b/sparse/src/sparse/csc_mat.rs @@ -1,17 +1,17 @@ //! Definition of CSC matrices. use crate::sparse::SparseMatType; -use rlst_common::types::RlstResult; +use rlst_dense::types::RlstResult; use crate::sparse::csr_mat::CsrMatrix; use crate::sparse::tools::normalize_aij; use itertools::Itertools; -use rlst_common::types::Scalar; use rlst_dense::traits::AijIterator; use rlst_dense::traits::Shape; +use rlst_dense::types::RlstScalar; #[derive(Clone)] -pub struct CscMatrix { +pub struct CscMatrix { mat_type: SparseMatType, shape: [usize; 2], indices: Vec, @@ -19,7 +19,7 @@ pub struct CscMatrix { data: Vec, } -impl CscMatrix { +impl CscMatrix { pub fn new( shape: [usize; 2], indices: Vec, @@ -130,13 +130,13 @@ impl CscMatrix { } } -pub struct CscAijIterator<'a, Item: Scalar> { +pub struct CscAijIterator<'a, Item: RlstScalar> { mat: &'a CscMatrix, col: usize, pos: usize, } -impl<'a, Item: Scalar> CscAijIterator<'a, Item> { +impl<'a, Item: RlstScalar> CscAijIterator<'a, Item> { pub fn new(mat: &'a CscMatrix) -> Self { // We need to move the col pointer to the first col that has at least one element. @@ -150,7 +150,7 @@ impl<'a, Item: Scalar> CscAijIterator<'a, Item> { } } -impl<'a, Item: Scalar> std::iter::Iterator for CscAijIterator<'a, Item> { +impl<'a, Item: RlstScalar> std::iter::Iterator for CscAijIterator<'a, Item> { type Item = (usize, usize, Item); fn next(&mut self) -> Option { @@ -184,7 +184,7 @@ impl<'a, Item: Scalar> std::iter::Iterator for CscAijIterator<'a, Item> { } } -impl AijIterator for CscMatrix { +impl AijIterator for CscMatrix { type Item = Item; type Iter<'a> = CscAijIterator<'a, Item> where Self: 'a; @@ -193,7 +193,7 @@ impl AijIterator for CscMatrix { } } -impl rlst_dense::traits::Shape<2> for CscMatrix { +impl rlst_dense::traits::Shape<2> for CscMatrix { fn shape(&self) -> [usize; 2] { self.shape } diff --git a/sparse/src/sparse/csr_mat.rs b/sparse/src/sparse/csr_mat.rs index 411fdf7b..5c59c995 100644 --- a/sparse/src/sparse/csr_mat.rs +++ b/sparse/src/sparse/csr_mat.rs @@ -1,17 +1,17 @@ //! Definition of CSR matrices. use crate::sparse::SparseMatType; -use rlst_common::types::RlstResult; use rlst_dense::traits::AijIterator; +use rlst_dense::types::RlstResult; use crate::sparse::tools::normalize_aij; -use rlst_common::types::Scalar; use rlst_dense::traits::Shape; +use rlst_dense::types::RlstScalar; use super::csc_mat::CscMatrix; #[derive(Clone)] -pub struct CsrMatrix { +pub struct CsrMatrix { mat_type: SparseMatType, shape: [usize; 2], indices: Vec, @@ -19,7 +19,7 @@ pub struct CsrMatrix { data: Vec, } -impl CsrMatrix { +impl CsrMatrix { pub fn new( shape: [usize; 2], indices: Vec, @@ -130,13 +130,13 @@ impl CsrMatrix { } } -pub struct CsrAijIterator<'a, Item: Scalar> { +pub struct CsrAijIterator<'a, Item: RlstScalar> { mat: &'a CsrMatrix, row: usize, pos: usize, } -impl<'a, Item: Scalar> CsrAijIterator<'a, Item> { +impl<'a, Item: RlstScalar> CsrAijIterator<'a, Item> { pub fn new(mat: &'a CsrMatrix) -> Self { // We need to move the row pointer to the first row that has at least one element. @@ -150,7 +150,7 @@ impl<'a, Item: Scalar> CsrAijIterator<'a, Item> { } } -impl<'a, Item: Scalar> std::iter::Iterator for CsrAijIterator<'a, Item> { +impl<'a, Item: RlstScalar> std::iter::Iterator for CsrAijIterator<'a, Item> { type Item = (usize, usize, Item); fn next(&mut self) -> Option { @@ -182,7 +182,7 @@ impl<'a, Item: Scalar> std::iter::Iterator for CsrAijIterator<'a, Item> { } } -impl AijIterator for CsrMatrix { +impl AijIterator for CsrMatrix { type Item = Item; type Iter<'a> = CsrAijIterator<'a, Item> where Self: 'a; @@ -191,7 +191,7 @@ impl AijIterator for CsrMatrix { } } -impl Shape<2> for CsrMatrix { +impl Shape<2> for CsrMatrix { fn shape(&self) -> [usize; 2] { self.shape } diff --git a/sparse/src/sparse/mpi_csr_mat.rs b/sparse/src/sparse/mpi_csr_mat.rs index 784fabfa..85b0bf1e 100644 --- a/sparse/src/sparse/mpi_csr_mat.rs +++ b/sparse/src/sparse/mpi_csr_mat.rs @@ -11,11 +11,11 @@ use crate::traits::index_layout::IndexLayout; use crate::distributed_vector::DistributedVector; use mpi::traits::{Communicator, Equivalence, Root}; -use rlst_common::types::Scalar; use rlst_dense::traits::Shape; use rlst_dense::traits::{RawAccess, RawAccessMut}; +use rlst_dense::types::RlstScalar; -pub struct MpiCsrMatrix<'a, T: Scalar + Equivalence, C: Communicator> { +pub struct MpiCsrMatrix<'a, T: RlstScalar + Equivalence, C: Communicator> { mat_type: SparseMatType, shape: [usize; 2], local_matrix: CsrMatrix, @@ -26,7 +26,7 @@ pub struct MpiCsrMatrix<'a, T: Scalar + Equivalence, C: Communicator> { domain_ghosts: crate::ghost_communicator::GhostCommunicator, } -impl<'a, T: Scalar + Equivalence, C: Communicator> MpiCsrMatrix<'a, T, C> { +impl<'a, T: RlstScalar + Equivalence, C: Communicator> MpiCsrMatrix<'a, T, C> { pub fn new( shape: [usize; 2], indices: Vec, @@ -283,7 +283,7 @@ impl<'a, T: Scalar + Equivalence, C: Communicator> MpiCsrMatrix<'a, T, C> { } } -impl<'a, T: Scalar + Equivalence, C: Communicator> Shape<2> for MpiCsrMatrix<'a, T, C> { +impl<'a, T: RlstScalar + Equivalence, C: Communicator> Shape<2> for MpiCsrMatrix<'a, T, C> { fn shape(&self) -> [usize; 2] { self.shape } diff --git a/sparse/src/sparse/tools.rs b/sparse/src/sparse/tools.rs index 559929f8..b98b90e1 100644 --- a/sparse/src/sparse/tools.rs +++ b/sparse/src/sparse/tools.rs @@ -1,6 +1,6 @@ //! Tools for sparse matrix handling use crate::sparse::SparseMatType; -use rlst_common::types::Scalar; +use rlst_dense::types::RlstScalar; /// Normalize an Aij matrix. /// @@ -9,7 +9,7 @@ use rlst_common::types::Scalar; /// deleted and the entries sorted row-wise with increasing column indicies /// (if `sort_mode` is [SparseMatType::Csr]) or column-wise with increasing /// row indices (if `sort_mode` is [SparseMatType::Csc]). -pub fn normalize_aij( +pub fn normalize_aij( rows: &[usize], cols: &[usize], data: &[T], diff --git a/sparse/src/sparse/umfpack.rs b/sparse/src/sparse/umfpack.rs index 2525b3a4..8fccc51a 100644 --- a/sparse/src/sparse/umfpack.rs +++ b/sparse/src/sparse/umfpack.rs @@ -3,11 +3,11 @@ use std::ffi::c_void; use crate::sparse::csc_mat::CscMatrix; -use rlst_common::types::Scalar; -use rlst_common::types::*; -use rlst_common::types::{RlstError, RlstResult}; use rlst_dense::array::Array; use rlst_dense::traits::{RawAccess, RawAccessMut, Shape, Stride}; +use rlst_dense::types::RlstScalar; +use rlst_dense::types::*; +use rlst_dense::types::{RlstError, RlstResult}; use rlst_umfpack as umfpack; use super::csr_mat::CsrMatrix; @@ -18,7 +18,7 @@ pub enum TransposeMode { ConjugateTrans, } -pub struct UmfpackLu { +pub struct UmfpackLu { shape: [usize; 2], indices: Vec, indptr: Vec, @@ -26,7 +26,7 @@ pub struct UmfpackLu { numeric: *mut c_void, } -impl Drop for UmfpackLu { +impl Drop for UmfpackLu { fn drop(&mut self) { unsafe { umfpack::umfpack_di_free_numeric(&mut self.numeric) }; } @@ -47,7 +47,7 @@ impl UmfpackLu { rhs: Array, mut x: Array, trans: TransposeMode, - ) -> rlst_common::types::RlstResult<()> { + ) -> rlst_dense::types::RlstResult<()> { assert_eq!(rhs.stride()[0], 1); assert_eq!(x.stride()[0], 1); @@ -98,7 +98,7 @@ impl UmfpackLu { rhs: Array, mut x: Array, trans: TransposeMode, - ) -> rlst_common::types::RlstResult<()> { + ) -> rlst_dense::types::RlstResult<()> { assert_eq!(rhs.stride()[0], 1); assert_eq!(x.stride()[0], 1); diff --git a/sparse/src/traits/index_layout.rs b/sparse/src/traits/index_layout.rs index c6d29e39..b4d14a12 100644 --- a/sparse/src/traits/index_layout.rs +++ b/sparse/src/traits/index_layout.rs @@ -3,7 +3,7 @@ //! An [IndexLayout] specified how degrees of freedom are distributed among processes. //! We always assume that a process has a contiguous set of degrees of freedom. -use rlst_common::types::RlstResult; +use rlst_dense::types::RlstResult; /// The Index Layout trait. It fully specifies how degrees of freedom are distributed /// among processes. Each process must hold a contiguous number of degrees of freedom (dofs). diff --git a/sparse/src/traits/indexable_vector.rs b/sparse/src/traits/indexable_vector.rs index 224777a3..092ea567 100644 --- a/sparse/src/traits/indexable_vector.rs +++ b/sparse/src/traits/indexable_vector.rs @@ -1,10 +1,10 @@ //! An indexable vector is the standard type for n-dimensional containers use crate::traits::index_layout::IndexLayout; -use rlst_common::types::Scalar; +use rlst_dense::types::RlstScalar; pub trait IndexableVector { - type Item: Scalar; + type Item: RlstScalar; type Ind: IndexLayout; type View<'a>: IndexableVectorView @@ -21,7 +21,7 @@ pub trait IndexableVector { } pub trait IndexableVectorView { - type T: Scalar; + type T: RlstScalar; type Iter<'a>: std::iter::Iterator where Self: 'a; @@ -61,37 +61,37 @@ pub trait IndexableVectorViewMut: IndexableVectorView { /// Inner product with another object. pub trait Inner: IndexableVector { - fn inner(&self, other: &Self) -> rlst_common::types::RlstResult; + fn inner(&self, other: &Self) -> rlst_dense::types::RlstResult; } /// Take the sum of the squares of the absolute values of the entries. pub trait AbsSquareSum: IndexableVector { - fn abs_square_sum(&self) -> ::Real; + fn abs_square_sum(&self) -> ::Real; } /// Return the 1-Norm (Sum of absolute values of the entries). pub trait Norm1: IndexableVector { - fn norm_1(&self) -> ::Real; + fn norm_1(&self) -> ::Real; } /// Return the 2-Norm (Sqrt of the sum of squares). pub trait Norm2: IndexableVector { - fn norm_2(&self) -> ::Real; + fn norm_2(&self) -> ::Real; } /// Return the supremum norm (largest absolute value of the entries). pub trait NormInfty: IndexableVector { - fn norm_infty(&self) -> ::Real; + fn norm_infty(&self) -> ::Real; } /// Swap entries with another vector. pub trait Swap: IndexableVector { - fn swap(&mut self, other: &mut Self) -> rlst_common::types::RlstResult<()>; + fn swap(&mut self, other: &mut Self) -> rlst_dense::types::RlstResult<()>; } /// Fill vector by copying from another vector. pub trait Fill: IndexableVector { - fn fill(&mut self, other: &Self) -> rlst_common::types::RlstResult<()>; + fn fill(&mut self, other: &Self) -> rlst_dense::types::RlstResult<()>; } /// Multiply entries with a scalar. @@ -105,5 +105,5 @@ pub trait MultSumInto: IndexableVector { &mut self, other: &Self, scalar: Self::Item, - ) -> rlst_common::types::RlstResult<()>; + ) -> rlst_dense::types::RlstResult<()>; }