3D geometry issues with specifyCoefficients (PDE Toolbox)

Hello everyone! I am trying to prepare a heat transfer simulation in an adiabatic pipe. There is a specific formula that I would like my students to compare results with by setting coefficients for the generalized form of the equations solved by PDE Toolbox. Please look at the code below - I would really appreciate your help!

I noticed this error:

Incorrect number or types of inputs or outputs for function specifyCoefficients.

When I checked the documentation, it turns out that my nonconstant coefficients are the correct format, but the model I am specifying is a femodel class as opposed to PDEmodel (requested by the first argument). I made the 3D geometry in house, and I was initially associating it to the model by using femodel(AnalysisType, Geometry), which you can see commented out below.

I then tried to add model_adiabatic.Geometry = g3 hoping that would keep model_adiabatic as a PDEmodel and not a femodel, but it still shows as a femodel when I use class(model) and I am getting another error on top.

I have seen examples with specifyCoefficients used on 3D geometries, but they use importGeometry from external files. I do not have that and I was hoping to keep the code simple and transferable.

Could anybody kindly suggest a fix or anything I have missed?

Thank you!

Michela

% Setting Adiabatic Model

R = 0.05; %m
rho = 1000; %kg/m3
cp = 4182; %J/kgC

model_2 = createpde("thermal","steadystate");

% Geometry

pdecirc(0,0,R);
g = decsg(gd,sf,ns);
g2 = fegeometry(g);
g3 = extrude(g2,1);

model_adiabatic.Geometry = g3
figure;

% model_adiabatic = femodel(AnalysisType="thermalSteady", ...
%     Geometry=g3); 

pdegplot(model_adiabatic,EdgeLabels="on",FaceLabels="on")
title('Pipe')

class(model_adiabatic)

model.MaterialProperties=materialProperties(ThermalConductivity=k,MassDensity=rho, SpecificHeat=cp);

% Boundary Conditions

model_adiabatic.FaceLoad(3) = faceLoad(Heat = 0); % adiabatic wall
model_adiabatic.FaceLoad(4) = faceLoad(Heat = 0); % adiabatic wall
model_adiabatic.FaceLoad(5) = faceLoad(Heat = 0); % adiabatic wall
model_adiabatic.FaceLoad(6) = faceLoad(Heat = 0); % adiabatic wall

T1 = 40;
model_adiabatic.FaceBC(1) = faceBC(Temperature = T1); % entrance temperature, z = 0

% Boundary Conditions

model_adiabatic.FaceLoad(3) = faceLoad(Heat = 0); % adiabatic wall
model_adiabatic.FaceLoad(4) = faceLoad(Heat = 0); % adiabatic wall
model_adiabatic.FaceLoad(5) = faceLoad(Heat = 0); % adiabatic wall
model_adiabatic.FaceLoad(6) = faceLoad(Heat = 0); % adiabatic wall

T1 = 40;
model_adiabatic.FaceBC(1) = faceBC(Temperature = T1); % entrance temperature, z = 0

Answer 1

[Ravi Kumar]

Try this. This should not throw any error. Note that I assumed k = 1 and created 2-D geometry by directly calling decsg.

% Setting Adiabatic Model
R = 0.05; %m
k = 1; % CHANGE THIS
rho = 1000; %kg/m3
cp = 4182; %J/kgC

% Geometry
C = [1,0,0,R]';
g = decsg(C);
g2 = fegeometry(g);
g3 = extrude(g2,1);
%model_adiabatic.Geometry = g3
figure;
model_adiabatic = femodel(AnalysisType="thermalSteady", ...
     Geometry=g3); 
pdegplot(model_adiabatic,EdgeLabels="on",FaceLabels="on")
title('Pipe')

class(model_adiabatic)

ans = 'femodel'

model.MaterialProperties=materialProperties(ThermalConductivity=k,MassDensity=rho, SpecificHeat=cp);
% Boundary Conditions
model_adiabatic.FaceLoad(3) = faceLoad(Heat = 0); % adiabatic wall
model_adiabatic.FaceLoad(4) = faceLoad(Heat = 0); % adiabatic wall
model_adiabatic.FaceLoad(5) = faceLoad(Heat = 0); % adiabatic wall
model_adiabatic.FaceLoad(6) = faceLoad(Heat = 0); % adiabatic wall
T1 = 40;
model_adiabatic.FaceBC(1) = faceBC(Temperature = T1); % entrance temperature, z = 0
% Boundary Conditions
model_adiabatic.FaceLoad(3) = faceLoad(Heat = 0); % adiabatic wall
model_adiabatic.FaceLoad(4) = faceLoad(Heat = 0); % adiabatic wall
model_adiabatic.FaceLoad(5) = faceLoad(Heat = 0); % adiabatic wall
model_adiabatic.FaceLoad(6) = faceLoad(Heat = 0); % adiabatic wall
T1 = 40;
model_adiabatic.FaceBC(1) = faceBC(Temperature = T1); % entrance temperature, z = 0

[Michela]

Hi Ravi, I really appreciate your help with geometry! Thank you very much.

However, maybe my title is a bit misleading: I am trying to set coefficients for PDE Toolbox, and that is the part (the latter part of my code) that I cannot quite figure out! Do you have any idea how to change that?

Thank you and have a lovely day!

[Ravi Kumar]

Hi Michela,

The specifyCoefficients function is for PDEModel, where you specify everything in equations format. The femodel you created is for thermal analysuis. You don't need to use specifyCoefficients function with femodel. Whatever you want to do with specifyCoefficients can be done by specifying appropriate material properites:

model.MaterialProperties=materialProperties(ThermalConductivity=k,MassDensity=rho, SpecificHeat=cp);

Note that you can specify ThermalConductivity as a function instead of constant. If this still does not answer your question, can you post what you tried and the error you are getting.

Ravi

[Michela]

Hello Ravi,

Before I test it out, I hope it is okay if I ask another question. Is there a way to specify a 3D geometry for a PDE model? I specifically wanted to try using specifyCoefficients.

THank you very much!

Michela

[Ravi Kumar]

Hi Michela,

Sure, you can use 3-D geometry with a PDEModel, the output of createpde. In this case you should use specifyCoefficients, applyBoundaryConditions, setInitialConditions to define your PDE and then call solvepde on the PDEModel object.

However, if you are solving a thermal PDE with typical material and BCs, you should be able to set it up using femodel. The PDEModel interface is useful when the PDEs don't fit into well pacakged analysis types that femodel doesn't support. You will be dealing with the equations directly if you use PDEModel. So, I suggest you try femodel if you are seting up a standard problem.

Regards,

Ravi