A Modeling analytical solution of 1D heat equation

I am trying to write code for analytical solution of 1D heat conduction equation in semi-infinite rod. The analytical solution is given by Carslaw and Jaeger 1959 (p305) as
$$
h(x,t) = \Delta H .erfc( \frac{x}{2 \sqrt[]{vt} } )
$$
where x is distance, v is diffusivity (material property) and t is time. I have written following code in MATLAB for to find heat distribution at 40 m from the input with time. The input is a sine wave at x=0.
Code:
clear; close all; clc;
x = 40;       % distance from input to observation point
v = 100;      % diffusivity
p = 300;        % number of steps
Ti = 0;       % Initial head
t =  1:p;    % time steps at which calculation is done, lenth(t) should be
equal to length(To)
delH = sin(linspace(0.1,60,p));

h = zeros(1,p);

for m = 1:p
    func = x/(2*sqrt(v*t(m)));
    h(m) = delH(m) * erfc(func);

end
plot(delH)
hold on;
plot(h)
legend('input','output')
xlabel('time(sec)')
ylabel('temperature')
I am not understanding why the output sinewave keeps its amplitude increasing with time? Should it just mimic the input sine wave with reduced amplitude. The yellow line in
this image shows result from finite difference and I expect result from Carslaw to mimic finite difference result.
 

mjc123

Science Advisor
697
306
I am not understanding why the output sinewave keeps its amplitude increasing with time?
Look at your equation. Is there a clue there? ΔH is your input sine wave; you are multiplying it by a factor that increases with time.
I am wondering whether the equation is the right one for the scenario - is it perhaps for constant ΔH? Someone correct me if I'm wrong. I remember doing something like this some years ago, and I think the equation was of the form
ΔT(x,t) = ΔT0ei(ωt-βx)e-αx
This (or the real part of it) gives at x a sine wave of reduced but constant amplitude, and phase shifted relative to x=0 (which your solution isn't).
 
Look at your equation. Is there a clue there? ΔH is your input sine wave; you are multiplying it by a factor that increases with time.
I am wondering whether the equation is the right one for the scenario - is it perhaps for constant ΔH? Someone correct me if I'm wrong. I remember doing something like this some years ago, and I think the equation was of the form
ΔT(x,t) = ΔT0ei(ωt-βx)e-αx
This (or the real part of it) gives at x a sine wave of reduced but constant amplitude, and phase shifted relative to x=0 (which your solution isn't).
Thank you for your response. You are exactly right that I should have a sine wave at x with constant but reduced amplitude and with a relative phas shift. I am using equation 15 on page 357 of Jan Taler's book. Page 15 is here. Can you please tell me if equation 15 is for this scenario or not?
 

mjc123

Science Advisor
697
306
It's a bit small, and I can't scroll up and down, but it looks to me as if ΔT(x=0) is supposed to be constant. From your table of contents, this section is dealing with "Formula Derivation for Temperature Distribution in a Half-Space with a Step Increase in Surface Temperature", i.e. T(x=0) increases at time t=0 from T0 to T0 + ΔT0 (constant). It looks as if you want p.366 "Formula Derivation for a Quasi-Steady State Temperature Field in a Half-Space when Surface Temperature Changes Periodically".
 
18,822
3,649
You're using the equation incorrectly. In your system, ##\Delta H## the temperature change at the boundary is not constant with time. You need to express the solution to your problem as a convolution integral with appropriate time delays in the erfc term. What would the solution look like if expressed as an integral in this way?
 
Thank you again for your insightful comments.
It's a bit small, and I can't scroll up and down, but it looks to me as if ΔT(x=0) is supposed to be constant. From your table of contents, this section is dealing with "Formula Derivation for Temperature Distribution in a Half-Space with a Step Increase in Surface Temperature", i.e. T(x=0) increases at time t=0 from T0 to T0 + ΔT0 (constant). It looks as if you want p.366 "Formula Derivation for a Quasi-Steady State Temperature Field in a Half-Space when Surface Temperature Changes Periodically".
The problem is that the input function in my real case is not completly sinosoidal but p.366 refers to cases where input is completly periodic. Here is complete Jan Taler's book. Actually in this paper, the author has shown how he used the solution of Carslaw and Jaeger to derive to the final equation. He uses equation equation of Carsla (2a in paper) to come to the conclusion (3a in paper). Equation 2a is same as I wrote above and as given by Carslaw on page 305 (equation 5).
 
You're using the equation incorrectly. In your system, ##\Delta H## the temperature change at the boundary is not constant with time. You need to express the solution to your problem as a convolution integral with appropriate time delays in the erfc term. What would the solution look like if expressed as an integral in this way?
Isn't ##\Delta H## changing with every step in my code? I am sorry I am not very good at math and its been quite a time I am looking at this problem. Can you tell what changes I have to make in my code?
 
18,822
3,649
The equation should be $$h(x,t)=\int_0^t{\frac{d\Delta H(t-\tau)}{d\tau}\left[erfc{\frac{x}{2\sqrt{\nu \tau}}}\right]d\tau}$$
 

Want to reply to this thread?

"Modeling analytical solution of 1D heat equation" You must log in or register to reply here.

Physics Forums Values

We Value Quality
• Topics based on mainstream science
• Proper English grammar and spelling
We Value Civility
• Positive and compassionate attitudes
• Patience while debating
We Value Productivity
• Disciplined to remain on-topic
• Recognition of own weaknesses
• Solo and co-op problem solving
Top