在 PDE 工具箱中使用时间相关(热)源

Using time dependent (heat) source in PDE Toolbox

我的目标是在解决传热问题时应用瞬态热源。

瞬态传导传热的偏微分方程为:

更多信息可以在这里找到:Solving a Heat Transfer Problem With Temperature-Dependent Properties

在我的例子中,所有参数都是常量,除了源项 f 需要随时间变化。

我在这里遵循示例代码:Nonlinear Heat Transfer In a Thin Plate 它提供了一种解决瞬态问题的方法,并且我能够绘制每个时间点的热数据。

将它应用于我的案例时的问题是,在示例中,源是一个常数值,贯穿整个区域和整个时间,并且与辐射和对流有关(在我的案例中它们应该全为零) ,但我需要提供一个随时间变化的源(通过时变电流产生的焦耳热)。来源可能具有以下格式之一:

  1. 解析:正值,例如 1 W/m^2,在时间 window 内,例如 0< t< 1 ns,否则为 0。
  2. 数值:数据由 1xN 向量提供,其中 N 是时间点的数量。

并且来源被限制在某个区域,例如。 0< x <1 毫米且 0< y<1 毫米。

我看到过类似的问题但没有回答:How to use a variable coefficient in PDE Toolbox to solve a parabolic equation (Matlab)

有没有办法用 PDE 工具箱实现这个?从头开始写代码会很复杂....

您可以使用 FEATool FEM Matlab Toolbox 非常轻松地定义和解决具有时间相关和非线性 PDE 系数的问题,如下面的 m 脚本代码片段所示。请注意,热源(汇)项 f 被缩放为 f*(t>2500) 这意味着它只会在 t=2500 之后激活(因为 switch 表达式的计算结果为 0 如果为假或 1 如果为真)。

% Coefficents and problem definition from https://www.mathworks.com/help/pde/examples/nonlinear-heat-transfer-in-a-thin-plate.html
k=400;rho=8960;specificHeat=386;thick=.01;stefanBoltz=5.670373e-8;hCoeff=1;ta=300;emiss=.5;

% Set up 2D fea struct with geometry and grid.
fea.sdim = {'x' 'y'};
fea.geom = { gobj_rectangle( 0, 1, 0, 1 ) };
fea.grid = rectgrid( 10 );

% Add heat transfer physics mode.
fea = addphys( fea, @heattransfer );

fea.phys.ht.eqn.coef{1,end}{1} = rho*thick;      % Density eqn coefficient.
fea.phys.ht.eqn.coef{2,end}{1} = specificHeat;   % C_p eqn coefficient.
fea.phys.ht.eqn.coef{3,end}{1} = k*thick;        % Thermal condictivity.
f = sprintf( '%g*( %g - T ) + %g*( %g^4 - T^4 )', ...
             2*hCoeff, ta, 2*emiss*stefanBoltz, ta );
fea.phys.ht.eqn.coef{6,end}{1} = ['(',f,')*(t>2500)'];   % Heat source term.

fea.phys.ht.bdr.sel(1) = 1;  % Set prescribed temperature for boundary 1.
fea.phys.ht.bdr.coef{1,end}{1} = 1000;  
fea.phys.ht.bdr.sel(2:4) = 3;   % Isolation BCs for boundaries 2-4.

% Check, parse, and solve fea problem.
fea = parsephys( fea );
fea = parseprob( fea );
[fea.sol.u,t] = solvetime( fea, 'tstep', 50, 'tmax', 5000, 'init', {ta} );

% Postprocessing and visualization.
for i=1:size(fea.sol.u,2)
  T_top(i) = evalexpr( 'T', [.5;1-sqrt(eps)], fea, i );
end

subplot(1,2,1)
postplot( fea, 'surfexpr', 'T', 'title', 'T @ t=5000' )
subplot(1,2,2)
plot( t, T_top, 'r-' )
xlabel( 't' )
ylabel( 'T(0.5,1) @ t=5000' )
grid on

在此处的解决方案中,您可以看到由于热量从下边界向上扩散,顶部边缘的温度呈线性上升,直到 t=2500,此时散热器被激活。

Matlab FEATool Nonlinear Time dependent Heat Transfer Solution

关于你的第二点数字源项,在这种情况下你可以创建并调用你自己的外部函数来制表和内插你的数据,在这种情况下它看起来像

fea.phys.ht.eqn.coef{6,end}{1} = 'my_fun( t )';

你有一个 Matlab 函数 my_fun.m 在表单的 Matlab 路径上可以访问的地方

function [ val ] = my_fun( t )
data  = [7 42 -100  0.1];   % Example data.
times = [0 10   99 5000];   % Time points.
val = interp1( data, times, t );   % Interpolate data.

最后,虽然模型是使用 m-script Matlab 代码定义的,以便在 Whosebug 上共享,但您可以轻松地使用 Matlab FEA GUI,如果需要,甚至可以将您的 GUI 模型导出为 m-script 代码。