Simple Model using ModelingToolkit and Clapeyron failing

I’ve been working on implementing some simple models utilizing both ModelingToolkit and Clapeyron. I put together this simple model of a generic compressor:

using ModelingToolkit, NonlinearSolve
using Clapeyron

model = PR(["oxygen", "nitrogen"])
x = [0.21, 0.79]
efficiency = 0.8

@variables T2, T2s
@parameters T1, P1 ,P2

eqs = [
    entropy(model, P1, T1, x) ~ entropy(model, P2, T2s, x),
    efficiency ~ (enthalpy(model, P2, T2s, x) - enthalpy(model, P1, T1, x))/
        (enthalpy(model, P2, T2, x) - enthalpy(model, P1, T1, x))
]

@mtkbuild ns = NonlinearSystem(eqs)

guesses = [T2 => T1 + 10.0, T2s => T1 + 10.0]
ps = [T1 => 293.0, P1 => 101325.0, P2 => 300000.0]

prob = NonlinearProblem(ns, guesses, ps)
sol = solve(prob, NewtonRaphson())

When I run this, I get the following error:

UndefVarError: `unknown` not defined

Stacktrace:
[1] macro expansion
  @ /opt/julia/packages/SymbolicUtils/cLLXc/src/code.jl:418 [inlined]
... (omitted a bit of the stactrace for clarity)
[20] solve(prob::NonlinearProblem{Vector{Float64}, true, ModelingToolkit.MTKParameters{Vector{Float64}, StaticArraysCore.SizedVector{0, Any, Vector{Any}}, Tuple{}, Tuple{}, Tuple{}, RuntimeGeneratedFunctions.RuntimeGeneratedFunction{(:ˍ₋out, :ˍ₋arg1), ModelingToolkit.var"#_RGF_ModTag", ModelingToolkit.var"#_RGF_ModTag", (0x8affda93, 0x1a5e37a7, 0xa8ce5e86, 0x0678db27, 0x99eb41c6), Nothing}, RuntimeGeneratedFunctions.RuntimeGeneratedFunction{(:ˍ₋arg1,), ModelingToolkit.var"#_RGF_ModTag", ModelingToolkit.var"#_RGF_ModTag", (0xdf1028c7, 0x42787cff, 0xadff9b07, 0x887df60f, 0x719e4156), Nothing}}, NonlinearFunction{true, SciMLBase.FullSpecialize, ModelingToolkit.var"#f#614"{RuntimeGeneratedFunctions.RuntimeGeneratedFunction{(:ˍ₋arg1, :ˍ₋arg2), ModelingToolkit.var"#_RGF_ModTag", ModelingToolkit.var"#_RGF_ModTag", (0xabe9932d, 0xd9194a65, 0x43ba96d3, 0x44399e97, 0xc7055b4f), Nothing}, RuntimeGeneratedFunctions.RuntimeGeneratedFunction{(:ˍ₋out, :ˍ₋arg1, :ˍ₋arg2), ModelingToolkit.var"#_RGF_ModTag", ModelingToolkit.var"#_RGF_ModTag", (0x313eaf66, 0xa195af1c, 0x1c98f348, 0x3bc97140, 0x6af6266b), Nothing}}, LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, ModelingToolkit.ObservedFunctionCache{NonlinearSystem}, Nothing, NonlinearSystem, Nothing}, @Kwargs{}, SciMLBase.StandardNonlinearProblem}, args::GeneralizedFirstOrderAlgorithm{nothing, :NewtonRaphson, NoLineSearch, Missing, NewtonDescent{Nothing, typeof(NonlinearSolve.DEFAULT_PRECS)}, Nothing, Nothing, Nothing})
    @ DiffEqBase /opt/julia/packages/DiffEqBase/DdIeW/src/solve.jl:1050
 [21] top-level scope
    @ In[2]:23

I’m newer to Julia, and the error and stacktrace is not that helpful. If I vectorize these equations and solve with NLsolve, it works fine. I was just hoping to utilize the organization and clarity that ModelingToolkit adds to the mix.

Any help is appreciated.

1 Like

Perhaps this issue is related: Integration with Symbolics.jl · Issue #186 · ClapeyronThermo/Clapeyron.jl · GitHub

If this would be the problem, then there should also be work-arounds using @register_symbolic.

This looks already wrong:

Did you try:

@variables T2 T2s
@parameters T1 P1 P2

???

I did try to remove the commas, and same error. Thanks for the suggestion.

I think it’s something related to that thread you linked to. I’m passing a model and mol fraction vector into the equation set, and not registering them as parameters or variables. In that link, they pass the model as I do, but separate out the mol fraction vector as a parameter in a way that won’t work for me.

Neither way seems to work for me.

  • @variables T2, T2s

Both of

@variables T2, T2s

and

@variables T2 T2s

are valid constructs, unless something has changed. It is also possible to put a comma separated list in parenthesis (tuple like), and break the list over two lines, e.g.,

@variables (T2,
          T2s)

Finally, it is possible to include the variables in a begin...end statement:

@variables begin
   T2
   T2s
end

I guess the choice depends on how vertically compact one wants the code to be. Using the begin ... end construct gives more space for adding various other information, such as textual description, units, etc.

I wanted to reply here, in case anyone searches this in the future. The issue arose from inside Clapeyron, and you need to explicitly specify the phase to get ModelingToolkit to work with it.

The code that works for this is:

using ModelingToolkit, NonlinearSolve
using Clapeyron

model = PR(["air", "water"])
x = [1., 0.]

@variables T2 T2s
@parameters T1 P1 P2 eff

eqs = [
    entropy(model, P1, T1, x, phase="unknown") ~ entropy(model, P2, T2s, x, phase="unknown"),
    eff ~ (enthalpy(model, P2, T2s, x, phase="unknown") - enthalpy(model, P1, T1, x, phase="unknown"))/
        (enthalpy(model, P2, T2, x, phase="unknown") - enthalpy(model, P1, T1, x, phase="unknown"))
]

@mtkbuild ns = NonlinearSystem(eqs)

guesses = [T2 => T1 + 10.0, T2s => T1 + 10.0]
ps = [T1 => 293.0, P1 => 101325.0, P2 => 300000.0, eff => 0.8]

prob = NonlinearProblem(ns, guesses, ps)
sol = solve(prob, NewtonRaphson())
3 Likes