Open juanmvenegas opened 1 year ago
xtalax
commented
1 year ago I am looking at this now, but I think I see the issue - you need to supply explicit initial conditions for all dependent variables, we don't have a setup for inferring these just yet.
All the variables that are missing don't have initial conditions supplied for them. Does this help?
Please ask away if you run in to any more problems :+1:
juanmvenegas
commented
1 year ago Thanks for your feedback on the issue! I included BCs and ICs for all variables in system and the discretization seems to work now. Unfortunately, the model does not solve. I get the following solver error message using Rodas4P or Rosenbrock23:
PDESolution:
Return Code:
Unstable
Dependent variables:
CMeOH(z, t): (20, 1) sized solution
CCH2Os(z, t): (20, 1) sized solution
CCH2O(z, t): (20, 1) sized solution
CCO(z, t): (20, 1) sized solution
CH2O(z, t): (20, 1) sized solution
k3(z, t): (20, 1) sized solution
CO2s(z, t): (20, 1) sized solution
CMeOHs(z, t): (20, 1) sized solution
T(z, t): (20, 1) sized solution
Ts(z, t): (20, 1) sized solution
CH2Os(z, t): (20, 1) sized solution
R2(z, t): (20, 1) sized solution
k1(z, t): (20, 1) sized solution
k2(z, t): (20, 1) sized solution
CCOs(z, t): (20, 1) sized solution
CO2(z, t): (20, 1) sized solution
R1(z, t): (20, 1) sized solution
Domain:
ERROR: BoundsError: attempt to access 1-element Vector{Float64} at index [2]I am not really sure what this is telling me, but I assume this is related to the boundary conditions of the system still not being correctly specified? it seems to show the discretized problem at the first time point.
I have tried a few things to no avail:
use_ODAE=true but this cause an error immediately after code execution that was different: ERROR: MethodError: no method matching (DAEProblem{true})(::ODESystem)jac=true when creating ODEProblemAny feedback on how to proceed debugging would be appreciated! If I should make this a new issue since the discxretization was resolved, let me know. I am including the latest code below for convenience:
using ModelingToolkit, MethodOfLines, DifferentialEquations, DomainSets, Latexify, LaTeXStrings, NonlinearSolve, ForwardDiff
#parameters of problem
@parameters z t dₜ dₚᵥ ϵ uₛ ρf cₚ Tin Tw ΔH₁ ΔH₂
@parameters Peₕᵣ Peₘᵣ Bi Uw kf hfs Deₚ λₚ av
#Parameters with a default value
const CO2ₒ=34.0 #mol/m3
const CMeOHₒ=1.74 #mol/m3
const L=0.7
const R = 8.314 #J/mol/K
#Variables, seprated in multiple lines for readability
@variables CO2(..) CMeOH(..) CCH2O(..) CCO(..) CH2O(..) P(..) T(..)
@variables CO2s(..) CMeOHs(..) CCH2Os(..) CCOs(..) CH2Os(..) Ts(..)
@variables R1(..) R2(..) k1(..) k2(..) k3(..)
#Differential terms
Dt = Differential(t)
Dz = Differential(z)
Dzz = Differential(z)^2
#PDE equations of problem
eqs = [
#Equations of gas phase species and bulk temperature
Dt(CO2(z,t)) ~ kf*av*(CO2s(z,t)-CO2(z,t)) - uₛ*Dz(CO2(z,t))
Dt(CMeOH(z,t))~ kf*av*(CMeOHs(z,t)-CMeOH(z,t))-uₛ*Dz(CMeOH(z,t))
Dt(CCH2O(z,t))~ kf*av*(CCH2Os(z,t)-CCH2O(z,t))-uₛ*Dz(CCH2O(z,t))
Dt(CCO(z,t)) ~ kf*av*(CCOs(z,t)-CCO(z,t))-uₛ*Dz(CCO(z,t))
Dt(CH2O(z,t)) ~ kf*av*(CH2Os(z,t)-CH2O(z,t))-uₛ*Dz(CH2O(z,t))
Dt(T(z,t)) ~ hfs*av*(Ts(z,t)-T(z,t)) - 4*Uw/dₜ * (T(z,t) - Tw)-uₛ*ρf*cₚ*Dz(T(z,t))
##Reaction terms (non isothermal)
k1(z,t) ~ 125E7*exp(-79496/R/Ts(z,t))
k2(z,t) ~ 1.12*exp(-8368/R/Ts(z,t))
k3(z,t) ~ 54E5*exp(-66944/R/Ts(z,t))
R1(z,t) ~ k1(z,t)*CMeOHs(z,t)^0.5/(1+k2(z,t)*CMeOHs(z,t)^0.5)
R2(z,t) ~ k3(z,t)*CCH2Os(z,t)^0.5/(1+0.2*CCH2Os(z,t)^0.5)
#Surface model terms, coupled to gas phase via algebraic equations
kf*av*(CO2s(z,t)-CO2(z,t)) ~ -0.5*R1(z,t)-0.5*R2(z,t)
kf*av*(CMeOHs(z,t)-CMeOH(z,t))~ -1.0*R1(z,t)
kf*av*(CCH2Os(z,t)-CCH2O(z,t))~ +1.0*R1(z,t)-1.0*R2(z,t)
kf*av*(CH2Os(z,t)-CH2O(z,t)) ~ +1.0*R1(z,t)+1.0*R2(z,t)
kf*av*(CCOs(z,t)-CCO(z,t)) ~ +1.0*R2(z,t)
hfs*av*(Ts(z,t)-T(z,t)) ~ ΔH₁*R1(z,t) + ΔH₂*R2(z,t) #heat transfer between bulk and surface
]
#Boundary conditions, Differential BCs currently commented out for troubleshooting
bcs =[
#z = 0 BCs
CO2(0.0,t)~CO2ₒ,
CMeOH(0.0,t)~CMeOHₒ,
CCH2O(0.0,t)~0.0,
CCO(0.0,t)~0.0,
CH2O(0.0,t)~0.0,
T(0.0,t)~517.0,
k1(0.0,t) ~ 0.0,
k2(0.0,t) ~ 0.0,
k3(0.0,t) ~ 0.0,
R1(0.0,t) ~ 0.0,
R2(0.0,t) ~ 0.0,
CO2s(0.0,t)~0.0,
CMeOHs(0.0,t)~0.0,
CCH2Os(0.0,t)~0.0,
CCOs(0.0,t)~0.0,
CH2Os(0.0,t)~0.0,
Ts(0.0,t)~517.0,
#z=L BCs
Dz(CO2(L,t))~0.0,
Dz(CMeOH(L,t))~0.0,
Dz(CCH2O(L,t))~0.0,
Dz(CCO(L,t))~0.0,
Dz(CH2O(L,t))~0.0,
Dz(T(L,t))~0.0,
#t=0 ICs
CO2(z,0.0)~0.0,
CMeOH(z,0.0)~0.0,
CCH2O(z,0.0)~0.0,
CCO(z,0.0)~0.0,
CH2O(z,0.0)~0.0,
T(z,0.0)~517.0,
k1(z,0.0) ~ 0.0,
k2(z,0.0) ~ 0.0,
k3(z,0.0) ~ 0.0,
R1(z,0.0) ~ 0.0,
R2(z,0.0) ~ 0.0,
CO2s(z,0.0)~0.0,
CMeOHs(z,0.0)~0.0,
CCH2Os(z,0.0)~0.0,
CCOs(z,0.0)~0.0,
CH2Os(z,0.0)~0.0,
Ts(z,0.0)~517.0
]
domains = [z ∈ Interval(0.0,L),
t ∈ Interval(0.0,100.0)]
#Parameter map
pars = [
dₜ=>0.0266,
dₚᵥ=> 0.0046,
ϵ=> 0.5,
uₛ=>2.47,
ρf=>1.018,
cₚ=> 952,
Tin=>517,
Tw=> 517,
ΔH₁=> -158700,
ΔH₂=> -158700,
Peₕᵣ=> 8.6,
Peₘᵣ=>6.6,
Bi=> 5.5,
Uw=>220,
kf=> 0.25,
hfs=> 400,
Deₚ=> 4.9E-6,
λₚ=>2,
av=>100
]
#Variable array
vars = [
CO2(z,t),
CMeOH(z,t),
CCH2O(z,t),
CCO(z,t),
CH2O(z,t),
T(z,t),
CO2s(z,t),
CMeOHs(z,t),
CCH2Os(z,t),
CCOs(z,t),
CH2Os(z,t),
Ts(z,t),
k1(z,t),
k2(z,t),
k3(z,t),
R1(z,t),
R2(z,t)]
@named pdesys = PDESystem(eqs,bcs,domains,[z,t],vars,pars)
# Method of lines discretization, keeping it small during troubleshooting
dz = 20
order = 2
discretization = MOLFiniteDifference([z=>dz],t,approx_order=order,grid_align=edge_align,use_ODAE=false)
prob = discretize(pdesys,discretization) # This gives an ODEProblem since it's time-dependent
sol = DifferentialEquations.solve(prob,Rodas4P(),reltol=1E-8,abstol=1E-8)
Hello,
I am attempting to solve a time-dependent packed bed reactor model with MethodOfLines. This model will have a length and time dependence (z, t domains). I am running into trouble during the ODEProblem creation step, as it seems that the discretized equations generated by MOLFiniteDifference contain variables that are not correctly carried over to the ODEPRoblem generation process. I assume there is some part of the structure of MOLFiniteDifference or the discretize function that I am missing. I will copy the code below, but before I'll add part of the error stack trace for reference:
Prior to the discretization, all relevant variables are correctly structured in the PDESystem as a function of z and t:
It seems that upon discretization and generating the corresponding X(t) arrays each variable for each z position, these are not carried over. I could not find information in the documentation that could help me troubleshoot this issue, but may have searched with improper terminology for the field.
The code of this problem is below:
I noticed that the error only shows variables that do not have a differential term in them. Are there limitations on how algebraic equations must be stated in MethodOfLines.jl?
I have tried to discretize with different number of points, different order and different boundary conditions (my code currently has the Neuman BCs commented out, but with them the system still runs into the same error). Is there a specific order on how to specify discretized and continuous domains in my variables (e.g. (z,t) vs t,z)? Are my equations with d/dt terms on the LHS and d/dz terms on the RHS correct? Finally, my BCs only cover the variables that have differential terms as I assume that the initial and boundary values of the other terms in my equations would be calculated directly from these. Do I need to explicitly state BCs for all variables in the system?
Thank you for your time and help in troubleshooting this model.