diff -r edb95d2b83cc -r d2acaaddd9af src/nalgebra_support.rs
--- a/src/nalgebra_support.rs	Sun Nov 10 09:02:57 2024 -0500
+++ b/src/nalgebra_support.rs	Tue Dec 31 09:12:43 2024 -0500
@@ -56,9 +56,9 @@
     }
 }
 
-impl<'a, SM,SV,N,M,K,E> Linear<Matrix<E,M,K,SV>> for Matrix<E,N,M,SM>
+impl<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>,
-        N : Dim, M : Dim, K : Dim, E : Scalar + ClosedMul + ClosedAdd + Zero + One,
+        N : Dim, M : Dim, K : Dim, E : Float + Scalar + ClosedMul + ClosedAdd + Zero + One,
         DefaultAllocator : Allocator<E,N,K>,
         DefaultAllocator : Allocator<E,M,K>,
         DefaultAllocator : Allocator<E,N,M>,
@@ -75,12 +75,31 @@
       DefaultAllocator : Allocator<E,M,N> {
 
     #[inline]
-    fn gemv(&self, y : &mut Matrix<E,N,K,SV2>, α : E, x : &Matrix<E,M,K,SV1>, β : E) {
+    fn gemv(&self, y : &mut Matrix<E,N,K,SV2>, α : E, x : Matrix<E,M,K,SV1>, β : E) {
+        Matrix::gemm(y, α, self, &x, β)
+    }
+
+    #[inline]
+    fn apply_mut(&self, y : &mut Matrix<E,N,K,SV2>, x : Matrix<E,M,K,SV1>) {
+        self.mul_to(&x, y)
+    }
+}
+
+impl<'a, SM,SV1,SV2,N,M,K,E> GEMV<E, &'a 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>,
+      N : Dim, M : Dim, K : Dim, E : Scalar + ClosedMul + ClosedAdd + Zero + One + Float,
+      DefaultAllocator : Allocator<E,N,K>,
+      DefaultAllocator : Allocator<E,M,K>,
+      DefaultAllocator : Allocator<E,N,M>,
+      DefaultAllocator : Allocator<E,M,N> {
+
+    #[inline]
+    fn gemv(&self, y : &mut Matrix<E,N,K,SV2>, α : E, x : &'a Matrix<E,M,K,SV1>, β : E) {
         Matrix::gemm(y, α, self, x, β)
     }
 
     #[inline]
-    fn apply_mut<'a>(&self, y : &mut Matrix<E,N,K,SV2>, x : &Matrix<E,M,K,SV1>) {
+    fn apply_mut(&self, y : &mut Matrix<E,N,K,SV2>, x : &'a Matrix<E,M,K,SV1>) {
         self.mul_to(x, y)
     }
 }
@@ -91,17 +110,33 @@
       DefaultAllocator : Allocator<E,M> {
 
     #[inline]
-    fn axpy(&mut self, α : E, x : &Vector<E,M,SV1>, β : E) {
+    fn axpy(&mut self, α : E, x : Vector<E,M,SV1>, β : E) {
+        Matrix::axpy(self, α, &x, β)
+    }
+
+    #[inline]
+    fn copy_from(&mut self, y : Vector<E,M,SV1>) {
+        Matrix::copy_from(self, &y)
+    }
+}
+
+impl<'a, SM,SV1,M,E> AXPY<E, &'a Vector<E,M,SV1>> for Vector<E,M,SM>
+where SM: StorageMut<E,M>, SV1: Storage<E,M>,
+      M : Dim, E : Scalar + ClosedMul + ClosedAdd + Zero + One + Float,
+      DefaultAllocator : Allocator<E,M> {
+
+    #[inline]
+    fn axpy(&mut self, α : E, x : &'a Vector<E,M,SV1>, β : E) {
         Matrix::axpy(self, α, x, β)
     }
 
     #[inline]
-    fn copy_from(&mut self, y : &Vector<E,M,SV1>) {
+    fn copy_from(&mut self, y : &'a Vector<E,M,SV1>) {
         Matrix::copy_from(self, y)
     }
 }
 
-impl<SM,M,E> Projection<E, Linfinity> for Vector<E,M,SM>
+impl<SM,M,E> Projection<Linfinity> for Vector<E,M,SM>
 where SM: StorageMut<E,M>,
       M : Dim, E : Scalar + ClosedMul + ClosedAdd + Zero + One + Float + RealField,
       DefaultAllocator : Allocator<E,M> {
@@ -114,7 +149,8 @@
 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>,
-      N : Dim, M : Dim, K : Dim, E : Scalar + ClosedMul + ClosedAdd + Zero + One + SimdComplexField,
+      N : Dim, M : Dim, K : Dim,
+      E : Float + Scalar + ClosedMul + ClosedAdd + Zero + One + SimdComplexField,
       DefaultAllocator : Allocator<E,N,K>,
       DefaultAllocator : Allocator<E,M,K>,
       DefaultAllocator : Allocator<E,N,M>,
@@ -128,7 +164,7 @@
     }
 }
 
-impl<E,M,S,Si> Dot<Vector<E,M,Si>,E>
+impl<E,M,S,Si> Dot<Vector<E,M,Si>>
 for Vector<E,M,S>
 where M : Dim,
       E : Float + Scalar + ClosedMul + ClosedAdd + Zero + One,
@@ -137,7 +173,21 @@
       DefaultAllocator : Allocator<E,M> {
 
     #[inline]
-    fn dot(&self, other : &Vector<E,M,Si>) -> E {
+    fn dot(&self, other : Vector<E,M,Si>) -> E {
+        Vector::<E,M,S>::dot(self, &other)
+    }
+}
+
+impl<'a, E,M,S,Si> Dot<&'a Vector<E,M,Si>>
+for Vector<E,M,S>
+where M : Dim,
+      E : Float + Scalar + ClosedMul + ClosedAdd + Zero + One,
+      S : Storage<E,M>,
+      Si : Storage<E,M>,
+      DefaultAllocator : Allocator<E,M> {
+
+    #[inline]
+    fn dot(&self, other : &'a Vector<E,M,Si>) -> E {
         Vector::<E,M,S>::dot(self, other)
     }
 }
@@ -167,7 +217,16 @@
 
 // TODO: should allow different input storages in `Euclidean`.
 
-impl<E,M,S> Euclidean<E>
+impl<E,M,K,S> HasScalarField
+for Matrix<E,M,K,S>
+where M : Dim, K : Dim,
+      S : Storage<E,M,K>,
+      E : Float + Scalar,
+      DefaultAllocator : Allocator<E,M,K> {
+    type Field = E;
+}
+
+impl<E,M,S> Euclidean
 for Vector<E,M,S>
 where M : Dim,
       S : StorageMut<E,M>,
@@ -192,7 +251,7 @@
     }
 }
 
-impl<E,M,S> StaticEuclidean<E>
+impl<E,M,S> StaticEuclidean
 for Vector<E,M,S>
 where M : DimName,
       S : StorageMut<E,M>,
@@ -205,7 +264,7 @@
     }
 }
 
-impl<E,M,S> Norm<E, L1>
+impl<E,M,S> Norm<L1>
 for Vector<E,M,S>
 where M : Dim,
       S : StorageMut<E,M>,
@@ -218,7 +277,7 @@
     }
 }
 
-impl<E,M,S> Dist<E, L1>
+impl<E,M,S> Dist<L1>
 for Vector<E,M,S>
 where M : Dim,
       S : StorageMut<E,M>,
@@ -230,7 +289,7 @@
     }
 }
 
-impl<E,M,S> Norm<E, L2>
+impl<E,M,S> Norm<L2>
 for Vector<E,M,S>
 where M : Dim,
       S : StorageMut<E,M>,
@@ -243,7 +302,7 @@
     }
 }
 
-impl<E,M,S> Dist<E, L2>
+impl<E,M,S> Dist<L2>
 for Vector<E,M,S>
 where M : Dim,
       S : StorageMut<E,M>,
@@ -255,7 +314,7 @@
     }
 }
 
-impl<E,M,S> Norm<E, Linfinity>
+impl<E,M,S> Norm<Linfinity>
 for Vector<E,M,S>
 where M : Dim,
       S : StorageMut<E,M>,
@@ -268,7 +327,7 @@
     }
 }
 
-impl<E,M,S> Dist<E, Linfinity>
+impl<E,M,S> Dist<Linfinity>
 for Vector<E,M,S>
 where M : Dim,
       S : StorageMut<E,M>,