Pde resolution with temperature-dependent material properties

3 次查看(过去 30 天)
Hi All,
I am tryng to solve a pde equation with temperature-dependent material properties and I was wondering if what I implemented was right. I have not managed to check if, for each time-step, Matlab calculates material properties using the current solution of temperature and the polynomial functions prescribed. Could you please give me a hint?
Here the main function is:
function [c,f,s] = pdex1pde(x,t,u,dudx,Tmelt,molarMass)
% Temperature-dependent material props
[ep,Kfr,cpfr,rhofr]=Kcprho(u,molarMass); %this is the function containing the polynomials for the material properties
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
diffr = Kfr/(rhofr*cpfr); % frozen thermal diffusivity
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
c = 1/diffr; % solid phase thermal diffusivity;
f = dudx;
s = 0
end
Boundary conditions:
function [pl,ql,pr,qr] = pdex1bc(xl,ul,xr,ur,t,sigma,Tenv,Boltz,molarMass)
[ep,Kfr,cpfr,rhofr] = Kcprho(ul,molarMass);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
diffr = Kfr/(rhofr*cpfr); % frozen thermal diffusivity
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
pl = (sigma)-(Boltz*ep*(ul^4-Tenv^4));
ql = Kfr;
pr = 0;
qr = 1;
end
Initial condition:
function u0 = pdex1ic(x,Tcl)
u0 = Tcl;
end
Material property file:
function [ep,Kfr,cpfr,rhofr] = Kcprho(u,molarMass)
% Temperature-dependent material props for W
ep=0.3; %constant emissivity
%thermal conductivity in the range [300 K < u < 3695 K]
Kfr = 149.441-45.466*10^(-3)*u+13.193*10^(-6)*u^2-1.484*10^(-9)*u^3+(3.866*10^(6))/u^2; % [W/m K]
%specific heat capacity in the range [300 K < u < 3695 K]
cpfr = (21.868372+8.068661*10^(-3)*u-3.756196*10^(-6)*u^2+1.075862*10^(-9)*u^3+(1.406637*10^(-4))/u^2)/molarMass; % [J/kg K]
%mass density in the range [300 K < u < 3695 K]
rhofr = 1000*(19.25-2.66207*10^(-4)*(u-273.15)-3.0595*10^(-9)*(u-273.15)^2-9.5185*10^(-12)*(u-273.15)^3); % [kg/m^3]
end

回答(0 个)

产品

Community Treasure Hunt

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

Start Hunting!

Translated by