Generating a 2D Composite Geometry with Subdomains for PDE Analysis in MATLAB

3 次查看（过去 30 天）
Darby Paez
Darby Paez 2024-11-29
编辑： Jiexian Ma 2025-4-5
Ran in:
I am trying to create a composite geometry for a PDE problem in MATLAB. My geometry consists of two regions:
  1. A general rectangular domain.
  2. A subdomain with an arbitrary shape (even more complex than sawtooth function), which is entirely contained within the rectangular domain.
I need to combine these two geometries into a PDEModel while preserving both domains as distinct subregions so that I can apply separate boundary conditions or specifyCoefficients to each "face". I also need to generate a high-quality finite element mesh that respects both regions. I attempted to use polyshape because I could generate complex geometries, not just squares, circles like decsg. I put an example generating both geometries and trying to combine using subtraction (same with union), but just one region is generate.
Here the example:
close all; clear; clc;
%% Step 1: Define Custom Geometry (Subdomain)
% here I am using sawtooth, but in my real problem is more complex geometry
T = 10 * (1 / 50); % Period of the sawtooth function
fs = 1000; % Sampling frequency
t = 0:1/fs:T-1/fs; % Time vector
fun = sawtooth(2 * pi * 50 * t); % Sawtooth function
ini_tickness = 1 * ones(size(fun)); % Base thickness
in_Domain = ini_tickness + normalize(fun, "range"); % Thickness profile normalized to range [1, 2]
%% Step 2: Define General Geometry (Domain)
Lx = 40; % Length of the rectangular domain [mm]
Ly = 20; % Height of the rectangular domain [mm]
pgon = polyshape([0, Lx, Lx, 0], [Ly, Ly, 0, 0]); % Rectangular domain
figure(1);
plot(pgon);
axis equal;
title('Rectangular Domain');
%% Step 3: Define the Subdomain Geometry
x_holo = linspace(0, Ly - 1, length(in_Domain)); % Adjust x-coordinates to fit domain
y_holo = in_Domain; % Variable thickness as a function of x
% Create the internal shape shape as a polyshape
holo_polygon = polyshape([x_holo, fliplr(x_holo)], [zeros(size(x_holo)), fliplr(y_holo)]);
Warning: Polyshape has duplicate vertices, intersections, or other inconsistencies that may produce inaccurate or unexpected results. Input data has been modified to create a well-defined polyshape.
holo_polygon = rotate(holo_polygon, -90, [20, 0]); % Rotate the polyshape
holo_polygon = translate(holo_polygon, -19, -0.5); % Translate to fit inside the rectangular domain
figure(2);
plot(holo_polygon);
axis equal;
title('internal Geometry');
%% Step 4: Combine the Geometries
% Subtract the hologram from the rectangular domain
% also with union, etc I'll obtain just one "face" and I need two.
polyout = subtract(pgon, holo_polygon);
figure(3);
plot(polyout);
daspect([1, 1, 1]);
title('Composite Geometry');
%% Step 5: Generate Mesh and Define PDEModel
% Attempt to create a triangulation from the composite geometry
tr = triangulation(polyout);
polyout = boundaryshape(tr); % Convert to a boundaryshape object
% Define PDEModel and import the geometry
pdem = createpde;
pdem.geometryFromMesh(tr.Points', tr.ConnectivityList');
figure(4);
subplot(2, 1, 1);
pdegplot(pdem, 'EdgeLabels', "on");
title('Edges of the Geometry');
subplot(2, 1, 2);
pdegplot(pdem, 'FaceLabels', "on");
title('Faces of the Geometry');
%% Step 6: Generate Mesh % Please notice that there is only one face
generateMesh(pdem, "GeometricOrder", "linear");
figure(5);
pdeplot(pdem);
title('Generated Mesh');

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

回答（2 个）

Jiexian Ma
Jiexian Ma 2025-4-5
编辑：Jiexian Ma 2025-4-5
Function generateMesh can handle 2d geometry with multi domains (defined by polyshapes), but it's not easy to play with.
To make life easier, you can use Im2mesh package to do that.
Matlab code
You could check demo14 to demo17 in Im2mesh package. I provide a few examples.
Divyam
Divyam 2024-12-10
Hi @Darby Paez,
The "polyshape" method doesn't directly support multiple regions in the form of subdomains; it typically combines them into a single region, and hence it isn't possible to define your polygons as distinct regions.
Since you already tried the "decsg" method and found it incompatible with your use case, you could try defining a function handle for the location of the sawtooth model that will add the necessary properties you need to the polygon based on its location while treating the rectangle and the sawtooth as a single face.
% Pseudocode
applyBoundaryCondition(model, "dirichlet", Edge = 1, r = @myboundaryfun);
function bcMatrix = myboundaryfun(location, state)
% add code here, can use the inpolygon method to check whether a point
% lies in the sawtooth polygon or the rectangle
end
For more information regarding the use of function handles to define properties, refer to this documentation: https://www.mathworks.com/help/pde/ug/steps-to-specify-a-boundary-conditions-object.html?#buftxz0-1:~:text=Q%2Cg%3DG)%3B-,Nonconstant%20Boundary%20Conditions,-Use%20functions%20to
As an alternative, you can try out the PDE Modeler app to interact with the polygons and then define boundary conditions specifically for the selected polygon: https://www.mathworks.com/help/pde/ug/pdemodeler-app.html

类别

在 Help Center 和 File Exchange 中查找有关 Geometry and Mesh 的更多信息

标签

产品


版本

R2024b

Community Treasure Hunt

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

Start Hunting!

Translated by