Hi, I’m starting to use ModelingToolkit.jl for chemical processes using an acausal approach. Here I’m showing an illustrative example of a simple mixing system with heating. Any feedback to make the model more efficient and general would be appreciated.
Two process streams are mixed to produce one of the feeds for our chemical reactor. After mixing, the blended stream is fed to a heating vessel before being sent to the reactor
Mathematical model
Acausal model
using ModelingToolkit, DifferentialEquations
using ModelingToolkit: t_nounits as t, D_nounits as D
# A port can be an inlet or an outlet
@connector Port begin
c(t) # concentration of component A
T(t) # temperature
q(t), [connect = Flow] # heat
v(t), [connect = Flow] # volumetric flow rate
end
@constants Tref = 25 # Reference temperature, °C
# Components
# Stream n
@mtkmodel Stream_n begin
@parameters begin
c # g/L
v # L/min
T # C
Cp = 1 # cal/g C
ρ = 1000 # g/L
end
@components begin
s = Port()
end
@equations begin
s.c ~ c
s.T ~ T
s.v ~ -v # flow is leaving
s.q ~ -ρ*v*Cp*(T - Tref)
end
end
# Mixing tee (two inlets, one outlet)
@mtkmodel MixT begin
@parameters begin
Cp = 1 # cal/g C
ρ = 1000 # g/L
end
@variables begin
c(t) # outlet concentration
v(t) # outlet volumetric flow rate
T(t) # outlet temperature
q(t) # outlet energy flow rate
end
@components begin
in₁ = Port() # inlet 1
in₂ = Port() # inlet 2
out = Port() # outlet
end
@equations begin
out.c ~ c
out.T ~ T
out.v ~ -v # flow is leaving
out.q ~ -q
q ~ ρ*v*Cp*(T - Tref)
v ~ in₁.v + in₂.v # overall mass balance
v*c ~ in₁.v*in₁.c + in₂.v*in₂.c # component mass balance
q ~ in₁.q + in₂.q # energy balance
end
end
# Tank
@mtkmodel Tank begin
@parameters begin
V = 150 # volume, L
Cp = 1 # cal/g C
ρ = 1000 # g/L
Q = 1.05e6 # cal/min
end
@variables begin
c(t)
v(t)
T(t)
q(t)
u(t)
end
@components begin
inlet = Port()
end
@equations begin
inlet.v ~ v # outlet flow rate = inlet flow rate
q ~ ρ*v*Cp*(T - Tref)
u ~ ρ*V*Cp*(T-Tref) # internal energy
# ODE system
V*D(c) ~ inlet.v*inlet.c - inlet.v*c
D(u) ~ inlet.q - q + Q
end
end
# Connecting and simulating the process
@mtkmodel MixProcess begin
@components begin
stream1 = Stream_n(c=1,v=20,T=25)
stream2 = Stream_n(c=4,v=10,T=55)
mixt = MixT()
tank = Tank(V=150,Q=1.05e6)
end
@equations begin
connect(stream1.s,mixt.in₁)
connect(stream2.s,mixt.in₂)
connect(mixt.out,tank.inlet)
end
end
@mtkbuild mix_model = MixProcess()
# Define the ODE (DAE) problem
c0 = 3.0 # initial concentration in the tank
T0 = 80.0 # initial temperature of tank
u0 = [mix_model.tank.c => c0, mix_model.tank.T => T0]
guess = [mix_model.mixt.c => c0,mix_model.mixt.T => T0]
prob = ODEProblem(mix_model,u0,(0.0,30.0);guesses = guess)
# Solve problem
sol = solve(prob)
using Plots
plot(sol;idxs=mix_model.tank.T)
plot(sol;idxs=mix_model.tank.c)