Showing posts with label Dirichlet. Show all posts
Showing posts with label Dirichlet. Show all posts

Sunday, 24 April 2022

1D Laplace's equation using Finite Difference Method

This post is about FDM for Laplace Equation with various boundary conditions.

MATLAB Code (1D, Dirichlet Boundary Conditions)

%% initialize the workspace, clear the command window


clear; clc


%% finite difference 1D laplace dirichlet boundary conditions %% d2u/dx2 = 0 %% u(o) = 10, u(L) = 4 %% Ax=b%%


N = 4 ; %number of grid points

L = 1; %length of domain

dx = L/(N-1); %element size


%% initialize variables %%


l = linspace(0,L,N); %independent

u=zeros(1,N); %dependent


%% boundary conditions %%


u(1)=10;

u(end)=4;


%% b vector %%


b=zeros(N-2,1);

b(1) = b(1) - u(1);

b(end) = b(end) - u(end);


%% A matrix


A = -2*eye(N-2,N-2);

for i=1:N-2

    if i<N-2

        A(i,i+1)=1;

    end

    if i>1

        A(i,i-1)=1;

    end

end


%% solve for unknowns %%


x = A\b;


%% fill the u vector with unknowns %%


u(2:end-1) = x;


%% plot the results %%


hold on; grid on; box on, grid minor

set(gca,'FontSize',40)

set(gca, 'FontName', 'Times New Roman')

ylabel('u','FontSize',44)

xlabel('l','FontSize',44)


plot(l,u,'-o','color',[0 0 0],'LineWidth',2,'MarkerSize',20)


MATLAB Code (1D, Mixed Boundary Conditions)

%% initialize the workspace, clear the command window

clear; clc

%% finite difference 1D laplace mixed boundary conditions %% d2u/dx2 = 0 %% u(o) = 10, du/dx(L) = 4 %% Ax=b%%

N = 5; %number of grid points
L = 1; %length of domain
dx = L/(N-1); %element size
a = 4;
%% initialize variables %%

l = linspace(0,L,N); %independent
u=zeros(1,N); %dependent

%% dirichlet boundary condition %%

u(1) = 10;

%% b vector %%

b=zeros(N-1,1);
b(1) = b(1) - u(1);
b(end) = b(end) + dx*a; %neumann boundary condition added to b vector

%% A matrix

A = -2*eye(N-1,N-1);
for i=1:N-1
    if i<N-1
        A(i,i+1)=1;
    end
    if i>1
        A(i,i-1)=1;
    end
end
A(N-1,N-1) = -1; %neumann boundary condition added to A matrix
%% solve for unknowns %%

x = A\b;

% fill the u vector with unknowns %%

u(2:end) = x;

%% plot the results %%

hold on; grid on; box on, grid minor
set(gca,'FontSize',40)
set(gca, 'FontName', 'Times New Roman')
ylabel('u','FontSize',44)
xlabel('l','FontSize',44)

plot(l,u','-o','color',[0 0 0],'LineWidth',2,'MarkerSize',20)