fix: handle I as mass matrix for CheckInit

Author: AayushSabharwal

src/initialization.jl

@@ -123,6 +123,7 @@ function get_initial_values(
     t = current_time(integrator)
     M = f.mass_matrix
 
+    M == I && return u0, p, true
     algebraic_vars = [all(iszero, x) for x in eachcol(M)]
     algebraic_eqs = [all(iszero, x) for x in eachrow(M)]
     (iszero(algebraic_vars) || iszero(algebraic_eqs)) && return u0, p, true

test/initialization.jl

@@ -1,4 +1,4 @@
-using StochasticDiffEq, OrdinaryDiffEq, NonlinearSolve, SymbolicIndexingInterface, Test
+using StochasticDiffEq, OrdinaryDiffEq, NonlinearSolve, SymbolicIndexingInterface, LinearAlgebra, Test
 
 @testset "CheckInit" begin
     @testset "ODEProblem" begin
@@ -27,6 +27,20 @@ using StochasticDiffEq, OrdinaryDiffEq, NonlinearSolve, SymbolicIndexingInterfac
                 prob, integ, f, SciMLBase.CheckInit(),
                 Val(SciMLBase.isinplace(f)); abstol = 1e-10)
         end
+
+        @testset "With I mass matrix" begin
+            function rhs(u, p, t)
+                return u
+            end
+            prob = ODEProblem(ODEFunction(rhs; mass_matrix = I), ones(2), (0.0, 1.0))
+            integ = init(prob)
+            u0, _, success = SciMLBase.get_initial_values(
+                prob, integ, prob.f, SciMLBase.CheckInit(),
+                Val(false); abstol = 1e-10
+            )
+            @test success
+            @test u0 == prob.u0
+        end
     end
 
     @testset "DAEProblem" begin