Media properties and component encapsulation in MTK

Specific Domains Modelling & Simulations
1

Hi,

Modelica has a Media library that captures the properties of various materials. It implements the equations given in IF97 (a version of was apparently implemented in Julia). It provides e.g. calculation of the Cp, Cv from temperature and pressure.

What’s the best practice to implement something similar in MTK (e.g. as part of a ThermalFluid library)?

In particular:

  • How do you define the “type” of a component? i.e. I’d like to have a Medium “abstract type” with properties Cp, Cv, ‘ρ’ etc. with “concrete” components such as Water and Oil implementing it (sure, I can document what a “Medium” is required to provide, but I’d like MTK to enforce these interfaces)
  • IF97 is rather complex, with dozens of parameters, which are of little interest to most using it. We could implement IF97 fully as part of MTK, but using that model would import all those parameters and equations to the system using it. I think it makes sense for sub-systems to be (optionally) “black-boxes”, without all of their “internals” exposed.

DD

1 Like
2

Hi,

I use CoolProp (low level interace) with MTK.
It requires some interfacing and it has some limitations (it is not thread safe).
I had to map functions with their derivatives and declare them to work with ForwardDiff.
I hope that Clapeyron will be able to replace CoolProp sooner or later. I think it would require less interfacing.

The components of the fluid (and the thermodynamic model) are declared when the MTK component are instantiated.

About your second point, the macro ModelingToolkit.@register allows to keep functions as black-boxes in the equation system.

1 Like
3

Hi, this is an old thread now, but do you have any more information or an example on how you did this? I now need to do the same thing.

4

I have example with CoolProp. I can provide that if you want.

5

Please do! That would be very helpful for me.

6

Would be great if you could share in public. I’m also looking for a solution for this.

7

Hello,
Here is a simple Refrigeration cycle (vapour compression). It is using CoolProp.

using ModelingToolkit, DifferentialEquations, CoolProp


@variables t
D = Differential(t)
const fluid = "SO2"

PropsSI(out::AbstractString, name1::AbstractString, value1::Real, name2::AbstractString, value2::Real, fluid::AbstractString) = CoolProp.PropsSI(out, name1, value1, name2, value2, fluid)
@register_symbolic PropsSI(out::AbstractString, name1::AbstractString, value1::Real, name2::AbstractString, value2::Real, fluid::AbstractString)

@connector  function CoolantPort(;name) 
vars = @variables begin 
    p(t), [input = true]
    h(t)
    mdot(t),[connect =Flow]
end
ODESystem(Equation[], t, vars, [];name=name)
end

function MassSource(;name,source_enthalpy,source_pressure,source_mdot) 
    @named port = CoolantPort()
    para = @parameters begin
        
    end
    vars = @variables begin

     end

    eqs = [
        port.mdot ~ source_mdot 
        port.p ~ source_pressure
        port.h ~ source_enthalpy
    ]
    compose(ODESystem(eqs, t, vars, para;name),port)
end

function compressor_h(πc, h_in, p_in)
    s = PropsSI("S", "H", h_in, "P", p_in, fluid)
    return PropsSI("H", "S", s, "P", πc * p_in, fluid)
end
    @register_symbolic compressor_h(πc, h_in, p_in)


function Compressor(;name,πc) 
    @named inport = CoolantPort()
    @named outport = CoolantPort()
    para = @parameters begin

    end
    vars = @variables begin

     end
   eqs = [   outport.p ~ πc * inport.p
            outport.h ~ compressor_h(πc, inport.h, inport.p)
            outport.mdot ~ inport.mdot
   ]    
   compose(ODESystem(eqs, t, vars, para;name), inport, outport)
end

function Condensor(;name)
    @named inport = CoolantPort()
    @named outport = CoolantPort()
    para = @parameters begin

    end
    vars = @variables begin

     end
   eqs = [   outport.p ~ inport.p
            outport.h ~ PropsSI("H","Q",0,"P",inport.p,fluid)
            outport.mdot ~ inport.mdot
   ]    
   compose(ODESystem(eqs, t, vars, para;name), inport, outport)
end

function ExpansionValve(;name,πc)
    @named inport = CoolantPort()
    @named outport = CoolantPort()
    para = @parameters begin

    end
    vars = @variables begin

     end
   eqs = [   outport.p ~ inport.p/πc
            outport.h ~ inport.h
            outport.mdot ~ inport.mdot
   ]    
   compose(ODESystem(eqs, t, vars, para;name), inport, outport)
end

function Evaporator(;name)
    @named inport = CoolantPort()
    @named outport = CoolantPort()
    para = @parameters begin

    end
    vars = @variables begin
     
     end
   eqs = [  
            outport.p ~ inport.p 
            outport.h ~ PropsSI("H","Q",1,"P",inport.p,fluid)
            outport.mdot ~ inport.mdot
   ]    
   compose(ODESystem(eqs, t, vars, para;name), inport, outport)
end

function MassSink(;name) 
    @named    port = CoolantPort()
    para = @parameters begin
        
    end
    vars = @variables begin
        p(t) 
        mdot(t) 
        h(t) 
     end

   eqs = [
        mdot ~ port.mdot
        p ~ port.p
        h ~ port.h
   ]
   compose(ODESystem(eqs, t, vars, para;name),port)
end

start_T = 300
start_p = 101325
start_h = PropsSI("H","T",start_T,"P",start_p,fluid)
mdot = 10 

@named source = MassSource(source_enthalpy = start_h,source_pressure = start_p,source_mdot = mdot)
@named compressor = Compressor(πc = 10.0)
@named condensor = Condensor()
@named expvalve = ExpansionValve(πc = 10.0)
@named evap = Evaporator()
@named sink = MassSink()
#@unpack p,outport = oneport
eqs = [
    connect(source.port,compressor.inport)
    connect(compressor.outport,condensor.inport)
    connect(condensor.outport,expvalve.inport)
    connect(expvalve.outport,evap.inport)
    connect(evap.outport,sink.port)
]
@named model = ODESystem(eqs, t, systems=[source,compressor,condensor,expvalve,evap,sink])
u0 = []
tspan = (0.0, 1.0)
sys = structural_simplify(model)
prob = ODEProblem(sys,u0,tspan)
 sol = solve(prob,saveat = 1.0)
3 Likes