熱モデルのGeome​tryの設定にimp​ortGeometr​y関数とmultis​phere関数とを使​うのでシミュレーショ​ンの挙動が異なるのは​なぜですか?

3 views (last 30 days)
私はPartial Differential Equation Toolboxを使って二層の球からなるジオメトリを持つ熱モデルのシミュレーションをしています. 内側の球を内部熱源として機能させています.その中で, importGeometry関数でGeometryを設定した熱モデルとmultisphere関数でGeometryを設定した熱モデルとでシミュレーションの挙動が異なることを確認しました. 内側の球の境界条件が変化しているように思えます. multisphere関数で設定した方は, 内側の球の熱をそのまま伝えているのに対し, importGeometry関数で設定した方は, 内側の球の境界で断熱が起きているような挙動になります. どうしてこのようなことが起きるのでしょうか.
%Using multisphere function
%Create thermalmodel
thermalmodel = createpde('thermal','transient');
gm = multisphere([3.15*10^-3,1*10^-2]);
thermalmodel.Geometry = gm;
%Visualize geometry
figure
subplot(1,2,1)
pdegplot(thermalmodel,'FaceAlpha',0.25,'CellLabel','on')
title('Geometry with Cell Labels')
subplot(1,2,2)
pdegplot(thermalmodel,'FaceAlpha',0.25,'FaceLabel','on')
title('Geometry with Face Labels')
%Set physical parameter
generateMesh(thermalmodel,'Hmax',1*10^-3);
thermalProperties(thermalmodel,'Cell',1,'ThermalConductivity',0.778,'MassDensity',1.66*10^3,'SpecificHeat',2.54*10^3);
thermalProperties(thermalmodel,'Cell',2,'ThermalConductivity',0.642,'MassDensity',1*10^3,'SpecificHeat',3.72*10^3);
internalHeatSource(thermalmodel,6.15*10^6,'cell',1);
thermalBC(thermalmodel,'Face',2,'ConvectionCoefficient',25,'AmbientTemperature',25);
thermalIC(thermalmodel,25);
%Calculate PDE
tlist=0:1:300;
result = solve(thermalmodel,tlist);
Tmin=0;
Tmax=max(result.Temperature(:,end));
%Visualize by contour
[YG,ZG] = meshgrid(linspace(-1*10^-2,1*10^-2,100),linspace(-1*10^-2,1*10^-2,100));
XG = zeros(size(YG));
tIndex = [6,51,101,151,201,301];
varNames = {'Time_index','Time_step'};
index_step = table(tIndex.',tlist(tIndex).','VariableNames',varNames);
disp(index_step);
TG = interpolateTemperature(result,XG,YG,ZG,tIndex);
t = linspace(0,2*pi);
ylayer1 = 3.15*10^-3*cos(t); zlayer1 = 3.15*10^-3*sin(t);
ylayer2 = 1*10^-2*cos(t); zlayer2 = 1*10^-2*sin(t);
figure('Position',[10,10,1000,550]);
for i = 1:numel(tIndex)
subplot(2,3,i)
contour(YG,ZG,reshape(TG(:,i),size(YG)),'ShowText','on')
colormap turbo
c=colorbar;
c.Label.String = 'Temperature(℃)';
title(['Temperature at Time ' num2str(tlist(tIndex(i)))]);
hold on
caxis([Tmin,Tmax])
axis equal
plot(ylayer1,zlayer1,'k','LineWidth',1.5)
plot(ylayer2,zlayer2,'k','LineWidth',1.5)
end
↓Result
%Using importGeometry function
%Create thermalmodel
thermalmodel = createpde('thermal','transient');
gm = importGeometry('ball_simu.stl');
thermalmodel.Geometry = gm;
%Visualize geometry
figure
subplot(1,2,1)
pdegplot(thermalmodel,'FaceAlpha',0.25,'CellLabel','on')
title('Geometry with Cell Labels')
subplot(1,2,2)
pdegplot(thermalmodel,'FaceAlpha',0.25,'FaceLabel','on')
title('Geometry with Face Labels')
%Set physical parameter
generateMesh(thermalmodel,'Hmax',1*10^-3);
thermalProperties(thermalmodel,'Cell',2,'ThermalConductivity',0.778,'MassDensity',1.66*10^3,'SpecificHeat',2.54*10^3);
thermalProperties(thermalmodel,'Cell',1,'ThermalConductivity',0.642,'MassDensity',1*10^3,'SpecificHeat',3.72*10^3);
internalHeatSource(thermalmodel,6.15*10^6,'cell',2);
thermalBC(thermalmodel,'Face',1,'ConvectionCoefficient',25,'AmbientTemperature',25);
thermalIC(thermalmodel,25);
%Calculate PDE
tlist=0:1:300;
result = solve(thermalmodel,tlist);
Tmin=0;
Tmax=max(result.Temperature(:,end));
%Visualize by contour
[YG,ZG] = meshgrid(linspace(-1*10^-2,1*10^-2,100),linspace(-1*10^-2,1*10^-2,100));
XG = zeros(size(YG));
tIndex = [6,51,101,151,201,301];
varNames = {'Time_index','Time_step'};
index_step = table(tIndex.',tlist(tIndex).','VariableNames',varNames);
disp(index_step);
TG = interpolateTemperature(result,XG,YG,ZG,tIndex);
t = linspace(0,2*pi);
ylayer1 = 3.15*10^-3*cos(t); zlayer1 = 3.15*10^-3*sin(t);
ylayer2 = 1*10^-2*cos(t); zlayer2 = 1*10^-2*sin(t);
figure('Position',[10,10,1000,550]);
for i = 1:numel(tIndex)
subplot(2,3,i)
contour(YG,ZG,reshape(TG(:,i),size(YG)),'ShowText','on')
colormap turbo
c=colorbar;
c.Label.String = 'Temperature(℃)';
title(['Temperature at Time ' num2str(tlist(tIndex(i)))]);
hold on
caxis([Tmin,Tmax])
axis equal
plot(ylayer1,zlayer1,'k','LineWidth',1.5)
plot(ylayer2,zlayer2,'k','LineWidth',1.5)
end
↓Result

Accepted Answer

Hiroyuki Hishida
Hiroyuki Hishida on 5 Sep 2022
佐藤様、
generateMeshで生成される形状が異なるためです。例えばthermalproperty設定の前あたりに、以下を差し込んで実行してみてください。すると、multiSphereとimportのときで、表示される形状(点群)が異なることが確認できると思います。
myMesh = generateMesh(thermalmodel,'Hmax',0.01);
figure
pdemesh(myMesh,'FaceAlpha',0)
title('generateMeshで生成される形状')
佐藤様のスクリプトでもそうですが、STLをimportしてもそれをそのまま解析に渡すことができず、generatemeshします。このgenerateMeshは、与えられた形状の主に頂点情報から有限要素を生成します。そのため、与える形状情報はstlでなくても、csg表現でも、alpha shapeにも対応します。
この与える形状情報ですが、multisphereの場合、球の中に球があるカタチにはなっておらず、あくまでモノとしての形状は一つの球です。一方、佐藤様のスクリプトにおいては元のSTLにおいても球の中に球を作成されています。ここが違いとしてあります。
シンプルな例題としては以下が参考になるかと思います。
では単一の球をSTLで与えて内部に球をCELLで追加したい場合どうするかですが、例えば以下はいかがでしょうか。
どうぞよろしくお願いします。
菱田
  2 Comments
和奈 佐藤
和奈 佐藤 on 6 Sep 2022
参考にさせて頂いたところ, うまくいきました.
ありがとうございます!
Hiroyuki Hishida
Hiroyuki Hishida on 6 Sep 2022
うまくいって良かったです。
菱田

Sign in to comment.

More Answers (0)

Community Treasure Hunt

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

Start Hunting!