Help with setting up coupled PDE system for adsorption column

Specific Domains Numerics
1

Hi all,

I’m trying to model an adsorption column with a system of PDEs, which have been discretized in space using a WENO finite volume scheme. I’m new to julia and this is my first time simulating a complex system like this. I am currently unsure what the best framework is in which to write out the equations and use solvers from DifferentialEquations.jl. I’ve attempted to write the model both as an ODEProblem() and as a DAEProblem(), but the solvers I’ve tried either crash by reducing the time step to an extremely small value, or run very slowly.

Specifically the PDEs that I’m solving are (after removing most constants):

  • Component Mass balance:
\frac{\partial y_i}{\partial t} = \frac{T}{P} \frac{\partial}{\partial z} \left(\frac{P}{T} \frac{\partial y_i}{\partial z} \right) - \frac{T}{P} \frac{\partial x_i}{\partial t} - \frac{y}{P} \frac{\partial P}{\partial t} + \frac{y}{T} \frac{\partial T}{\partial t}

and (hopefully self-enforced) y_1 + y_2 + y_3 = 1.

  • Total Mass balance:
\frac{\partial P}{\partial t} = - T \frac{\partial }{\partial Z} (P v / T) - T \sum_{i=1}^2 \frac{\partial x_i}{\partial t} + \frac{P}{T} \frac{\partial T}{\partial t}
  • Solid phase balance:
\frac{\partial x_i}{\partial t} = \alpha_i (x_i^*(T, P) - x_i)
  • Column energy balance:
\frac{\partial T}{\partial t} = \frac{\partial^2 T}{\partial Z^2} - \frac{\partial}{\partial Z} (v P) - T \sum_{i=1}^2 \frac{\partial x_i}{\partial t} - \Omega(x_1, x_2) \frac{dP}{dt}
  • Pressure drop
- \frac{dP}{dZ} = a v + b v |v|

These equations and boundary/initial conditions are taken directly from the paper Multiobjective Optimization of a Four-Step Adsorption Process for Postcombustion CO2 Capture Via Finite Volume Simulation.

The only purely algebraic equation is the pressure drop relation, which in theory could be solved separately at each time step. However, because the system is so tightly coupled, I’m not sure whether I should be formulating this as an ODEProblem or a DAEProblem. I also tried writing the system in mass matrix form, but it seems that isn’t supported when the coefficients in front of the time derivatives depend on the state variables.

The other issue I’m having is that the simulations without the solid phase balance equation (dx_i/dt = 0) are using up an enormous amount of memory, e.g. 45.900171 seconds (3.93 G allocations: 62.034 GiB, 12.26% gc time, 6.25% compilation time). Even though I’m using preallocated buffers for intermediate arrays, this seems excessive. I’m not sure what’s causing such high memory usage.

My project can be found at GitHub - shivang57721/adsorption_process: Numerical simulation of adsorption. I apologize that it’s not a ‘minimal example’, but I’m not really sure how to best isolate the problem. I would greatly appreciate any advice on debugging memory/stability issues, or maybe more general suggestions for organizing coupled PDEs in julia. For example, I see that the package ModelingToolkit.jl might be useful, but I’m not really sure.

2

Interesting problem. For the mass matrix formulation you could use an internal energy formulation so that you have \frac{\partial U}{\partial t} instead of \frac{P}{T}\frac{\partial T}{\partial t} and introduce an algebraic equation for U(P,T). In principle then it should be possible to write the full system as in the paper as a DAE system and to use some stiff DAE capable solver.

That said, we tackled a slightly similar problem in this paper: https://doi.org/10.1016/j.cej.2025.162027 with VoronoiFVM.jl (Disclaimer: I am the main author). This could work for the upwind discretization, but probably not for WENO et al.

1 Like