Create time dependent internal heat source for heat transfer problem

Hello,

I am writing a heat source function that is supposed to be time dependent, but it is not produccing the right results. I wrote a internal heat source function "HeatSourceFunc" that will be used as follows:

internalHeatSource(model,@HeatSourceFunc);

I have written two different heat source functions that basically do the same thing, except that the first case is not time-dependent, and the second case is time-dependent and is supposed to start at time t = 5 seconds.

Here is case 1 (not time dependent):

function Q = HeatSourceFunc(location,~) % case 1
    xx = location.x;
    yy = location.y;
    
    Q = zeros(size(xx));         
    idx = (xx-0.006).^2 + (yy-0.004).^2 <= 0.002^2; % circular heat source
    Q(idx) = 50;
end

Here is case 2 (time dependent):

function Q = HeatSourceFunc(location,state) % case 2
    xx = location.x;
    yy = location.y;
    tt = state.time;
    
    Q = zeros(size(xx));
    
    if tt < 5     % no heat sources before t = 5 seconds
        return;
    else          % circular heat source starting at t = 5 seconds
        idx = (xx-0.006).^2 + (yy-0.004).^2 <= 0.002^2;
        Q(idx) = 20 - tt;
    end
end

I ran the simulation from t = 0 seconds to t = 100 seconds as follows:

time = 0:100; % times at which to sample solution
sol = solve(model,time); % solve problem

Then I obtained the following graphs. Why am I not getting a similar picture in case 2 as in case 1? 100 seconds is plenty of time for a temperature change to occur. Any help would be appreciated, thank you.

Answer 1

[Alan Weiss]

I am not sure, but I believe that the answer might be that the solver checks if anything is going on at small times, and then decides that the problem is not dependent on any variables. This is a fallacious decision, of course, but that might be the problem.

The solution (if I am correct) is to solve two problems. Solve one for times . Then take the solution at time 5 as the initial condition for the problem for . This workaround should enable you to obtain correct answers.

Also, there are good reasons to solve problems this way even in the absence of this bug. When the coefficients of your problem depend nonsmoothly on time, the ODE solver can skip over the time of the nonsmooth point. You should help the solver by explicitly restarting it at such points.

Good luck,

Alan Weiss

MATLAB mathematical toolbox documentation

[Ekaterina Startseva]

For this problem you can use analytic approximation of Heaviside function. If internal heat source is on for times and off for times Heaviside function approximation is , where a larger k corresponds to a steeper rise of the function at the point . Given the required width of the transition region of the Heaviside function , the value of k can be estimated as .

So for example heat sourse is . In this case, the discontinuous function is replaced by a smooth one and you can use

q = @(location,state)2000*location.x*(0.5+0.5*tanh(k*(state.time-5)));

or use a separate function:

internalHeatSource(thermalModel,@heatsource, 'Face',1)
function q = heatsource(location,state)
t  = state.time;
x = location.x;
k = 100;
q = 2000*x.*(0.5-0.5*tanh(k*(t-5)));
end