solving 2D Poisson equation having already generated a mesh

So i need to resolve the 2D poisson equation using Finite Element Method, i have an L shaped domain . So first of all i generated Triangular free mesh on the domain. but i dont know how to solve the eqution using my script. The conditions at the boundary is : u(x,t)=0

thank you.

% mesh generation

N = 1;

Visualize_Plot = true; % se vero visualizza il grafico

% Input Boundary :

% mi fa scegliere il file text dalla cartella

% [file,path] = uigetfile('*.txt');

% boundary = importdata(fullfile(path,file));

% oppure importo direttamente conoscendo il nome del file

boundary = importdata(fullfile('coordinate.txt'));

% dentro al file .txt c'è una matrice con prima colonna le x e la seconda colonna le y

tic

% if Visualize_Plot

% figure

% end

[T]= CreateMesh(boundary,N,Visualize_Plot);

min = min(min(boundary(:,1)),min(boundary(:,2))); % calcolo il valore minimo delle coordinate

max = max(max(boundary(:,1)),max(boundary(:,2))); % il valore max delle coordinate

margXY = 0.1*(max - min); % margine del 10%

if Visualize_Plot

figure(1)

xlim([min-margXY , max+margXY]) % imposto il limite dell'asse delle x

ylim(xlim) % stessa cosa per l'asse delle y

xlabel('x')

ylabel('y')

title('Mesh generation')

end

toc

function [T] = CreateMesh(boundary,N,PLOT)

% matrice che ha al suo interno la descrizione del poligono

Geom = [2 % ci indica che è un poligono

size(boundary,1)-1 % Numero di segmenti

boundary(1:end-1,1) % x-coordinate

boundary(1:end-1,2)]; % y-coordinate

[dl,~] = decsg(Geom);

r = min(sqrt(sum((boundary(1:end-1,1:2) - boundary(2:end,1:2)).^2,2))); % calcolo la distanza minima tra 2 vertici consecutivi del contorno

% genero una mesh iniziale con

[p,~,t] = initmesh(dl,'Hgrad',1.3,'Hmax',r/N);

% Hmax : la dimensione massima della mesh

% Hgrad : Mesh growth rate

T = zeros(size(t,2),6);

hold on

for i=1:size(t,2) % su ciascun triangolo nella mesh (t è il num di trangoli)

T(i,:) = [p(1,t(1,i)),p(2,t(1,i)),p(1,t(2,i)),p(2,t(2,i)),p(1,t(3,i)),p(2,t(3,i))];

% memorizza i vertici del traingolo

if PLOT

plot([p(1,t(1,i)),p(1,t(2,i)),p(1,t(3,i)),p(1,t(1,i))],...

[p(2,t(1,i)),p(2,t(2,i)),p(2,t(3,i)),p(2,t(1,i))],...

'color',[0.5 0.5 0.5])

end

end

plot(boundary(:,1),boundary(:,2),'b')

hold off

Answer 1

[Shivansh]

Hi Malak,

It seems like you have generated a mesh without using the MATLAB inbuilt functions. To solve the Poisson equation using the Finite Element Method (FEM), you'll need to assemble the system of equations that arises from the weak form of the Poisson equation and apply the boundary conditions.

You can refer to the following steps to proceed further:

1. Assemble the Stiffness Matrix and Load Vector:
2. Apply Boundary Conditions.
3. Solve the Linear System.

Here's a simplified version of the MATLAB code that you could use to perform these steps. You will need to define your own functions for computing the element stiffness matrix and load vector based on the specific form of your Poisson equation.

% Assuming you have the mesh points p and the triangulation t from your mesh generation
% p: 2xNp matrix containing the coordinates of the nodes
% t: 4xNt matrix containing the indices of the nodes for each triangle

% Number of points
Np = size(p, 2);

% Initialize the global stiffness matrix K and load vector F
K = sparse(Np, Np);  % Using sparse for efficiency
F = zeros(Np, 1);

% Assemble the stiffness matrix and load vector
for i = 1:size(t, 2)
    nodes = t(1:3, i);  % Nodes of the i-th triangle
    [Ke, Fe] = elementStiffnessMatrixAndLoadVector(p(:, nodes));
    K(nodes, nodes) = K(nodes, nodes) + Ke;
    F(nodes) = F(nodes) + Fe;
end

% Apply boundary conditions
% You will need to identify the indices of the nodes on the boundary
boundaryNodes = findBoundaryNodes(p, boundary);
K(boundaryNodes, :) = 0;
K(boundaryNodes, boundaryNodes) = speye(length(boundaryNodes), length(boundaryNodes));
F(boundaryNodes) = 0;

% Solve the linear system
u = K \ F;

% u now contains the solution at the mesh points

In this code, "elementStiffnessMatrixAndLoadVector" is a placeholder for a function that you need to implement which computes the local stiffness matrix "Ke" and load vector "Fe" for each element. The function "findBoundaryNodes" is also a placeholder for a function that identifies the nodes on the boundary of the domain.

You can also use to the PDE toolbox and refer to the inbuilt functions to generate mesh and solve the 2D Poisson equation. You can use the "generatemesh" function for mesh creation and "solvepde" tpo solve the eqyation. You can refer to the following documentation for more information:

1. generatemesh: https://www.mathworks.com/help/pde/ug/generate-mesh.html.
2. solvepde: https://www.mathworks.com/help/pde/ug/pde.pdemodel.solvepde.html.

I hope it helps!