Merge pull request #395 from SciML/myb/wfact_fix

Fix Wfact(_t) calculation
SciML · May 25, 2020 · ee79407 · ee79407 · YingboMa · May 25, 2020
2 parents 470610e + 910b6a4
commit ee79407
Show file tree

Hide file tree

Showing 3 changed files with 5 additions and 5 deletions.
diff --git a/Project.toml b/Project.toml
@@ -1,7 +1,7 @@
 name = "ModelingToolkit"
 uuid = "961ee093-0014-501f-94e3-6117800e7a78"
 authors = ["Chris Rackauckas <[email protected]>"]
-version = "3.6.3"
+version = "3.6.4"
 
 [deps]
 ArrayInterface = "4fba245c-0d91-5ea0-9b3e-6abc04ee57a9"

diff --git a/src/systems/diffeqs/abstractodesystem.jl b/src/systems/diffeqs/abstractodesystem.jl
@@ -63,16 +63,17 @@ function calculate_factorized_W(sys::AbstractODESystem, simplify=true)
     isempty(sys.Wfact[]) || return (sys.Wfact[],sys.Wfact_t[])
 
     jac = calculate_jacobian(sys)
+    M = calculate_massmatrix(sys)
     gam = Variable(:__MTKWgamma)()
 
-    W = - LinearAlgebra.I + gam*jac
+    W = - M + gam*jac
     Wfact = lu(W, Val(false), check=false).factors
 
     if simplify
         Wfact = ModelingToolkit.simplify.(Wfact)
     end
 
-    W_t = - LinearAlgebra.I/gam + jac
+    W_t = - M/gam + jac
     Wfact_t = lu(W_t, Val(false), check=false).factors
     if simplify
         Wfact_t = ModelingToolkit.simplify.(Wfact_t)

diff --git a/test/odesystem.jl b/test/odesystem.jl
@@ -212,9 +212,8 @@ p  = [k₁ => 0.04,
 tspan = (0.0,100000.0)
 prob1 = ODEProblem(sys,u0,tspan,p)
 prob2 = ODEProblem(sys,u0,tspan,p,jac=true)
-# Wfact version is not very stable because of the lack of pivoting
 prob3 = ODEProblem(sys,u0,tspan,p,Wfact=true,Wfact_t=true)
-for (prob, atol) in [(prob1, 1e-12), (prob2, 1e-12), (prob3, 0.1)]
+for (prob, atol) in [(prob1, 1e-12), (prob2, 1e-12), (prob3, 1e-12)]
     sol = solve(prob, Rodas5())
     @test all(x->≈(sum(x), 1.0, atol=atol), sol.u)
 end