Media properties and component encapsulation in MTK

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