# 1D Convection Diffusion Equation - Inlet Mixing Effect

• A
• HumanistEngineer
In summary, the conversation is about solving a 1D convection-diffusion equation using a Matlab code to model a sensible stratified storage tank. The issue being discussed is the abnormal temperature increase that occurs when including an additional mixing effect in the simulation. The person is seeking help in understanding why this problem is occurring and if there is a solution to it. The conversation also includes some details about the code and the parameters used in the simulation.
HumanistEngineer
I have a working Matlab code solving the 1D convection-diffusion equation to model sensible stratified storage tank by use of Crank-Nicolson scheme (without εeff in the below equation).

As indicated by Zurigat et al; there is an additional mixing effect having a hyperbolic decaying form from the top of the tank to the bottom (at the inlet we have a high rate of mixing). So I try to involve εeff as a magnifying effect on the thermal diffusivity as shown in the equation above. Please see the section Inlet Mixing in the Matlab code below.

The problem is that: there occurs an abnormal temperature increase more than the coming inlet temperature whenever I involve the mixing effect as a hyperbolic decaying form (please see pictures below).

Would you please tell me why I have this abnormal temperature increase problem? I tried also Forward-Time Central Space method and CN with upwind treatment but at all approaches, I encountered this abnormal issue. Any solution please? I am not so experienced with the finite difference methods!

With εeff involved, the problem with abnormal temperature increase encountered:

The case without εeff, the solution of convection-diffusion equation for a stratified sensible storage tank by use of Crank-Nicolson scheme:

Matlab Code:

Code:
function [T,dh,dt] = CDR_CrankNicolson_Exp_InletMixing(t_final,n)
%% 1D Convection-Diffusion-Reaction (CDR) equation - Crank-Nicolson scheme
%   solves a d2u/dx2 + b du/dx + c u = du/dt

%% References
% pg. 28 - http://www.ehu.eus/aitor/irakas/fin/apuntes/pde.pdf
% pg. 03 - https://www.sfu.ca/~rjones/bus864/notes/notes2.pdf
% http://nptel.ac.in/courses/111107063/24

tStart = tic; % Measuring computational time

%% INPUT
% t_final   : overall simulation time [s]
% n         : Node number [-]

%% Tank Dimensions
h_tank=1;       % [m]       Tank height
d_tank=0.58;     % [m]       Tank diameter
vFR=11.75;        % [l/min]   Volumetric Flow Rate
v=((vFR*0.001)/60)/((pi*d_tank^2)/4);    % [m/s]     Water velocity

dh=h_tank/n;       % [m] mesh size

dt=1;             % [s] time step
maxk=round(t_final/dt,0);    % Number of time steps

% Constant CDR coefficients | a d2u/dx2 + b du/dx + c u = du/dt
a_const=0.15e-6;     % [m2/s]Thermal diffusivity
b=-v;                        % Water velocity
c=0;                         %! Coefficient for heat loss calculation

%% Inlet Mixing - Decreasing hyperbolic function through tank height
eps_in=2;             % Inlet mixing magnification on diffusion
A_hyp=(eps_in-1)/(1-1/n);
B_hyp=eps_in-A_hyp;

Nsl=1:1:n;

eps_eff=A_hyp./Nsl+B_hyp;

a=a_const.*eps_eff;         % Case 0
% a=ones(n,1)*a_const;      % Case 1
% a=ones(n,1)*a_const*20;    % Case 2

%% Initial condition - Tank water at 20 degC
T=zeros(n,maxk);

T(:,1)=20;

%% Boundary Condition - Inlet temperature at 60 degC
T(1,:)=60;

%% Formation of Tridiagonal Matrices
%   Tridiagonal Matrix @Left-hand side
subL(1:n-1)=-(2*dt*a(1:n-1)-dt*dh*b);             % Sub diagonal - Coefficient of u_i-1,n+1
maiL(1:n-0)=4*dh^2+4*dt*a-2*dh^2*dt*c;          % Main diagonal - Coefficient of u_i,n+1
supL(1:n-1)=-(2*dt*a(1:n-1)+dt*dh*b);             % Super diagonal - Coefficient of u_i+1,n+1
tdmL=diag(maiL,0)+diag(subL,-1)+diag(supL,1);
%   Tridiagonal Matrix @Right-hand side
subR(1:n-1)=2*dt*a(1:n-1)-dt*dh*b;                % Sub diagonal - Coefficient of u_i-1,n
maiR(1:n-0)=4*dh^2-4*dt*a-2*dh^2*dt*c;          % Main diagonal - Coefficient of u_i,n
supR(1:n-1)=2*dt*a(1:n-1)+dt*dh*b;               % Super diagonal - Coefficient of u_i+1,n
tdmR=diag(maiR,0)+diag(subR,-1)+diag(supR,1);

%% Boundary Condition - Matrices
tdmL(1,1)=1; tdmL(1,2)=0;
tdmL(end,end-1)=0; tdmL(end,end)=1;

tdmR(1,1)=1; tdmR(1,2)=0;
tdmR(end,end-1)=1; tdmR(end,end)=0;

MMtr=tdmL\tdmR;

%% Solution - System of Equations
for j=2:maxk % Time Loop
Tpre=T(:,j-1);

T(:,j)=MMtr*Tpre;

if T(end,j)>=89.9
T(:,j+1:end)=[];
finishedat=j*dt;
ChargingTime=sprintf('Charging time is %f [s]', finishedat)
tElapsed = toc(tStart);
SimulationTime=sprintf('Simulation time is %f [s]',tElapsed)
return
end
end

end

#### Attachments

• 2018-07-24_14_43_59-_Zurigat_et_al_AE_1988.pdf.png
2.3 KB · Views: 1,411
• 2018-07-24_14_58_12-_Figures_-_Figure_1.png
19.9 KB · Views: 1,750
• 2018-07-24_14_19_45-_Figures_-_Figure_1.png
18.6 KB · Views: 1,877
It's not clear what you did. What is the difference between the solid curves and the dotted curves?

Sorry Chestermiller, my mistake. The dotted curves are from an experiment while solid lines are from the simulation (the Matlab code given). So the aim of the simulation is to obtain the temperature distribution through the tank height at one dimension by solving the convection-diffusion equation dT/dt = -v dT/dx + alpha d2T/dx2. I could obtain my solution when Crank-Nicolson scheme is applied if eps_eff (inlet mixing effect) not included but far different than the experimental data (please see the second graph). After applying eps_eff as an hyperbolic decaying form as from the top of the tank to the bottom, there occurs an abnormal temperature increase in the solution (please see the first graph). I try to understand/get rid off these abnormal temperature increase so asking you.

I don't want to slog through your code to try to find an error. I'm a little confused about your spatial finite difference scheme. First you mention central differencing, and then you mention upwind differencing. Are you saying that you use central differencing on the diffusion term and upwind differencing on the advection term?

Chestermiller said:
I don't want to slog through your code to try to find an error. I'm a little confused about your spatial finite difference scheme. First you mention central differencing, and then you mention upwind differencing. Are you saying that you use central differencing on the diffusion term and upwind differencing on the advection term?

For the graphs given, I used the Crank-Nicolson scheme (the graphs given are from that). However, I also mentioned that in other attempts (results not given here) I also tried different schemes (i.e. Forward time central space, Crank-Nicolson with Upwind Treatment (the hybrid scheme with the Upwind scheme on the convection/advection term), and Backward Euler) and encountered the same problem. So at each of these difference schemes applied, the solution was obtained with abnormal temperature increase (at each without the involvement of εeff, the results, the temperature distribution through the tank height is reasonable).

My recent discovery on this issue is that: when the 'effective diffusivity factor - εeff ' has a hyperbolic decaying form (high at the top unity at the bottom) this issue with abnormal temperature increase occurs (the first graph). This is because I just give a try of multiplying a constant εeff with the thermal diffusivity α. Even at high values of constant εeff there was no abnormal temperature increase obtained but a high diffusion through the tank height, as physically expected.

HumanistEngineer said:
For the graphs given, I used the Crank-Nicolson scheme (the graphs given are from that). However, I also mentioned that in other attempts (results not given here) I also tried different schemes (i.e. Forward time central space, Crank-Nicolson with Upwind Treatment (the hybrid scheme with the Upwind scheme on the convection/advection term), and Backward Euler) and encountered the same problem. So at each of these difference schemes applied, the solution was obtained with abnormal temperature increase (at each without the involvement of εeff, the results, the temperature distribution through the tank height is reasonable).

My recent discovery on this issue is that: when the 'effective diffusivity factor - εeff ' has a hyperbolic decaying form (high at the top unity at the bottom) this issue with abnormal temperature increase occurs (the first graph). This is because I just give a try of multiplying a constant εeff with the thermal diffusivity α. Even at high values of constant εeff there was no abnormal temperature increase obtained but a high diffusion through the tank height, as physically expected.
Are you saying that you are making ##\epsilon## vary spatially from top to bottom? If so, shouldn't the diffusion term be ##\frac{\partial}{\partial x}\left(\alpha \epsilon\frac{\partial T}{\partial x}\right)##?

If I take your equation and add ##\frac{\partial (\alpha \epsilon)}{dx}\frac{\partial T}{\partial x}## to both sides, I get: $$\frac{\partial T}{\partial t}+\left[V+\frac{\partial (\alpha \epsilon)}{dx}\right]\frac{dT}{dx}=\frac{\partial}{\partial x}\left(\alpha \epsilon\frac{\partial T}{\partial x}\right)$$

HumanistEngineer
Chestermiller said:
Are you saying that you are making ##\epsilon## vary spatially from top to bottom? If so, shouldn't the diffusion term be ##\frac{\partial}{\partial x}\left(\alpha \epsilon\frac{\partial T}{\partial x}\right)##?

If I take your equation and add ##\frac{\partial (\alpha \epsilon)}{dx}\frac{\partial T}{\partial x}## to both sides, I get: $$\frac{\partial T}{\partial t}+\left[V+\frac{\partial (\alpha \epsilon)}{dx}\right]\frac{dT}{dx}=\frac{\partial}{\partial x}\left(\alpha \epsilon\frac{\partial T}{\partial x}\right)$$

That's right. I shall update the solution for the diffusion term. Thanks a lot, Chestermiller.

There is a much better way of setting up the difference equation for this problem to get much better accuracy, by eliminating numerical dispersion in the solution. In an advection-diffusion problem like this, numerical dispersion (false dispersion due to the numerical scheme) is always an issue.

The CN version of this scheme is as follows:
$$\frac{T(x+\Delta x, t+\Delta t)-T(t,x)}{\Delta t}=\frac{1}{2}\left[\frac{\partial }{\partial x}\left(\alpha \epsilon\frac{\partial T}{\partial x}\right)\right]_{(x-\Delta x),t}+\frac{1}{2}\left[\frac{\partial }{\partial x}\left(\alpha \epsilon\frac{\partial T}{\partial x}\right)\right]_{x,(t+\Delta t)}$$with, specifically, ##\Delta t=\Delta x/V##

## 1. What is the 1D convection-diffusion equation and how does it relate to inlet mixing effect?

The 1D convection-diffusion equation is a mathematical model used to describe the transport of a substance in a fluid through a one-dimensional space. It combines the effects of convection (movement of the fluid) and diffusion (random movement of the substance) to determine the overall transport. Inlet mixing effect refers to the impact of introducing a substance at the inlet of the system on the overall transport.

## 2. How is the inlet mixing effect represented in the 1D convection-diffusion equation?

The inlet mixing effect is represented in the 1D convection-diffusion equation through the use of boundary conditions. These conditions specify the concentration of the substance at the inlet, which affects the overall transport and mixing within the system. The strength of the inlet mixing effect can also be adjusted by varying the inlet velocity and diffusion coefficient in the equation.

## 3. What factors influence the inlet mixing effect in the 1D convection-diffusion equation?

The strength of the inlet mixing effect depends on several factors, including the inlet velocity and diffusion coefficient, as mentioned before. Other important factors include the geometry of the system, the concentration and properties of the substance being introduced at the inlet, and any other flow conditions that may be present (such as turbulence).

## 4. How does the inlet mixing effect impact the accuracy of the 1D convection-diffusion equation?

The inlet mixing effect can significantly impact the accuracy of the 1D convection-diffusion equation, especially in cases where the substance being introduced at the inlet has a high concentration or is highly reactive. In these situations, the inlet mixing effect can cause significant changes in the overall transport and mixing within the system, leading to deviations from the predicted results.

## 5. Are there any limitations to using the 1D convection-diffusion equation for studying inlet mixing effect?

While the 1D convection-diffusion equation is a useful tool for studying inlet mixing effect, it does have some limitations. For example, it assumes that the system is one-dimensional and that the substance being transported is well-mixed. In reality, many systems are more complex and may not meet these assumptions, which can affect the accuracy of the equation's predictions.

• Differential Equations
Replies
2
Views
2K
• Chemistry
Replies
13
Views
1K
• Differential Equations
Replies
1
Views
6K
• Programming and Computer Science
Replies
1
Views
1K
• MATLAB, Maple, Mathematica, LaTeX
Replies
41
Views
8K
• Differential Equations
Replies
1
Views
7K
• MATLAB, Maple, Mathematica, LaTeX
Replies
3
Views
3K
• Differential Equations
Replies
5
Views
5K