PDE Toolbox: Can I use a time dependent ambient temperature?

I'm trying to simulate 1D semi-infinite conduction with a convective boundary condition on one surface. I'd like to apply a time varying ambient temperature to the convective boundary in the form of of a square wave. When running with a constant ambient temperature my code works perfectly and gives the expected solution, however when I use a function to apply the time varying condition I get

Error using pde.ThermalBC (line 83)

The value of 'AmbientTemperature' is invalid. Must be specified as a numeric value.

Is there any way around this constraint? It seems a bit strange that other boundary conditions can be set as a function, but ambient temperature can't. If it isn't possible to vary the ambient temperature then I think I could get around it by setting a heat flux boundary condition at every time step as and vary the ambient temperature this way, but this isn't ideal. This is my code

% Material properties
k = 0.36; % Thermal conductivity, W/(m*K)
rho = 1610; % Density, kg/m^3
cp = 958; % Specific heat, J/(kg*K)

% Thermal conditions
htc = 500; % Heat transfer coefficient, W/(m^2*K)
T_i = 300; % Initial material temperature, K
T_aw = 325; % Adiabatic wall temperature, K

% Set a square wave ambient temperature
frequency = 0.25; % Heating frequency, Hz
dT = @(location,state)...
    T_i+(T_aw-T_i)*0.5*(square(state.time*2*pi*frequency)+1);

% Set time steps
tfinal = 300; % End time, s
tlist = 0:1:tfinal; % Time array

% Set the maximum cell size
Hmax = 0.001;

% Rectangular domain
half_height = 2*Hmax; % Height either side of x axis, m
width = 0.1; % Depth into x axis, m

% Set-up the PDE solver
thermalModelT = createpde('thermal','transient');

% Create the domain geometry
g = decsg([3 4 0 0 width width -half_height half_height half_height -half_height]');
geometryFromEdges(thermalModelT,g);

% Generate the mesh for the model
generateMesh(thermalModelT,'Hmax',Hmax);

% Set material properties
thermalProperties(thermalModelT,'ThermalConductivity',k,...
                                'MassDensity',rho,...
                                'SpecificHeat',cp);

% Set thermal boundary conditions

% Convective on the surface
thermalBC(thermalModelT,'Edge',1,...
                        'ConvectionCoefficient',htc,...
                        'AmbientTemperature',dT);

% Heat flux BC on other 3 edges
thermalBC(thermalModelT,'Edge',2,'HeatFlux',0);
thermalBC(thermalModelT,'Edge',3,'HeatFlux',0);
thermalBC(thermalModelT,'Edge',4,'HeatFlux',0);

% Isothermal initial condition
thermalIC(thermalModelT,T_i);

% Solve
result = solve(thermalModelT,tlist);

Thanks for any help you can give!

[Mattia Battaglia]

I have the same problem. The ability to use time dependent ambient temperature is mandatory for my task and for the time being keeps me from using the PDE toolbox. Is there any chance this could be implemented?

[Ravi Kumar]

You can use 'HeatFlux' parameter to code the function which could accommodate both convection coefficient and ambient temperature as non-constant quantities.

[Philippe]

hello,

not sure, because heat flux function is not dependent of temperature variable on face , as convective BC is ??

Philippe

Answer 1

[Bradley]

If your ambient temperature is variable using a square wave, you could easily break up the call to solve based on the time increments of your square wave. So, you'll end up with several smaller time increments which you just stitch together at the end.

To get around the issue of carrying the final values from one time interval to the next, you can re-assign the initial conditions to the next from the final result of the last by using the normal thermalIC function:

thermalIC(model,previous_results)

I've done this in a model looking at the thermal response of a battery cell going from discharge to charge at certain time intervals. It's worked well, with little fuss.

[William Davis]

This is a great suggestion!