periodic boundary conditions - finite differences

Question

Hi I have a code below that solves non linear coupled PDE's. However I need to implement periodic boundary conditions. The periodic boundary conditions are troubling me, what should I add into my code to enforce periodic boundary conditions? Updated based on modular arithmetic suggestions below.

Note, t>=0 and x is in the interval [0,1]. Here are the coupled equations, below that I provide my code

where a, b > 0.

Those are the initial conditions, but now I need to impose periodic boundary conditions. These can be mathematically written as u(0,t)=u(1,t) and du(0,t)/dx=du(1,t)/dx, the same holds for f(x,t). The du/dx I have for the boundary conditions are really meant to be partial derivatives.

My code is below

program coupledPDE 

integer, parameter :: n = 10, A = 20 
real, parameter :: h = 0.1, k = 0.05 
real, dimension(0:n-1) :: u,v,w,f,g,d 
integer:: i,m 
real:: t, R, x,c1,c2,c3,aa,b 

R=(k/h)**2.
aa=2.0
b=1.0
c1=(2.+aa*k**2.-2.0*R)/(1+k/2.)
c2=R/(1.+k/2.)
c3=(1.0-k/2.)/(1.0+k/2.)
c4=b*k**2./(1+k/2.)


do i = 0,n-1 !loop over all space points except 0 and n
  x = real(i)*h    
  w(i) = z(x)  !u(x,0)
  d(i) = z(x)  !f(x,0)
end do


do i=0,n-1
  ip=mod(i+1,n)
  il=modulo(i-1,n)
  v(i) = (c1/(1.+c3))*w(i) + (c2/(1.+c3))*(w(ip)+w(il)) -(c4/(1.+c3))*w(i)*((w(i))**2.+(d(i))**2.)    !\partial_t u(x,0)=0
  g(i) = (c1/(1.+c3))*d(i) + (c2/(1.+c3))*(d(ip)+d(il)) -(c4/(1.+c3))*d(i)*((w(i))**2.+(d(i))**2.)    !\partial_t f(x,0)=0
end do

do m=1,A 

   do i=0,n-1
       ip=mod(i+1,n)
       il=modulo(i-1,n)
       u(i)=c1*v(i)+c2*(v(ip)+v(il))-c3*w(i)-c4*v(i)*((v(i))**2.+(g(i))**2.) 
       f(i)=c1*g(i)+c2*(g(ip)+g(il))-c3*d(i)-c4*g(i)*((v(i))**2.+(g(i))**2.) 
   end do 
     print*, "the values of u(x,t+k) for all m=",m
   print "(//3x,i5,//(3(3x,e22.14)))",m,u   

  do i=0,n-1
   w(i)=v(i)
   v(i)=u(i)
   d(i)=g(i)
   t=real(m)*k
   x=real(i)*k
  end do

end do


end program coupledPDE

function z(x)
real, intent(in) :: x
real :: pi
pi=4.0*atan(1.0)
z = sin(pi*x)
end function z

Thanks for reading, if I should reformat my question in a more proper way please let me know.

Answer 1

One option to boundary conditions in PDE discretization is to use ghost (halo) cells (gridpoints). It may be not the most clever one for periodic BC, but it can be used for all other boundary condition types.

So you declare your arrays as

real, dimension(-1:n) :: u,v,w,f,g,d

but you solve your PDE only in points 0..n-1 (point n is identical with point 0). You could also do from 1..n and declare arrays form 0..n+1.

Then you set

 u(-1) = u(n-1)

and

 u(n) = u(0)

and the same for all other arrays.

At each time-step you set this again for u and f or all other fields that are modified during the solution:

do m=1,A 
   u(-1) = u(n-1)
   u(n) = u(0)
   f(-1) = f(n-1)
   f(n) = f(0)

   do i=0,n-1 !Discretization equation for all times after the 1st step
       u(i)=...
       f(i)=...
   end do 
end do

All above assumed explicit temporal discretization and spatial discretization with finite differences and it assumed that x(0) = 0 and x(n) = 1 are your boundary points.