Which boundary condition is correct or what will be the correct alternative?

7 views (last 30 days)
% Constants
Ly = 0.001; % Thickness of the sheet (meters)
T = 15; % Total simulation time (seconds)
Ny = 10; % Number of spatial grid points
Nt = 300; % Number of time steps
k1 = 0.15; % Thermal conductivity of PVC W/m-K
rho = 1250; % density in kg/m3
cp = 1350 ; % specific heat in J/kg-K
alpha = k1/(rho*cp); % Thermal diffusivity (m^2/s)
T_ambient = 25; % Ambient temperature (°C)
velocity_x = 0.1667; % Velocity in x-direction (m/s)
%contact with hot roller
k2 = 50; %Thermal conductivity of steel roller W/m-K
h_roller = 500; %convective heat transfer coefficient of roller W/m2-K
d_roller = 0.4; %diameter of roller in m
contact_angle = 180; %in °
contact_length = pi*d_roller*contact_angle/360; %contact length with roller
T_roller = 100; % Roller temperature
%Contact in air
distance_x = 0.5; % Distance traveled in x-direction (meters)
h_air = 12; % Convective heat transfer coefficient of air W/m2-K
%Contact in IR field
distance_ir = 0.3 ; % in m
radiative_flux = 50e3; %W/m2
absorption = 50; % in %
performance = 80; % in %
Net_radiative_intensity = radiative_flux*(absorption/100)*(performance/100);
%Contact with air again
distance_v = 1 ; % in m
%Naming the distance variables
y1 = contact_length;
y2 = contact_length + distance_x;
y3 = contact_length + distance_x+ distance_ir;
y4 = contact_length + distance_x+ distance_ir+distance_v;
% Discretization
dy = Ly / (Ny - 1);
dt = T / Nt;
y = linspace(0, Ly, Ny);
t = linspace(0, T, Nt);
% Initial temperature distribution
initial_temperature = 25;
distance_travelled = zeros(1,Nt); %initial distance travelled
% Initialize temperature matrix
u = zeros(Ny, Nt);
u(:, 1) = initial_temperature;
% Forward euler method
r = alpha * dt / (dy^2); % as 0 < r < 0.5
% Time-stepping loop (explicit method)
for n = 1:Nt - 1
% Update position in x-direction
distance_travelled(n+1) = distance_travelled(n) + velocity_x * dt;
% Check if the total distance is within the contact length
if distance_travelled(n) <= y1
for i = 2:Ny - 1
u(i, n + 1) = u(i, n) + r * (u(i + 1, n) - 2 * u(i, n) + u(i - 1, n));
end
% Apply boundary conditions for contact with the roller
u(1, n + 1) = ((h_roller * T_roller * dy) + (k1 * u(2, n + 1))) / (k1 + h_roller * dy);
u(Ny, n + 1) = ((h_air * T_ambient * dy) + (k1 * u(Ny - 1, n + 1))) / (k1 + h_air * dy);
% For the next part in contact with air
elseif distance_travelled(n) <= y2
for i = 2:Ny - 1
u(i, n + 1) = u(i, n) + r * (u(i + 1, n) - 2 * u(i, n) + u(i - 1, n));
end
% Apply boundary conditions for contact with air
u(1, n + 1) = ((h_air * T_ambient * dy) + (k1 * u(2, n + 1))) / (k1 + h_air * dy);
u(Ny, n + 1) = ((h_air * T_ambient * dy) + (k1 * u(Ny - 1, n + 1))) / (k1 + h_air * dy);
%In contact with IR field radiation
elseif distance_travelled(n) <= y3
for i = 2:Ny - 1
u(i, n + 1) = u(i, n) + r * (u(i + 1, n) - 2 * u(i, n) + u(i - 1, n));
end
% Apply boundary conditions for contact with IR field
u(1, n + 1) = ((h_air * T_ambient * dy) + (k1 * u(2, n + 1))) / (k1 + h_air * dy);
u(Ny, n + 1) = u(Ny - 1, n + 1) + ((Net_radiative_intensity * dy) / k1);
% Distance in contact with air again
else
for i = 2:Ny - 1
u(i, n + 1) = u(i, n) + r * (u(i + 1, n) - 2 * u(i, n) + u(i - 1, n));
end
% Apply boundary conditions for contact with air
u(1, n + 1) = ((h_air * T_ambient * dy) + (k1 * u(2, n + 1))) / (k1 + h_air * dy);
u(Ny, n + 1) = ((h_air * T_ambient * dy) + (k1 * u(Ny - 1, n + 1))) / (k1 + h_air * dy);
end
% Break the loop when the total distance is reached
if distance_travelled(n + 1) >= y4
break;
end
end
% Plot the temperature distribution over time
figure;
for n = 1:min(Nt,numel(distance_travelled))
plot(y, u(:, n));
xlabel('Thickness (m)');
ylabel('Temperature (°C)');
title(['Time = ', num2str(t(n)), ' seconds']);
axis([0 Ly 20 150]); % Adjust the axis limits if needed
pause(0.1); % Pause to create animation effect
% Break the plotting loop when the loop breaks
if distance_travelled(n) >= y4
break;
end
end
which boundary condition is correct? The PVC material's 1st element is in contact with the roller whose temperature is maintained at 100°C so ideally the mode of heat transfer should be conduction. In that case, either dirichlet boundary condition or neumann boundary condition should be implemented
If i apply dirichlet boundary condition u(1,n+1) = T_roller so it is fixing the temperature to 100°C until it is contact with roller which will not be the case . In reality, the temperature of the first element will rise more and try to reach 100°C as compared to other elements but won't be 100°C
For Neumann boundary condition, I need the information of temperature gradient or flux but unfortunately, I don't have neither of that.
I am not sure the boundary condition which I have applied here is correct or not? Can someone please clarify this or help me in this issue?

Accepted Answer

Torsten
Torsten on 26 Sep 2023
Edited: Torsten on 26 Sep 2023
If you can estimate the heat transfer coefficient k between PVC material and roller, you could use
-lambda_PVC * dT_PVC/dn = k * (T_PVD - T_roller)
where d/dn is the spatial derivative in the direction of the outer normal. Look up "Robin boundary condition".
If k is large, this means that the PVC temperature follows the roller temperature very fast. In this case, it is justified to use the Dirichlet boundary condition for T_PVC.
But if the two materials are in direct contact, they should share the same temperature at this location, shouldn't they ?
  9 Comments
Joydeep Bagchi
Joydeep Bagchi on 27 Sep 2023
Once agian thanks Torsten for the reply.
There is no transition zone used in the last code I posted.
Ok, but I need a continuous transition not like the temperature distribution at each time step rather a continuous temperature distribution through the whole time span.
Which question ? You mean the boundary condition you have to use on both ends of the web ?
The link of the pdf which I have posted has the answer in it so I have implemented that condition which is right.
Your boundary condition terms use values for u at time t(n+1) on the right-hand sides. Thus you have a mixture of explicit and implicit Euler for which I'm not able to predict whether it will work out correctly.
I need to give the term for gradient thatswhy I have written like that which is working fine for me.
Thanks a lot for your effort and time.
Torsten
Torsten on 27 Sep 2023
Ok, but I need a continuous transition not like the temperature distribution at each time step rather a continuous temperature distribution through the whole time span.
The temperature distribution over the web will be continuous if you choose more curves to plot.
The link of the pdf which I have posted has the answer in it so I have implemented that condition which is right.
I cannot find a boundary condition for the roller on page 31 of the document.
I need to give the term for gradient thatswhy I have written like that which is working fine for me.
It's your code, but your advancement in time is not Explicit Euler.

Sign in to comment.

More Answers (0)

Products


Release

R2023a

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!