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)