How to create a geometry consisting of a prism, with an inner inclined plane, for a diffusion pde problem?

Question

Gabriel Madirolas 2018-9-4

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/417429-how-to-create-a-geometry-consisting-of-a-prism-with-an-inner-inclined-plane-for-a-diffusion-pde-pr

回答： Gabriel Madirolas 2019-9-2

For example, the nodes could be placed at (0, 0 ,0), (1, 0 ,0), (1, 1 ,0), (0, 1 ,0), (0, 0 ,1), (1, 0 ,1), (1, 1 ,1), (0, 1 ,1). The plane could go from the upper left edge to the lower right edge, so its nodes would be (0, 0 ,1), (0, 1 ,0), (1, 0, 0), (1, 1, 0).

Everytime I tried

model = createpde(); geometryFromMesh(model,nodes,elements); or importGeometry(model,'PrismAndPlane.stl');

I have gotten all kinds of error messages (no matter how I generated the nodes and elements), including

'The input mesh does not form a closed volume. This may be due to missing faces or gaps in the input.' and 'Failed to create geometry. The stl file is invalid, more than two facets share an edge.'

I cannot include specific code, because I have literally tried a dozen different solutions (alpha shapes, geometry from meshes, create and import stl files, etc). The only one that worked provided a double inner plane, which is not ok to solve a diffusion pde.

So take it as a quiz. Glory for the one who solves it!

Thanks in advance,

Gabriel

5 个评论
显示 3更早的评论隐藏 3更早的评论

Gabriel Madirolas 2018-9-5

在 MATLAB Online 中打开

Thanks Alan,

The code below, for example, generates two contiguous regions:

% Define the two volumes (wedges) that will form the geometry
xUpperWedge = [0 0 1 1 1 1]';
yUpperWedge = [0 1 0 1 0 1]';
zUpperWedge = [1 1 1 1 0 0]';
xLowerWedge = [0 0 0 0 1 1]';
yLowerWedge = [0 1 0 1 0 1]';
zLowerWedge = [1 1 0 0 0 0]';
% Create an alpha shape for each wedge,
% and obtain a surface mesh from each
upperWedge = alphaShape(xUpperWedge,yUpperWedge,zUpperWedge);
[upperElements,upperNodes] = boundaryFacets(upperWedge);
lowerWedge = alphaShape(xLowerWedge,yLowerWedge,zLowerWedge);
[lowerElements,lowerNodes] = boundaryFacets(lowerWedge);
lowerElements = lowerElements +6; % this is to concatenate with upper 
% Put the data in the correct shape for geometryFromMesh.
upperElements = upperElements';
upperNodes = upperNodes';
lowerElements = lowerElements';
lowerNodes = lowerNodes';
elements = [upperElements lowerElements];
nodes = [upperNodes lowerNodes];
% Create a PDE model and import the surface mesh.
plate = createpde();
geometryFromMesh(plate,nodes,elements);
% Plot for visualzation purposes
figure
plot(upperWedge,'facecolor','g','facealpha',0.3)
hold on
plot(lowerWedge,'facecolor','y','facealpha',0.3)
% View the geometry and face numbers.
figure
pdegplot(plate,'FaceLabels','on','FaceAlpha',0.5)

I want to study diffusion from one to the other. But the problem in this case is that I get a double inner face (labeled F3 and F9, at least when I run it in my computer). I would like the boundary between the two regions to be a single face, because otherwise the pde solver won't consider both regions to be connected (as far as I know).

I found an example of a multiple region in the documentation of geometryFromMesh:

load MultidomainMesh3D
model = createpde;
geometryFromMesh(model,nodes,elements,ElementIdToRegionId);
pdegplot(model,'CellLabels','on','facelabels','on','facealpha',0.5)

This looks similar, but I would want the inner plane to be inclined. Unfortunately, it is not specified how the MultidomainMesh3D geometry was created, so I cannot custom my own geometry inspired by this one.

Again, I have tried many things, including create and import an stl file, but nothing worked.

Thanks in advance,

Gabriel Madirolas

Simone Galante 2018-10-11

Immagine1.jpg

Hello Gabriel,

Have you been able to solve your problem? Have you found the method to merge your inner face F3 and F9? I have found a similar problem with the edges about two cyilindical solids that I show you in the attached immage file. Thanks in advance, Simone

Gabriel Madirolas 2018-10-15

Hi Simone,

Unfortunately, I couldn't manage a solution to my problem. I think Matlab's pde is still not working smoothly enough with multi-region 3D geometries.

What I finally did is to write my own pde solver code, using finite differences method.

Gabriel Madirolas

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Konstantin Kovalev 2019-6-18

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/417429-how-to-create-a-geometry-consisting-of-a-prism-with-an-inner-inclined-plane-for-a-diffusion-pde-pr#answer_379780

在 MATLAB Online 中打开

Hello Gabriel,

Please find below a script that does the job. Feel free to comment back should you have any questions.

Best Regards,

Konstantin

close all
clear variables
% columns 5 and 6 belong to upper wedge only, 7 and 8 - lower
x = [0 0 1 1 1 1 0 0];
y = [0 1 0 1 0 1 0 1];
z = [1 1 0 0 1 1 0 0];
% create delaunay triangulation for upper wedge
w1 = delaunayTriangulation(x(1:6)',y(1:6)',z(1:6)');
% copy elements (the copy will be lower wedge)
e1 = w1.ConnectivityList;
e2 = e1;
% update nodes 5-6 with nodes 7-8 (lower elements are inverted, V<0)
e2(e2==5) = 7;
e2(e2==6) = 8;
% combine both wedges' elements
% (1st two nodes of lower elements are swapped to undo inversion)
e = [e1; e2(:,[2 1 3:end])]';
% domain IDs - elements 1-3 are upper wedge, 4-6 lower
d = [1 1 1 2 2 2]';
pdem = createpde;
g = geometryFromMesh(pdem, [x; y; z], e, d);
pdegplot(pdem,'FaceLabels','on','CellLabels','on','FaceAlpha',0.75);
return