diff -r 7b2ee3e84c5f -r 5e3c1874797d src/nalgebra_support.rs
--- a/src/nalgebra_support.rs	Fri Dec 06 16:14:41 2024 -0500
+++ b/src/nalgebra_support.rs	Wed Dec 11 20:45:17 2024 -0500
@@ -10,7 +10,7 @@
 
 use nalgebra::{
     Matrix, Storage, StorageMut, OMatrix, Dim, DefaultAllocator, Scalar,
-    ClosedMul, ClosedAdd, SimdComplexField, Vector, OVector, RealField,
+    ClosedAddAssign, ClosedMulAssign, SimdComplexField, Vector, OVector, RealField,
     LpNorm, UniformNorm
 };
 use nalgebra::Norm as NalgebraNorm;
@@ -28,11 +28,11 @@
 
 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 + ClosedMul + ClosedAdd + Zero + One,
-        DefaultAllocator : Allocator<E,N,K>,
-        DefaultAllocator : Allocator<E,M,K>,
-        DefaultAllocator : Allocator<E,N,M>,
-        DefaultAllocator : Allocator<E,M,N> {
+        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]
@@ -43,11 +43,11 @@
 
 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 + ClosedMul + ClosedAdd + Zero + One,
-        DefaultAllocator : Allocator<E,N,K>,
-        DefaultAllocator : Allocator<E,M,K>,
-        DefaultAllocator : Allocator<E,N,M>,
-        DefaultAllocator : Allocator<E,M,N> {
+        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]
@@ -58,21 +58,21 @@
 
 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>,
-        N : Dim, M : Dim, K : Dim, E : Scalar + ClosedMul + ClosedAdd + Zero + One,
-        DefaultAllocator : Allocator<E,N,K>,
-        DefaultAllocator : Allocator<E,M,K>,
-        DefaultAllocator : Allocator<E,N,M>,
-        DefaultAllocator : Allocator<E,M,N> {
+        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>;
 }
 
 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>,
-      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> {
+      N : Dim, M : Dim, K : Dim, E : Scalar + Zero + One + Float,
+      DefaultAllocator : Allocator<N,K>,
+      DefaultAllocator : Allocator<M,K>,
+      DefaultAllocator : Allocator<N,M>,
+      DefaultAllocator : Allocator<M,N> {
 
     #[inline]
     fn gemv(&self, y : &mut Matrix<E,N,K,SV2>, α : E, x : &Matrix<E,M,K,SV1>, β : E) {
@@ -87,8 +87,8 @@
 
 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>,
-      M : Dim, E : Scalar + ClosedMul + ClosedAdd + Zero + One + Float,
-      DefaultAllocator : Allocator<E,M> {
+      M : Dim, E : Scalar + Zero + One + Float,
+      DefaultAllocator : Allocator<M> {
 
     #[inline]
     fn axpy(&mut self, α : E, x : &Vector<E,M,SV1>, β : E) {
@@ -103,8 +103,8 @@
 
 impl<SM,M,E> Projection<E, 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> {
+      M : Dim, E : Scalar + Zero + One + Float + RealField,
+      DefaultAllocator : Allocator<M> {
     #[inline]
     fn proj_ball_mut(&mut self, ρ : E, _ : Linfinity) {
         self.iter_mut().for_each(|v| *v = num_traits::clamp(*v, -ρ, ρ))
@@ -114,11 +114,11 @@
 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,
-      DefaultAllocator : Allocator<E,N,K>,
-      DefaultAllocator : Allocator<E,M,K>,
-      DefaultAllocator : Allocator<E,N,M>,
-      DefaultAllocator : Allocator<E,M,N> {
+      N : Dim, M : Dim, K : Dim, E : Scalar + Zero + One + SimdComplexField,
+      DefaultAllocator : Allocator<N,K>,
+      DefaultAllocator : Allocator<M,K>,
+      DefaultAllocator : Allocator<N,M>,
+      DefaultAllocator : Allocator<M,N> {
     type AdjointCodomain = OMatrix<E,M,K>;
     type Adjoint<'a> = OMatrix<E,M,N> where SM : 'a;
 
@@ -131,10 +131,10 @@
 impl<E,M,S,Si> Dot<Vector<E,M,Si>,E>
 for Vector<E,M,S>
 where M : Dim,
-      E : Float + Scalar + ClosedMul + ClosedAdd + Zero + One,
+      E : Float + Scalar + Zero + One,
       S : Storage<E,M>,
       Si : Storage<E,M>,
-      DefaultAllocator : Allocator<E,M> {
+      DefaultAllocator : Allocator<M> {
 
     #[inline]
     fn dot(&self, other : &Vector<E,M,Si>) -> E {
@@ -171,8 +171,8 @@
 for Vector<E,M,S>
 where M : Dim,
       S : StorageMut<E,M>,
-      E : Float + Scalar + ClosedMul + ClosedAdd + Zero + One + RealField,
-      DefaultAllocator : Allocator<E,M> {
+      E : Float + Scalar + Zero + One + RealField,
+      DefaultAllocator : Allocator<M> {
 
     type Output = OVector<E, M>;
     
@@ -196,8 +196,8 @@
 for Vector<E,M,S>
 where M : DimName,
       S : StorageMut<E,M>,
-      E : Float + Scalar + ClosedMul + ClosedAdd + Zero + One + RealField,
-      DefaultAllocator : Allocator<E,M> {
+      E : Float + Scalar + Zero + One + RealField,
+      DefaultAllocator : Allocator<M> {
 
     #[inline]
     fn origin() -> OVector<E, M> {
@@ -209,8 +209,8 @@
 for Vector<E,M,S>
 where M : Dim,
       S : StorageMut<E,M>,
-      E : Float + Scalar + ClosedMul + ClosedAdd + Zero + One + RealField,
-      DefaultAllocator : Allocator<E,M> {
+      E : Float + Scalar + Zero + One + RealField,
+      DefaultAllocator : Allocator<M> {
 
     #[inline]
     fn norm(&self, _ : L1) -> E {
@@ -222,8 +222,8 @@
 for Vector<E,M,S>
 where M : Dim,
       S : StorageMut<E,M>,
-      E : Float + Scalar + ClosedMul + ClosedAdd + Zero + One + RealField,
-      DefaultAllocator : Allocator<E,M> {
+      E : Float + Scalar + Zero + One + RealField,
+      DefaultAllocator : Allocator<M> {
     #[inline]
     fn dist(&self, other : &Self, _ : L1) -> E {
         LpNorm(1).metric_distance(self, other)
@@ -234,8 +234,8 @@
 for Vector<E,M,S>
 where M : Dim,
       S : StorageMut<E,M>,
-      E : Float + Scalar + ClosedMul + ClosedAdd + Zero + One + RealField,
-      DefaultAllocator : Allocator<E,M> {
+      E : Float + Scalar + Zero + One + RealField,
+      DefaultAllocator : Allocator<M> {
 
     #[inline]
     fn norm(&self, _ : L2) -> E {
@@ -247,8 +247,8 @@
 for Vector<E,M,S>
 where M : Dim,
       S : StorageMut<E,M>,
-      E : Float + Scalar + ClosedMul + ClosedAdd + Zero + One + RealField,
-      DefaultAllocator : Allocator<E,M> {
+      E : Float + Scalar + Zero + One + RealField,
+      DefaultAllocator : Allocator<M> {
     #[inline]
     fn dist(&self, other : &Self, _ : L2) -> E {
         LpNorm(2).metric_distance(self, other)
@@ -259,8 +259,8 @@
 for Vector<E,M,S>
 where M : Dim,
       S : StorageMut<E,M>,
-      E : Float + Scalar + ClosedMul + ClosedAdd + Zero + One + RealField,
-      DefaultAllocator : Allocator<E,M> {
+      E : Float + Scalar + Zero + One + RealField,
+      DefaultAllocator : Allocator<M> {
 
     #[inline]
     fn norm(&self, _ : Linfinity) -> E {
@@ -272,8 +272,8 @@
 for Vector<E,M,S>
 where M : Dim,
       S : StorageMut<E,M>,
-      E : Float + Scalar + ClosedMul + ClosedAdd + Zero + One + RealField,
-      DefaultAllocator : Allocator<E,M> {
+      E : Float + Scalar + Zero + One + RealField,
+      DefaultAllocator : Allocator<M> {
     #[inline]
     fn dist(&self, other : &Self, _ : Linfinity) -> E {
         UniformNorm.metric_distance(self, other)