diff --git a/stim/math/matrix.h b/stim/math/matrix.h
index f4cbe42..513b588 100644
--- a/stim/math/matrix.h
+++ b/stim/math/matrix.h
@@ -64,6 +64,16 @@ struct matrix
 
 		return *this;
 	}
+	
+	// M - rhs*I
+	CUDA_CALLABLE matrix<T, N> operator-(T rhs)
+	{
+		for(int i=0; i<N; i++)
+			for(int j=0 ; j<N; j++)
+				if(i == j)
+					M[i*N+j] -= rhs;
+		return *this;
+	}
 
 	template<typename Y>
 	vec<Y> operator*(vec<Y> rhs){
diff --git a/stim/math/tensor3.h b/stim/math/tensor3.h
index f0ed778..b18a67b 100644
--- a/stim/math/tensor3.h
+++ b/stim/math/tensor3.h
@@ -80,11 +80,11 @@ namespace stim {
 			stim::matrix<T, 3> M2 = matrix_sym<T, 3>::mat();
 			stim::matrix<T, 3> M3 = matrix_sym<T, 3>::mat();	// fill a tensor with symmetric values
 
-			M1.operator_minus(M1, lambda[0]);					// M1 = A - lambda[0] * I
+			M1 = M1 - lambda[0];					// M1 = A - lambda[0] * I
 	
-			M2.operator_minus(M2, lambda[1]);					// M2 = A - lambda[1] * I
+			M2 = M2 - lambda[1];					// M2 = A - lambda[1] * I
 
-			M3.operator_minus(M3, lambda[2]);					// M3 = A - lambda[2] * I
+			M3 = M3 - lambda[2];					// M3 = A - lambda[2] * I
 
 			T Mod = 0;											// module of one column
 
--
libgit2 0.21.4