# How I write periodic BCs

3 views (last 30 days)
lulu on 1 Dec 2020
Commented: lulu on 1 Dec 2020
clear all;
xl=0; xr=1; %domain[xl,xr]
J=40; % J: number of dividion for x
dx=(xr-xl)/J; % dx: mesh size
tf=0.5; % final simulation time
Nt=50; %Nt:number of time steps
dt=tf/Nt;
c=50; % paremeter for the solution
mu=dt/dx;
if mu>1.0 %make sure dt satisfy stability condition
error('mu should<1.0!')
end
%Evaluate the initial conditions
x=xl:dx:xr; %generate the grid point
f=exp(-c*(x-0.2).^2); %dimension f(1:J+1)
%store the solution at all grid points for all time steps
u=zeros(J+1,Nt);
%find the approximate solution at each time step
for n=1:Nt
t=n*dt; % current time
gl = exp(-c*(xl-t-0.2)^2); % BC at left side
gr = exp(-c*(xr-t-0.2)^2); % BC at right side
if n==1 % first time step
for j=2:J % interior nodes
u(j,n)=f(j)-0.5*mu*(f(j+1)-f(j-1))+0.5*mu^2*(f(j+1)-2*f(j)+f(j-1));
end
u(1,n) = gl; % the left-end point
u(J+1,n) = gr; % the right-end point
else
for j=2:J % interior nodes
u(j,n)=u(j,n-1)-0.5*mu*(u(j+1,n-1)-u(j-1,n-1))+0.5*mu^2*(u(j+1,n-1)-2*u(j,n-1)+u(j-1,n-1));
end
u(1,n)=gl; % the left-end point
u(J+1,n)=gr; % the right-end point
end
%calculate the analytical solution
for j=1:J+1
xj=xl+(j-1)*dx;
u_ex(j,n)=exp(-c*(xj-t-0.2)^2);
end
end
%plot the analytical and numerical solution at different times
figure;
hold on;
n=10;
plot(x,u(:,n),'r.',x, u_ex(:,n),'r-');
n=30;
plot(x,u(:,n),'g.',x, u_ex(:,n),'g-');
n=50;
plot(x,u(:,n),'b.',x, u_ex(:,n),'b-');
legend('Numerical t=0.1','Analytical t=0.1','Numerical t=0.3','Analytical t=0.3','Numerical t=0.5','Analytical t=0.5');
title('Numerical and Analytical Solutions at t=0.1,0.3,0.5');
lulu on 1 Dec 2020
No, f is the initial condition
I want to change gl and gr to periodic boundary conditions. If the initial conditions leaves the domain at one end, and re-enters it at the other end then the Periodic BCs are correct.

### Categories

Find more on Partial Differential Equation Toolbox in Help Center and File Exchange

### Community Treasure Hunt

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

Start Hunting!

Translated by