Modelling of two phase flow in packed bed (continued)

[Chestermiller]

Just in addition - in order to credibly use this model, I suppose I would need to add the $H_2O$ ODEs back in and tune the model with these included? Then remove them after tuning

This is a matter of your best judgment.

 

[casualguitar]

This is a matter of your best judgment.

I can add these back in to see if changing the porosity 0.7 in the above plots was just a lucky coincidence. If it wasn't a coincidence then we should see that the $H_2O$ buildup plots are also approximately correct

In regards to fine tuning the above plots (post #398), can we similarly use the number of tanks here to tune to Tuinier et al? Every other variable seems to be tied up in correlations

I'll work towards the above now unless you think otherwise

 

[Chestermiller]

I can add these back in to see if changing the porosity 0.7 in the above plots was just a lucky coincidence. If it wasn't a coincidence then we should see that the $H_2O$ buildup plots are also approximately correct

oIsn't the porosity 0.3?

In regards to fine tuning the above plots (post #398), can we similarly use the number of tanks here to tune to Tuinier et al? Every other variable seems to be tied up in correlations

I think so. This all is not exact representation of reality, right? It's a model approximation.

 

[casualguitar]

I think so. This all is not exact representation of reality, right? It's a model approximation.

Thats true, I'll use the tank number then

oIsn't the porosity 0.3?

Hmm. The paper doesn't seem to provide a porosity value. I found the PhD thesis that this publication came from online. In the section about this model, they don't seem to give the porosity value, but they do give the porosity for a number of other very similar simulations (using the same model as is in the paper), and the porosity value they use is $0.7$. In addition, when I use their value of 0.7 I get those plots in post #398 that are very similar to theirs. Here's the thesis: https://pure.tue.nl/ws/files/3606277/719418.pdf

If 0.7 is reasonable, then all that remains is to increase n until we (hopefully) get the same output as Tuinier

 

[casualguitar]

Here are the CO2 deposition and temperature profile plots with n=100, not including the $H_2O$ ODEs:

CO2 deposition plot comparisons (our model, Tuinier model):

Temperature profiles:

The deposition profiles are very similar. The temperature profile from our model is missing that first almost constant temperature zone at around 40C. I'm not sure what that is. Presumably its the water freezing but the temperature seems to be a bit high. Besides this it looks pretty good

P.s. My apologies for not providing the same units as Tuinier. I'll get this included in the next version that will have the $H_2O$ ODEs included

 

[casualguitar]

I think so. This all is not exact representation of reality, right? It's a model approximation.

Here is the output after adding $H_2O$ back in. It looks good in my opinion. I ran for n=50 because n=100 was going to take too long (I'll rerun now at n=100 to get smoother plots)

Gas temperature:

CO2 deposition:

H2O Deposition:

Besides plotting the deposition of CO2 and H2O on the same plot, and fixing the scale units, is there anything you see here as being required in order to consider this model now "fit for use"?

I'll start defining some performance parameters for a parametric analysis if not

 

[Chestermiller]

Here is the output after adding $H_2O$ back in. It looks good in my opinion. I ran for n=50 because n=100 was going to take too long (I'll rerun now at n=100 to get smoother plots)

Gas temperature:
View attachment 326962

CO2 deposition:
View attachment 326963

H2O Deposition:
View attachment 326965

Besides plotting the deposition of CO2 and H2O on the same plot, and fixing the scale units, is there anything you see here as being required in order to consider this model now "fit for use"?

I'll start defining some performance parameters for a parametric analysis if not

Check of overall mass balances on CO2 and H2O:: cumulative mass in = deposited mass + mass out?

 

[casualguitar]

Check of overall mass balances on CO2 and H2O:: cumulative mass in = deposited mass + mass out?

My apologies for the delay. It was quite difficult to extract the values of the intermediate variables (molar holdup, molar flow, molar desublimation rate, etc) from the integrator.

So we have -

Cumulative mass in:
$M_{IN,TOTAL} = \dot{m_{in}}*y_{CO_2, IN} * t\tag{1}$

The total molar amount of CO2 leaving the column up to time t is:
$\dot{m_T}$ = $\int_0^t{\dot{m}_n(t')y_{CO2}(t')dt'}$

Total molar holdup of CO2 in the gas phase:
$M_{CO2} = \sum_{j=1}^{n}(P/RT_j)y_{j,CO2}(A_cdz*\epsilon)\tag{4}$

Total molar holdup of CO2 in the solid phase:
$M_{SOLID} = \sum_{j=1}^{n}(M_j)\tag{5}$

(And the same for H2O)

Using these relationships, we get this output:

At long times CO2 in = CO2 out, and at long times CO2 Gas phase molar holdup = CO2 in * number of tanks

And the same trend for H2O, except the times involved are longer as it takes H2O longer to move through the bed:

If the above looks reasonable, I was hoping to define some parameters that could assess the performance of this system, before running the simulation for suitable ranges of initial/boundary conditions Possibly:

Captured Fraction:
$\eta = \frac{M_{SOLID}}{M_{IN,TOTAL}}$

In addition I was interested in defining a parameter that measures the total amount of energy required to 'capture' a given amount of CO2. This parameter would take into account the energy required to cool the bed (I should have this from the other model), and also the energy required to send combustion gas to the packed bed (via compressor). The idea being to assess different conditions versus this "energy required per unit of CO2 captured" parameter. If you think this is reasonable I'll make an initial attempt at defining that performance parameter more accurately

 

[casualguitar]

If you think this is reasonable I'll make an initial attempt at defining that performance parameter more accurately

Here's an inital attempt -

Assuming here that these variables are the ones we have operational control over:
1. Initial bed temperature
2. Inlet gas temperature
3. Bed pressure
4. Inlet CO2 concentration
5. Bed material

It would be useful to know which values of each of the above variables lead to the largest amount of CO2 captured. We could simply run a parametric analysis here to see that

I think it would be more useful though to know which values of each of the above variables lead to the largest amount of CO2 captured $\textbf{per unit of energy input}$ i.e. how much total energy input is required per unit mass of CO2 captured. In this packed bed system, we've got energy input in two forms:
1. Cooling of the packed bed
2. 'Pressurising' of the hot N2/CO2/H2O stream to pass it through the packed bed

Note here we have been running the model at atmospheric pressure, but we could in theory run it at any pressure up to the triple point pressure (about 5atm), which is why I included energy input #2.

So the total energy would be the sum of these (not defined properly yet).

The 'total amount of CO2 captured' could be defined as the maximum amount of CO2 captured. We have data for the amount of CO2 deposited at each point in the bed at each point in time, so we could simply take the maximum value.

Do you think the above is reasonable? If so I'll start defining the energy inputs properly. Very much open to the possibility of using other performance parameters at this stage if you think of any!

 

[Chestermiller]

My apologies for the delay. It was quite difficult to extract the values of the intermediate variables (molar holdup, molar flow, molar desublimation rate, etc) from the integrator.

So we have -

Cumulative mass in:
$M_{IN,TOTAL} = \dot{m_{in}}*y_{CO_2, IN} * t\tag{1}$

The total molar amount of CO2 leaving the column up to time t is:
$\dot{m_T}$ = $\int_0^t{\dot{m}_n(t')y_{CO2}(t')dt'}$

Total molar holdup of CO2 in the gas phase:
$M_{CO2} = \sum_{j=1}^{n}(P/RT_j)y_{j,CO2}(A_cdz*\epsilon)\tag{4}$

Total molar holdup of CO2 in the solid phase:
$M_{SOLID} = \sum_{j=1}^{n}(M_j)\tag{5}$

(And the same for H2O)

Using these relationships, we get this output:
View attachment 327171
At long times CO2 in = CO2 out, and at long times CO2 Gas phase molar holdup = CO2 in * number of tanks

And the same trend for H2O, except the times involved are longer as it takes H2O longer to move through the bed:
View attachment 327172

You should also show a graph of yellow + red + green in comparison to blue to show that CO2 and H2O are conserved. Im surprised that the holdup in the gas within the column is so low, even though all the CO2 and H20 exiting the column exits in the gas streams from the final tank.

If the above looks reasonable, I was hoping to define some parameters that could assess the performance of this system, before running the simulation for suitable ranges of initial/boundary conditions Possibly:

Captured Fraction:
$\eta = \frac{M_{SOLID}}{M_{IN,TOTAL}}$

In addition I was interested in defining a parameter that measures the total amount of energy required to 'capture' a given amount of CO2. This parameter would take into account the energy required to cool the bed (I should have this from the other model), and also the energy required to send combustion gas to the packed bed (via compressor). The idea being to assess different conditions versus this "energy required per unit of CO2 captured" parameter. If you think this is reasonable I'll make an initial attempt at defining that performance parameter more accurately

The fraction captured should pass through a maximum, after which, at some time, the flow should stop.

 

[Chestermiller]

Here's an inital attempt -

Assuming here that these variables are the ones we have operational control over:
1. Initial bed temperature
2. Inlet gas temperature
3. Bed pressure
4. Inlet CO2 concentration
5. Bed material

It would be useful to know which values of each of the above variables lead to the largest amount of CO2 captured. We could simply run a parametric analysis here to see that

I think it would be more useful though to know which values of each of the above variables lead to the largest amount of CO2 captured $\textbf{per unit of energy input}$ i.e. how much total energy input is required per unit mass of CO2 captured. In this packed bed system, we've got energy input in two forms:
1. Cooling of the packed bed
2. 'Pressurising' of the hot N2/CO2/H2O stream to pass it through the packed bed

Note here we have been running the model at atmospheric pressure, but we could in theory run it at any pressure up to the triple point pressure (about 5atm), which is why I included energy input #2.

So the total energy would be the sum of these (not defined properly yet).

The 'total amount of CO2 captured' could be defined as the maximum amount of CO2 captured. We have data for the amount of CO2 deposited at each point in the bed at each point in time, so we could simply take the maximum value.

Do you think the above is reasonable? If so I'll start defining the energy inputs properly. Very much open to the possibility of using other performance parameters at this stage if you think of any!

These all sound like good possibilities. Once you start playing with the model as if it is a real system, you will get other ideas.

 

[casualguitar]

These all sound like good possibilities. Once you start playing with the model as if it is a real system, you will get other ideas.

Thanks Chet. Working on this. Will update with progress in due course

 

[casualguitar]

Hi Chet - just posting my update here to keep this timeline intact -

My thoughts on the axial dispersion coefficient for model 1 -

My thoughts -

We defined $l$, the dispersion length, as:
$$l = \frac{\Delta x}{2}$$

We know:
$$\Delta x = \frac{L}{n}$$

Subbling $\Delta x$ in:
$$l = \frac{L}{2n}$$

We also know $D = u*l$, where D is the dispersion coefficient, so:
$$\frac{D}{u} = \frac{L}{2n}$$

So n, the number of tanks, as a function of the dispersion coefficient (and the equation we could use to find the number of tanks) is:
$$n = \frac{L*u}{2D}$$

Therefore, if I find a correlation for the dispersion coefficient in literature then I can sub this in above. You mentioned we couldnt use a correlation that was a function of the Reynolds number though I dont understand what the problem with using one that is a function of the Reynolds number would be?

Or have I made an incorrect assumption about what it is I'm doing?

 

[casualguitar]

You should also show a graph of yellow + red + green in comparison to blue to show that CO2 and H2O are conserved. Im surprised that the holdup in the gas within the column is so low, even though all the CO2 and H20 exiting the column exits in the gas streams from the final tank.

Just updating with the mole balance check. Plot of CO2 in vs out+gas holdup + solid holdup:

It seems to be ok based on this plot. Also it matches the Tuinier data well:

 

[Chestermiller]

Just updating with the mole balance check. Plot of CO2 in vs out+gas holdup + solid holdup:
View attachment 327496
It seems to be ok based on this plot. Also it matches the Tuinier data well:
View attachment 327497
View attachment 327498

Considering the differences between the two models, these are excellent matches.

 

[Chestermiller]

\Hi Chet - just posting my update here to keep this timeline intact -

My thoughts on the axial dispersion coefficient for model 1 -

My thoughts -

We defined $l$, the dispersion length, as:
$$l = \frac{\Delta x}{2}$$

We don't define the diffusion length as $$l=\Delta x/2$$We choose $\Delta x$ such that $$\Delta x=2l$$This enables us to use backward differencing to simulate the actual amount of dispersion without explicitly including the dispersion.

We know:
$$\Delta x = \frac{L}{n}$$

Subbling $\Delta x$ in:
$$l = \frac{L}{2n}$$

We also know $D = u*l$, where D is the dispersion coefficient, so:
$$\frac{D}{u} = \frac{L}{2n}$$

So n, the number of tanks, as a function of the dispersion coefficient (and the equation we could use to find the number of tanks) is:
$$n = \frac{L*u}{2D}$$

Therefore, if I find a correlation for the dispersion coefficient in literature then I can sub this in above. You mentioned we couldnt use a correlation that was a function of the Reynolds number though I dont understand what the problem with using one that is a function of the Reynolds number would be?

Or have I made an incorrect assumption about what it is I'm doing?

According to Freeze and Cherry, l is a characteristic of the porous medium independent of the flow rate.

 

[casualguitar]

We don't define the diffusion length as $$l=\Delta x/2$$We choose $\Delta x$ such that $$\Delta x=2l$$This enables us to use backward differencing to simulate the actual amount of dispersion without explicitly including the dispersion.

According to Freeze and Cherry, l is a characteristic of the porous medium independent of the flow rate.

I saw this also if you're referring to page 390. It defines it as:
$$D = \alpha * \bar{v} + D^*$$
where $\alpha$ is the thermal diffusivity and $D^*$ is the molecular diffusivity.

They extend this to say that $D = \alpha * \bar{v}$ is a reasonable estimate (as you previously defined it). However they don't seem to give $\alpha$ as a function of anything in particular (or reference anyone that might).

Are you saying here that it needs to be determined experimentally, or am I not understanding how we can approximate it?

 

[casualguitar]

At a bit of a loss for what to do about the dispersivity - but I did find something possibly useful in this paper by Rastegar et al. They cite some correlations for the axial dispersion coefficient that are not dependent on the Reynolds number : https://www.sciencedirect.com/scien...aLlLLqg7a40zkvohwbhpnLD7a6uZqchExe8QqNC-NCgiI

I was thinking that I could use one of these, and divide by the superficial velocity to get the value of $l$?

Here are the correlations they mention, where $D_{ab}$ is the axial dispersion coefficient:

First one:

Second one:

Do you think either of these are valid? If so I'll plot both for a range of superficial velocities

This question of quantifying the dispersion is the last roadblock in the model 1 simulation. After that its just parametric analysis

 

[casualguitar]

According to Freeze and Cherry, l is a characteristic of the porous medium independent of the flow rate.

Ah here we go - I found this in "A Review of Diffusion and Dispersion in Porous Media" by Perkins et al:

This is on the order of the particle diameter as you said, and is not a function of the Reynolds number. Looks promising?

So since we know $D=u*l$, I guess we can say:
$$D = u*(1.75*d_p)$$
So:
$$l = 1.75*d_p$$

And using our previously derived $n = \frac{L}{2*l}$, particle diameter of 5mm and column length of 1m, we get:
$$n \approx 57$$

I'm semi optimistic about the above. Have I done anything outlandish?

 

[casualguitar]

Hey Chet,

Just updating here to say the parametric analysis for the CO2 model is in progress. As of now its set up like the following:
Parametric variables: Bed temperature (-170C to -90C), CO2 inlet mole fraction (0.1, 0.15, 0.2), and material (alumina, stainless steel, and silica). The inlet temperature is 100C.
Performance parameters: Specific Energy (cooling energy required per mole of co2 captured), captured fraction (solid co2 divided by total co2 in)

In addition, I'm doing a similar analysis for standard air (not flue gas). So for the same bed temperatures and materials I'll run a CO2 mole fraction of 0.000414 and an inlet temperature of 15C

 

[casualguitar]

Updating on an initial attempt to reformulate the sharp front model on a mole basis:

Converting their capture step sharp front equations to a mole basis -

Component mass balance for defrost front:
\begin{equation}
\Phi'_1 \omega'_{i1} = \Phi'_2 \omega'_{i2} + A m'_i v_d
\end{equation}
$\Phi'_1$ and $\Phi'_2$ are the molar flow rates after and before the defrost front, respectively (in moles per second), $\omega'_{i1}$ and $\omega'_{i2}$ are the molar fractions of component 'i' after and before the defrost front, respectively. $A$ is the cross-sectional area of the bed through which the gas flows (in square meters). $m_i$ is the number of moles of component 'i' deposited per unit of bed volume (in moles per cubic meter). $v_d$ is the velocity of the defrost front (in meters per second), calculated as the distance the defrost front travels divided by the time taken.

Component mole balance for frost front:
\begin{equation}
\Phi'_0 \omega'_{i0} = \Phi'_1 \omega'_{i1} + A m'_i v_f
\end{equation}

So we can see here that the molar flow out of the defrost front is equal to the molar flow into the frost front.

The frost and defrost front velocities dont change for our model. I am guessing that we will use the $\Delta z$ value to track the position of the front:
\begin{equation}
v_d = \frac{z_{d,2} - z_{d,1}}{\Delta t}
\end{equation}
\begin{equation}
v_f= \frac{z_{f,2} - z_{f,1}}{\Delta t}
\end{equation}

Energy balances on a mole basis:
\begin{equation}
A_{vd} \left[ \rho_s C_{p,s}' (T_2 - T_1) + m_i' \Delta H' \right] = \Phi_2' (T_2 - T_1) \left( \omega_{i2}' C_{p,i}' + \omega_{j2}' C_{p,j}' \right)
\end{equation}
\begin{equation}
A_{vf} \left[ \rho_s C_{p,s}' (T_1 - T_0) - m_i' \Delta H' \right] = \Phi_0' (T_1 - T_0) \left( \omega_{i0}' C_{p,i}' + \omega_{j0}' C_{p,j}' \right)
\end{equation}

where the heat capacities are in $J/mol.K$

Overall mole balance for each front:
\begin{equation}
\Phi'_0 \omega'_{i0} = \Phi'_1 \omega'_{i1} + A m'_i v_f
\end{equation}

\begin{equation}
\Phi'_0 = \Phi'_1 + A m'_i v_f
\end{equation}

We also have the molar desublimation rate equation which might be useful to calculate $m_i$:
\begin{equation}
M_i'' = k_i\frac{Py_i - p_T}{RT}
\end{equation}

Just one question - have I missed a relation here? It looks like we have 7 equations and 8 unknowns. I'm defining the unknowns as: $T$, the molar fluxes out of the fronts $\Phi_0$ and $\Phi_1$, the mole fractions out of the fronts $\omega_0$ and $\omega_1$, the front velocities $v_d$ and $v_f$, and the amount of solid buildup $m_i$ in mol/m3.

If this is not the intended reformulation just let me know and I can reform these equations

 

[casualguitar]

Initial attempt to generate the sharp front heat balances from our equations:

The heat balance for the gas phase is given by:
\begin{equation}
\epsilon \rho_m C_{p,g,m} \frac{\partial T_g}{\partial t} = - \phi_m C_{p,g,m}\frac{\partial T_g}{\partial z} + \frac{\partial}{\partial z}(k_{eff}\frac{\partial T_g}{\partial z}) - q_{g,I}a_s
\end{equation}
And the heat balance for the packed bed:
\begin{equation}

\rho_s(1-\epsilon_g)C_{p,s}\frac{\partial T_b}{\partial t} = q_{I,b}a_s

\end{equation}
Assuming an infinite heat transfer coefficient means that Tg = Tb. Also assuming axial dispersion = 0 we can discount the dispersion term.

If we make these assumptions, and add the two heat balances we get:
\begin{equation}

\epsilon \rho_m C_{p,g,m} \frac{\partial T}{\partial t} + \rho_s(1-\epsilon_g)C_{p,s}\frac{\partial T}{\partial t} = - \phi_m C_{p,g,m}\frac{\partial T_g}{\partial z} + q_{I,b}a_s - q_{g,I}a_s

\end{equation}
We know the relationship between $q_{I,b}$ and $q_{g,I}$ so we can sub out $q_{I,b}$ using:
\begin{equation}

q_{g,I} + \sum\limits_{i=1}^{n_c} M_j''\Delta h_j = q_{I,b}

\end{equation}
Subbing in and cancelling:
\begin{equation}

\epsilon \rho_m C_{p,g,m} \frac{\partial T}{\partial t} + \rho_s(1-\epsilon_g)C_{p,s}\frac{\partial T}{\partial t} = - \phi_m C_{p,g,m}\frac{\partial T_g}{\partial z} + \sum\limits_{i=1}^{n_c} M_j''\Delta h_ja_s

\end{equation}

Assuming an infinite heat transfer coefficient means that all co2 in a spatial element is either sublimated or desublimated depending on the temperature (this happens instantly). So $\sum M''_i a_s$ is actually replaced with $m_i Av$, where A is the CSA and v is the velocity (of the front):
\begin{equation}

\epsilon \rho_m C_{p,g,m} \frac{\partial T}{\partial t} + \rho_s(1-\epsilon_g)C_{p,s}\frac{\partial T}{\partial t} = - \phi_m C_{p,g,m}\frac{\partial T_g}{\partial z} + m_i Av\Delta h_j

\end{equation}
At this point I'm not sure how to finish this out. If we let subscript 0 denote before the defrost front, subscript 1 denote between the defrost and frost fronts, and subscript 2 denote after the frost front, we can say that $\frac{\partial T_g}{\partial z}$ = $T_1 - T_0$ for defrost and $\frac{\partial T_g}{\partial z}$ = $T_2 - T_1$ for frost.

I can recreate the RHS of their heat balance for the defrost front for example knowing that:
\begin{equation}

- \phi_1 C_{p,g,m} = \phi_m (C_{p,i}y_i + C_{p,j}y_j)(T_1 - T_0)

\end{equation}
That leaves the LHS as below, which I am not sure how to simplify to what Tuinier has:
\begin{equation}

\epsilon \rho_m C_{p,g,m} \frac{\partial T}{\partial t} + \rho_s(1-\epsilon_g)C_{p,s}\frac{\partial T}{\partial t} + m_i Av_d\Delta h_i = \phi_m (C_{p,i}y_i + C_{p,j}y_j)(T_1 - T_0)

\end{equation}
They seem to have factored out $Av_d$ but we only have one velocity term on the LHS. Are we going to integrate the time partial derivatives here or do we know the value of this term?

For reference, their LHS is = $Av_d(\rho_s C_{p,s}(T_1-T_0) - m_i\Delta h_i)$, and their RHS is the same as ours

 

[Chestermiller]

$$\epsilon \rho_m C_{p,g,m} \frac{\partial T}{\partial t} + \rho_s(1-\epsilon_g)C_{p,s}\frac{\partial T}{\partial t} = - \phi_m C_{p,g,m}\frac{\partial T_g}{\partial z} + m_i Av\Delta h_j$$

I need to think about this some more. If the deposition term were not there, we Ould have $$\frac{\partial T}{\partial t}+\frac{\phi_mC_{p,g,m}}{\epsilon \rho C_{pig,m}+\rho_s(1-\epsilon)C_{p,s}}\frac{\partial T}{\partial x}=0$$

 

[casualguitar]

I need to think about this some more. If the deposition term were not there, we Ould have $$\frac{\partial T}{\partial t}+\frac{\phi_mC_{p,g,m}}{\epsilon \rho C_{pig,m}+\rho_s(1-\epsilon)C_{p,s}}\frac{\partial T}{\partial x}=0$$

Interesting. I assume this form is useful because it is a form of the advection equation, where everything before the spatial partial derivative is the velocity of the front. This equation seems to apply when we're outside the defrost/frost zone (as $m_i$ is zero)

If we did have the deposition term (which is the case in between the frost and defrost fronts):
\begin{equation}
\frac{\partial T}{\partial t}+\frac{\phi_mC_{p,g,m}}{\epsilon \rho C_{pig,m}+\rho_s(1-\epsilon)C_{p,s}}\frac{\partial T}{\partial x}= m_i Av\Delta h_j
\end{equation}
We get a first-order linear non-homogeneous PDE, which is the same as you had with an additional source term:
\begin{equation}
\frac{\partial T}{\partial t}+v\frac{\partial T}{\partial x}= m_i Av\Delta h_j
\end{equation}
If we're going the road of using a step function then I suppose the general solution is below but it doesn't seem useful:
\begin{equation}
T(x,t) = f(x-vt) + \frac{m_i Av\Delta h_j}{v}
\end{equation}
Can we equate the velocity in the deposition term and the velocity term in the advection equation defined by you above?

 

[casualguitar]

Hi Chet, I just wanted to clear up the wording on the constant dispersion length assumption - are we saying the first or second statement below?
1) the dispersion length $l$ is assumed to be equal to $\frac{\Delta x}{2}$
2) the spatial increment $\Delta x$ is assumed to be equal to $2l$, where $l$ is the dispersion length

I realise that these are very similar. In our discussion it seems that you're saying its "let $\Delta x$ = $2l$", but when we get to this point in our discretisation scheme:
\begin{flalign*}
\phi\frac{\partial y}{\partial z} &= \frac{1}{2}[\phi_{z+\Delta z/2}(\frac{y_{z+\Delta z - y_z}}{\Delta z}) &\\
&\phantom{=}+ \phi_{z-\Delta z/2}(\frac{y_z-y_{z-\Delta z}}{\Delta z})] &\\
\phi C_p\frac{\partial T}{\partial z} &= \frac{1}{2}[\phi_{z+\Delta z/2}C_{p,z+\Delta z/2}(\frac{T_{z+\Delta z} - T_z}{\Delta z}) &\\
&\phantom{=}+ \phi_{z-\Delta z/2}C_{p,z-\Delta z/2}(\frac{T_z-T_{z-\Delta z}}{\Delta z})] &\\
l\frac{\partial}{\partial z}(\phi_m\frac{\partial y}{\partial z}) &= \frac{l}{\Delta z}[\phi_{z+\Delta z/2}(\frac{y_{z+\Delta z} - y_z}{\Delta z}) &\\
&\phantom{=}- \phi_{z-\Delta z/2}(\frac{y_z-y_{z-\Delta z}}{\Delta z})] &\\
l\frac{\partial}{\partial z}(\phi_m C_{p,g,m}\frac{\partial T}{\partial z}) &= \frac{l}{\Delta z}[\phi_{z+\Delta z/2}C_{p,z+\Delta z/2}(\frac{T_{z+\Delta z}-T_z}{\Delta z}) &\\
&\phantom{=}- \phi_{z-\Delta z/2}C_{p,z-\Delta z/2}(\frac{T_z-T_{z-\Delta z}}{\Delta z})] &
\end{flalign*}

We sub in $l$ = $\frac{\Delta x}{2}$ and we dont sub in $2l$ for $\Delta x$:
\begin{equation*}
\epsilon \rho_m C_{p,g,m}\frac{\partial T_g}{\partial t} = \phi_{z - \Delta z/2}C_{p,z - \Delta z/2}(\frac{T_{z-\Delta z} - T_z}{\Delta z}) - q_{g,I,z}a_s
\end{equation*}
\begin{equation*}
\epsilon\frac{\partial \rho_m}{\partial t} = \frac{\phi_{z-\Delta z/2} - \phi_{z+\Delta z/2}}{\Delta z} - \sum\limits_{i=1}^{n_c} M_{j,z}'a_s
\end{equation*}
\begin{equation*}
\epsilon_g \rho_m \frac{\partial y_i}{\partial t} = \phi_{z-\Delta z/2}(\frac{y_{z-\Delta z} - y_z}{\Delta z}) - M_i''a_s + y_{i,z}\sum\limits_{i=1}^{n_c} M_{j,z}''a_s
\end{equation*}

Is it true then to say that we are actually saying "assume $l$ = $\Delta x$/2", and not "assume $\Delta x$ = 2*l?

Apologies for the verbose question

 

[Chestermiller]

Interesting. I assume this form is useful because it is a form of the advection equation, where everything before the spatial partial derivative is the velocity of the front. This equation seems to apply when we're outside the defrost/frost zone (as $m_i$ is zero)

If we did have the deposition term (which is the case in between the frost and defrost fronts):
\begin{equation}
\frac{\partial T}{\partial t}+\frac{\phi_mC_{p,g,m}}{\epsilon \rho C_{pig,m}+\rho_s(1-\epsilon)C_{p,s}}\frac{\partial T}{\partial x}= m_i Av\Delta h_j
\end{equation}
We get a first-order linear non-homogeneous PDE, which is the same as you had with an additional source term:
\begin{equation}
\frac{\partial T}{\partial t}+v\frac{\partial T}{\partial x}= m_i Av\Delta h_j
\end{equation}
If we're going the road of using a step function then I suppose the general solution is below but it doesn't seem useful:
\begin{equation}
T(x,t) = f(x-vt) + \frac{m_i Av\Delta h_j}{v}
\end{equation}
Can we equate the velocity in the deposition term and the velocity term in the advection equation defined by you above?

I'm not so sure this is exactly right. I'm struggling with this, and am having trouble getting the mental concentration to do this right. I'm usually much better at formulation than this.

 

[casualguitar]

I'm not so sure this is exactly right. I'm struggling with this, and am having trouble getting the mental concentration to do this right. I'm usually much better at formulation than this.

As a side point, I'm not exactly sure if the time required here is worth what it would give us? It strengthens the case for saying the model is 'validated' yes, however we have already validated this model at least somewhat with the Tuinier experimental data

 

[Chestermiller]

As a side point, I'm not exactly sure if the time required here is worth what it would give us? It strengthens the case for saying the model is 'validated' yes, however we have already validated this model at least somewhat with the Tuinier experimental data

I have a different perspective. It looks like the simple sharp front model tells 90 % of the story, except for a small amount of dispersion at the fronts. It certainly would be much simpler to do the real modeling and design with the sharp front model, and the results would be easier to interpret physically.