Ensure the output of addition is sorted.

src/sparse/cs_matrix_ops.rs

@@ -241,7 +241,11 @@ where
                 &mut res,
             );
 
-            for p in res.data.p[j]..nz {
+            // Keep the output sorted.
+            let range = res.data.p[j]..nz;
+            res.data.i[range.clone()].sort();
+
+            for p in range {
                 res.data.vals[p] = workspace[res.data.i[p]]
             }
         }

tests/sparse/cs_ops.rs

@@ -47,8 +47,8 @@ fn cs_mat_mul() {
 #[test]
 fn cs_mat_add() {
     let m1 = Matrix4x5::new(
-        4.0, 1.0, 4.0, 0.0, 9.0,
-        5.0, 6.0, 0.0, 8.0, 10.0,
+        4.0, 1.0, 4.0, 0.0, 0.0,
+        5.0, 6.0, 0.0, 8.0, 0.0,
         9.0, 10.0, 11.0, 12.0, 0.0,
         0.0, 0.0, 1.0, 0.0, 10.0
     );