# Simple Model using ModelingToolkit and Clapeyron failing

**URL:** https://discourse.julialang.org/t/simple-model-using-modelingtoolkit-and-clapeyron-failing/120874
**Category:** Modelling & Simulations
**Tags:** question, package, modelingtoolkit, nonlinearsolve
**Created:** [October 3, 2024, 9:44pm UTC](https://discourse.julialang.org/t/simple-model-using-modelingtoolkit-and-clapeyron-failing/120874 "2024-10-03T21:44:22Z")
**Posts on this page:** 6
**Page:** 1

<div class="post-metadata">

### Author: ![JarodMcCormick](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jarodmccormick/32/212396_2.png) [@JarodMcCormick](https://discourse.julialang.org/u/JarodMcCormick)
#### Post date: [October 3, 2024, 9:44pm UTC](https://discourse.julialang.org/t/simple-model-using-modelingtoolkit-and-clapeyron-failing/120874/1 "2024-10-03T21:44:22Z")

</div>

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

```julia
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:

```julia
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.

---

<div class="post-metadata">

### Author: ![ufechner7](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/ufechner7/32/51363_2.png) [@ufechner7](https://discourse.julialang.org/u/ufechner7)
#### Post date: [October 4, 2024, 12:08am UTC](https://discourse.julialang.org/t/simple-model-using-modelingtoolkit-and-clapeyron-failing/120874/2 "2024-10-04T00:08:06Z")

</div>

Perhaps this issue is related: [Integration with Symbolics.jl · Issue #186 · ClapeyronThermo/Clapeyron.jl · GitHub](https://github.com/ClapeyronThermo/Clapeyron.jl/issues/186)

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

---

<div class="post-metadata">

### Author: ![ufechner7](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/ufechner7/32/51363_2.png) [@ufechner7](https://discourse.julialang.org/u/ufechner7)
#### Post date: [October 4, 2024, 12:17am UTC](https://discourse.julialang.org/t/simple-model-using-modelingtoolkit-and-clapeyron-failing/120874/3 "2024-10-04T00:17:25Z")

</div>

This looks already wrong:

> [@JarodMcCormick](#):
>
> ```julia
> @variables T2, T2s
> @parameters T1, P1 ,P2
> 
> ```

Did you try:

```julia
@variables T2 T2s
@parameters T1 P1 P2

```

???

---

<div class="post-metadata">

### Author: ![JarodMcCormick](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jarodmccormick/32/212396_2.png) [@JarodMcCormick](https://discourse.julialang.org/u/JarodMcCormick)
#### Post date: [October 4, 2024, 3:05am UTC](https://discourse.julialang.org/t/simple-model-using-modelingtoolkit-and-clapeyron-failing/120874/4 "2024-10-04T03:05:03Z")

</div>

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.

---

<div class="post-metadata">

### Author: ![BLI](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/bli/32/37206_2.png) [@BLI](https://discourse.julialang.org/u/BLI)
#### Post date: [October 4, 2024, 7:00am UTC](https://discourse.julialang.org/t/simple-model-using-modelingtoolkit-and-clapeyron-failing/120874/5 "2024-10-04T07:00:13Z")

</div>

- `@variables T2, T2s`

Both of

```julia
@variables T2, T2s

```

and

```julia
@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.,

```julia
@variables (T2,
          T2s)

```

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

```julia
@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.

---

<div class="post-metadata">

### Author: ![JarodMcCormick](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jarodmccormick/32/212396_2.png) [@JarodMcCormick](https://discourse.julialang.org/u/JarodMcCormick)
#### Post date: [October 7, 2024, 11:57pm UTC](https://discourse.julialang.org/t/simple-model-using-modelingtoolkit-and-clapeyron-failing/120874/6 "2024-10-07T23:57:01Z")

</div>

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:

```julia
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())

```
