Adding Differential Equations to User Defined Files

Hello MFiX Developers and Users,

I hope this message finds you well! Just to note, I’m using MFiX-24.4.

I’m posting here to inquire about (potentially) adding an ordinary different equation (ODE) to a user file - specifically, usr_rates.f (I haven’t swapped to the new way of implementing rate expressions as of yet). Nonetheless, the ODE I want to add is apart of some rate expressions for my system and currently I’m using an algebraic expression as an approximation. However, I believe this approximation is leading to slightly inaccurate results and was hoping to implement the ODE in some way. From other posts I’ve seen, I’m under the impression that I would need be re-write source code to accomplish this. Therein, I’m posting here to see if there’s any work around.

As some insight (I can upload a .zip file of my project to this thread if needed), I have four (4) non-catalytic reactions that are functions of a conversion term. Represented as an ODE, this conversion term is a function of these four (4) non-catalytic reactions. Specifically,

dX/dt = (2r_1(X) + r_2(X) + r_3(X) + r_4(X))/C,

where X is a form of solids conversion, r_i are the rate expressions for the i ∈ {1,2,3,4} non-catalytic reactions in the system, and C is a constant. I first thought of applying an approximation to the ODE to transform it into an algebraic expression (different then what I’m currently simulating with). While possible, it would be difficult (or impossible) to isolate X.

I’m not sure if it’s possible, but I wanted to post here to see if there are any work around outside of re-writing source code. I greatly appreciate any insight.

Can you give a bit more detail about what the expressions r_i are, what X represents (is it a mass fraction?), and what algebraic expression you have tried?

Hi @cgw, thanks for the reply. You can find more information below.

  1. Here’s a photo of the 4 non-catalytic reactions (don’t mind that they’re numbered 2-5. There’s a reaction that comes before these but it’s not considered in the reactor that I’m modeling):

Here, a0 is the initial specific surface area of the oxygen carrier (OC) (a constant), k_i (i ∈ {2,3,4,5} represents the chemical reaction rate constant, M_NiO and M_Ni are the molecular weights of nickel oxide (NiO) and nickel (Ni), respectively (constants), W_NiO,0 is the initial NiO concentration (a constant), Ci is gas concentration for i ∈ {CH4, H2,CO}, and X, and (1-X), are related to conversion which I’ll explain in (2.).

  1. X denotes NiO conversion.

  2. The algebraic expression I’m currently simulating with is an expression for (1-X) (this term has a few different names - e.g., solid conversion and active surface area). Here’s a photo of the expression I was using, derived by Levenspiel:

rc, rp, and Z_NiO aren’t used in the model, so I won’t discuss them here. Y_NiO and Y_Ni are the mass fractions of NiO and Ni, respectively, while p_NiO and p_Ni are the densities of NiO and Ni, respectively. Therein, 1 - (expression on the right) = X.

If there’s any more information I can provide, please let me know. Thanks for the assistance!

Just coming back to this with a potential idea (I haven’t tested this as of yet, but is something that I thought of): using one of the other user-definable files (e.g., usr1 or usr2), I’d be able to represent the ODE as an algebraic expression (since the files are called every timestep and/or iteration). The problem I’m still facing is that the resulting expression is nonlinear wherein the variable of interest (X) might be hard to isolate. Not sure if this is a viable solution, but I wanted to document my thought process in case it spirals into a solution.

I’ve come up with another idea for the time being: treating this implicit ODE as an explicit ODE. Wherein, the dependent variable on the left is different than the dependent variable on the right (specifically, I’m using an older value (time wise) to compute a new one, then storing the new on as the old one, and rinse and repeat until simulation completion). While this is a VERY poor approximation (might lead to instability, large error, etc.), I wanted to conduct this experiment to see if it would bridge any gaps for me (i.e., get closer to a plausible solution). When implementing the code, it appears that the initial value I set is used in all the calculations rather than just the first which is causing very low consumption and generation rates. This likely might be a syntax problem - specifically, storing a new value as an old value and using it in a future computation. Thus, I’m wondering if someone could provide some insight on how to get around this.

I will note that I’ve tried “tricking” the software into only running the initial condition at t = 0 (if statement), but the same problem arises. I’ve tried printing the current values of the arrays I’ve created but I get an error regarding memory. The attached files are what I’ve created as an attempt to approximate this ODE using explicit euler (EE).

It’s also important to mention that the older value is what I’m using as my approximation for the 4 noncatalytic reactions as well. Therein, the old value is what’s used for everything - computing a new conversion at the next time step and computing the rate expressions at the current time step. Any insight into this matter would be much appreciated until a new solution is found. I’m not looking for anyone to run the simulation (if need be, let me know), just seeing if there’s a problem with my syntax. Thanks for the help!

usr0.f (2.2 KB)
usr1.f (12.9 KB)
usr_mod.f (4.8 KB)
usr_rates.f (48.3 KB)