MethodOfLines.jl
MethodOfLines.jl copied to clipboard
PDESystem & arguments
Trying out MethodOfLines.jl for industrial problem
I'm trying to make MOL work for gas-lifted oil production, where I have a liquid and gas entering the distributed volume via valves, and the mixture leaving via a valve. Essentially, one could say that I consider a PDAE = Partial and Differential Algebraic Equation. And in my case, the AEs can be found both in the eqs
part and in the bcs
part.
Thus, I have boundary conditions typically of form:
In other words: some of my BC equations relate distributed variables (form: u(..)
) to boundaries, e.g.,
u(t,0) ~ 0
while others involve temporal functions (v(t)
), e.g.,
v ~ sqrt(u(t,L) - 1)
[as abstracted examples].
By trial and error, I have found that MethodOfLines, does not allow me to define equations for such intermediate variables (e.g., v(t)
) in the boundary conditions.
To get around this problem, I have to eliminate these intermediate variables (e.g., v(t)
) from the BCs, leading to horribly complicated BCs (see my real examples above) which are error prone to express.
In the core form of using PDESystem
:
@named gas_lift = PDESystem(eqs,bcs,domains,ivs,dvs,params)
it seems that I can not use the vector of variables (dvs
above) in the same form as when I define the variables. What I mean is:
vars = @variables begin
# Distributed variables
m̌_ℓ(..), [description = "Liquid bulk density, kg/m^3"]
ṁ_ℓ(..), [description = "Liquid mass flow rate, kg/s"]
m̌_g(..), [description = "Gas bulk density, kg/m^3"]
...
can not be used as the dvs
argument. So I have to recreate the same information (more or less) in a new vector, where it seems like have to:
- strip off the description
- replace
(..)
with(t,x)
- change e.g.
x(t)
tox
Not a horrible problem, but somewhat inconvenient.
It is not clear to me what form params
must/should have. Following the more modern versions of ModelingToolkit, I have defined parameters as follows:
params = @parameters begin
g = 9.81, [description = "Acceleration of gravity, m/s^2"]
R = 8.31, [description = "Ideal gas constant, J/(K*mol)"]
T = 280, [description = "Temperature, assumed constant, K"]
...
It seems like this doesn't work. Is this because...
- I need to strip of the
description
etc.? - The examples in the MethodOfLines documentation uses the form
[g => 9.81, ...]
. Is this the cause for the problem?
Describe the solution you’d like
In ModelingToolkit
, a typical working call to ODESystem
is:
@named model = ODESystem(eqs, iv, dvs, params)
where iv
typically is t
, while dependent variables dvs
can be described in the "modern" MTK form, and the same for params
.
It would be very nice if PDESystem
could use the same basic form. Currently:
@named gas_lift = PDESystem(eqs,bcs,domains,ivs,dvs,params)
The call structure is very similar to that of ODESystem
, but where eqs
is expanded with bcs
. It would be very nice if:
-
params
inPDESystem
could be used similarly as inODESystem
, i.e., with description, macro@parameters
, etc. Maybe it can?? But it doesn't work in my case. -
dvs
inPDESystem
coulduse the variables defined via@variables
. There might be a problem with that, mainly due to the way dependent variables with more than one independent variable are specified: asu(..)
-- while thedvs
in thePDESystem
argument seem to require that one explicitly specify the independent variablesu(t,x)
.
It is possible to add purely temporal variables in the dvs
. It is not clear to me where the initial conditions for those should be given if they describe the solution of an ODE -- in the @variables
macro, or in the bcs
?. I have used temporal variables in a much simpler example, and equated them with external forcing functions -- that works.
Finally, does the functions unknowns
(or formerly states
), observed
, etc. work with the structure created by PDESystem
?
Describe alternatives you’ve considered
Some AEs that naturally relate to the BCs can be moved to eqs
. But not all (at least not in a simple way) ... consider:
Here, Y_a2t and Dp_at2 are only defined at x=0. So I need to describe these at the BC using temporal variables. OK -- a possible alternative could be to define them as distributed variables, i.e.,
Y_a2t(..)
and Dp_a2t(..)
and compute the variables in the eqs
section, where they are never really used, and then just apply the values Y_a2t(t,0)
and Dp_a2t(t,0)
in the BCs.
But this seems like a waste in that one would artificially compute Y_a2t
and Dp_a2t
over the entire interior (x=(0,L)), but only need the values at the boundaries.
Additional context
Add any other context or screenshots about the feature request here.