diff -r 1a38447a89fa -r 9226980e45a7 src/nalgebra_support.rs
--- a/src/nalgebra_support.rs	Sat Dec 14 09:31:27 2024 -0500
+++ b/src/nalgebra_support.rs	Tue Dec 31 08:30:02 2024 -0500
@@ -23,51 +23,49 @@
 use num_traits::identities::{Zero, One};
 use crate::linops::*;
 use crate::euclidean::*;
+use crate::mapping::{Space, BasicDecomposition};
 use crate::types::Float;
 use crate::norms::*;
+use crate::instance::Instance;
 
-impl<SM,SV,N,M,K,E> Apply<Matrix<E,M,K,SV>> for Matrix<E,N,M,SM>
-where SM: Storage<E,N,M>, SV: Storage<E,M,K>,
-        N : Dim, M : Dim, K : Dim, E : Scalar + Zero + One + ClosedAddAssign + ClosedMulAssign,
-        DefaultAllocator : Allocator<N,K>,
-        DefaultAllocator : Allocator<M,K>,
-        DefaultAllocator : Allocator<N,M>,
-        DefaultAllocator : Allocator<M,N> {
-    type Output = OMatrix<E,N,K>;
-
-    #[inline]
-    fn apply(&self, x : Matrix<E,M,K,SV>) -> Self::Output {
-        self.mul(x)
-    }
+impl<SM,N,M,E> Space for Matrix<E,N,M,SM>
+where
+    SM: Storage<E,N,M> + Clone,
+    N : Dim, M : Dim, E : Scalar + Zero + One + ClosedAddAssign + ClosedMulAssign,
+    DefaultAllocator : Allocator<N,M>,
+{
+    type Decomp = BasicDecomposition;
 }
 
-impl<'a, SM,SV,N,M,K,E> Apply<&'a Matrix<E,M,K,SV>> for Matrix<E,N,M,SM>
-where SM: Storage<E,N,M>, SV: Storage<E,M,K>,
-        N : Dim, M : Dim, K : Dim, E : Scalar + Zero + One + ClosedAddAssign + ClosedMulAssign,
-        DefaultAllocator : Allocator<N,K>,
-        DefaultAllocator : Allocator<M,K>,
-        DefaultAllocator : Allocator<N,M>,
-        DefaultAllocator : Allocator<M,N> {
-    type Output = OMatrix<E,N,K>;
-
-    #[inline]
-    fn apply(&self, x : &'a Matrix<E,M,K,SV>) -> Self::Output {
-        self.mul(x)
-    }
-}
-
-impl<'a, SM,SV,N,M,K,E> Linear<Matrix<E,M,K,SV>> for Matrix<E,N,M,SM>
-where SM: Storage<E,N,M>, SV: Storage<E,M,K>,
+impl<SM,SV,N,M,K,E> Mapping<Matrix<E,M,K,SV>> for Matrix<E,N,M,SM>
+where SM: Storage<E,N,M>, SV: Storage<E,M,K> + Clone,
         N : Dim, M : Dim, K : Dim, E : Scalar + Zero + One + ClosedAddAssign + ClosedMulAssign,
         DefaultAllocator : Allocator<N,K>,
         DefaultAllocator : Allocator<M,K>,
         DefaultAllocator : Allocator<N,M>,
         DefaultAllocator : Allocator<M,N> {
     type Codomain = OMatrix<E,N,K>;
+
+    #[inline]
+    fn apply<I : Instance<Matrix<E,M,K,SV>>>(
+        &self, x : I
+    ) -> Self::Codomain {
+        x.either(|owned| self.mul(owned), |refr| self.mul(refr))
+    }
+}
+
+
+impl<'a, SM,SV,N,M,K,E> Linear<Matrix<E,M,K,SV>> for Matrix<E,N,M,SM>
+where SM: Storage<E,N,M>, SV: Storage<E,M,K> + Clone,
+        N : Dim, M : Dim, K : Dim, E : Scalar + Zero + One + ClosedAddAssign + ClosedMulAssign,
+        DefaultAllocator : Allocator<N,K>,
+        DefaultAllocator : Allocator<M,K>,
+        DefaultAllocator : Allocator<N,M>,
+        DefaultAllocator : Allocator<M,N> {
 }
 
 impl<SM,SV1,SV2,N,M,K,E> GEMV<E, Matrix<E,M,K,SV1>, Matrix<E,N,K,SV2>> for Matrix<E,N,M,SM>
-where SM: Storage<E,N,M>, SV1: Storage<E,M,K>, SV2: StorageMut<E,N,K>,
+where SM: Storage<E,N,M>, SV1: Storage<E,M,K> + Clone, SV2: StorageMut<E,N,K>,
       N : Dim, M : Dim, K : Dim, E : Scalar + Zero + One + Float,
       DefaultAllocator : Allocator<N,K>,
       DefaultAllocator : Allocator<M,K>,
@@ -75,34 +73,36 @@
       DefaultAllocator : Allocator<M,N> {
 
     #[inline]
-    fn gemv(&self, y : &mut Matrix<E,N,K,SV2>, α : E, x : &Matrix<E,M,K,SV1>, β : E) {
-        Matrix::gemm(y, α, self, x, β)
+    fn gemv<I : Instance<Matrix<E,M,K,SV1>>>(
+        &self, y : &mut Matrix<E,N,K,SV2>, α : E, x : I, β : E
+    ) {
+        x.eval(|x̃| Matrix::gemm(y, α, self, x̃, β))
     }
 
     #[inline]
-    fn apply_mut<'a>(&self, y : &mut Matrix<E,N,K,SV2>, x : &Matrix<E,M,K,SV1>) {
-        self.mul_to(x, y)
+    fn apply_mut<'a, I : Instance<Matrix<E,M,K,SV1>>>(&self, y : &mut Matrix<E,N,K,SV2>, x : I) {
+        x.eval(|x̃| self.mul_to(x̃, y))
     }
 }
 
 impl<SM,SV1,M,E> AXPY<E, Vector<E,M,SV1>> for Vector<E,M,SM>
-where SM: StorageMut<E,M>, SV1: Storage<E,M>,
+where SM: StorageMut<E,M> + Clone, SV1: Storage<E,M> + Clone,
       M : Dim, E : Scalar + Zero + One + Float,
       DefaultAllocator : Allocator<M> {
 
     #[inline]
-    fn axpy(&mut self, α : E, x : &Vector<E,M,SV1>, β : E) {
-        Matrix::axpy(self, α, x, β)
+    fn axpy<I : Instance<Vector<E,M,SV1>>>(&mut self, α : E, x : I, β : E) {
+        x.eval(|x̃| Matrix::axpy(self, α, x̃, β))
     }
 
     #[inline]
-    fn copy_from(&mut self, y : &Vector<E,M,SV1>) {
-        Matrix::copy_from(self, y)
+    fn copy_from<I : Instance<Vector<E,M,SV1>>>(&mut self, y : I) {
+        y.eval(|ỹ| Matrix::copy_from(self, ỹ))
     }
 }
 
 impl<SM,M,E> Projection<E, Linfinity> for Vector<E,M,SM>
-where SM: StorageMut<E,M>,
+where SM: StorageMut<E,M> + Clone,
       M : Dim, E : Scalar + Zero + One + Float + RealField,
       DefaultAllocator : Allocator<M> {
     #[inline]
@@ -111,9 +111,9 @@
     }
 }
 
-impl<'own,SV1,SV2,SM,N,M,K,E> Adjointable<Matrix<E,M,K,SV1>,Matrix<E,N,K,SV2>>
+impl<'own,SV1,SV2,SM,N,M,K,E> Adjointable<Matrix<E,M,K,SV1>, Matrix<E,N,K,SV2>>
 for Matrix<E,N,M,SM>
-where SM: Storage<E,N,M>, SV1: Storage<E,M,K>, SV2: Storage<E,N,K>,
+where SM: Storage<E,N,M>, SV1: Storage<E,M,K> + Clone, SV2: Storage<E,N,K> + Clone,
       N : Dim, M : Dim, K : Dim, E : Scalar + Zero + One + SimdComplexField,
       DefaultAllocator : Allocator<N,K>,
       DefaultAllocator : Allocator<M,K>,
@@ -170,7 +170,7 @@
 impl<E,M,S> Euclidean<E>
 for Vector<E,M,S>
 where M : Dim,
-      S : StorageMut<E,M>,
+      S : StorageMut<E,M> + Clone,
       E : Float + Scalar + Zero + One + RealField,
       DefaultAllocator : Allocator<M> {
 
@@ -195,7 +195,7 @@
 impl<E,M,S> StaticEuclidean<E>
 for Vector<E,M,S>
 where M : DimName,
-      S : StorageMut<E,M>,
+      S : StorageMut<E,M> + Clone,
       E : Float + Scalar + Zero + One + RealField,
       DefaultAllocator : Allocator<M> {