diff --git a/NDTensors/src/diag/tensoralgebra/outer.jl b/NDTensors/src/diag/tensoralgebra/outer.jl
index 7e6c5ef6b2..4170a98bcb 100644
--- a/NDTensors/src/diag/tensoralgebra/outer.jl
+++ b/NDTensors/src/diag/tensoralgebra/outer.jl
@@ -1,10 +1,9 @@
 function outer!(
         R::DenseTensor{<:Number, NR}, T1::DiagTensor{<:Number, N1}, T2::DiagTensor{<:Number, N2}
     ) where {NR, N1, N2}
-    for i1 in 1:diaglength(T1), i2 in 1:diaglength(T2)
-        indsR = CartesianIndex{NR}(ntuple(r -> r ≤ N1 ? i1 : i2, Val(NR)))
-        R[indsR] = getdiagindex(T1, i1) * getdiagindex(T2, i2)
-    end
+    t1 = T1.storage.data
+    t2 = T2.storage.data
+    array(R) .= t1 .* t2'
     return R
 end
 
diff --git a/NDTensors/test/test_diag.jl b/NDTensors/test/test_diag.jl
index 50b922dc6d..ac54e023fe 100644
--- a/NDTensors/test/test_diag.jl
+++ b/NDTensors/test/test_diag.jl
@@ -123,6 +123,11 @@ end
 
     ## Test dot on GPU
     @test dot(t, A) ≈ dot(dev(array(t)), array(A)) rtol = sqrt(eps(elt))
+    
+    NDTensors.outer!(A, t,t);
+    for i in NDTensors.cpu(A)
+       @test i == one(elt)
+    end
 end
 nothing
 end