3
$\begingroup$

Hi I have a code below that solves non linear coupled PDE's given Dirichlet boundary conditions. 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? Thanks for the help.

Note, $t\geq 0$ and $x\in [0,1]$. Here are the coupled equations, below that I provide my code

\begin{eqnarray} \partial^2_tu(x,t) +\partial_t u(x,t) -a u(x,t)-\partial^2_x u(x,t)+b u(u^2+f^2)=0\\ \partial^2_tf(x,t) +\partial_t f(x,t) -a f(x,t)-\partial^2_x f(x,t)+bf(u^2+f^2)=0\\ u(x,0)=f(x,0)=z(x),\ \partial_t u(x,0)=\partial_t f(x,0)=0 \end{eqnarray} 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 $\partial_x u(0,t)=\partial_x u(1,t)$, the same holds for $f(x,t)$.

My code is below

program coupledPDE

integer, parameter :: n = 10, A = 20 !total number of position and time steps, respectively
real, parameter :: h = 0.1, k = 0.05 !step size for position and time respectively , x goes from 0 to 1, t goes from 0 to 1 also.
real, dimension(0:n) :: u,v,w,f,g,d !these variables store my solutions.
integer:: i,m !loop indices for space and time respectively
real:: t,R,x,c1,c2,c3,eta,a,b  

R=(k/h)**2.
a=1.0
b=1.0
c1=(2.+a*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.)


u=0.0 !Boundary conditions, need to make periodic..
f=0.0 !boundary conditions, but not correct.

do i = 0,n+1 
  x = real(i)*h    
  w(i) = z(x)  !u(x,0)
  d(i) = z(x)  !f(x,0)
end do

do i=1,n !discretization for 1st time step only, as usual for hyperbolic PDE's
  v(i) = (c1/(1.+c3))*w(i) + (c2/(1.+c3))*(w(i+1)+w(i-1)) -(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(i+1)+d(i-1)) -(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=1,n-1 !Discretization equation for all times after the 1st step
       u(i)=c1*v(i)+c2*(v(i+1)+v(i-1))-c3*w(i)-c4*v(i)*((v(i))**2.+(g(i))**2.) 
       f(i)=c1*g(i)+c2*(g(i+1)+g(i-1))-c3*d(i)-c4*g(i)*((v(i))**2.+(g(i))**2.) 
   end do 
   write(*,*)'check if u-f=0', u-f !this should be zero
   print*, "the values of f(x,t+k) for all m=",m
   print "(//3x,i5,//(3(3x,e22.14)))",m,f,u

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, I am new to here, so if I should reformat my question in a more proper way please let me know.

$\endgroup$
1
  • 6
    $\begingroup$ You need to map the left and right nodes. Have you checked this questions? 1 2 $\endgroup$
    – nicoguaro
    Mar 10 '16 at 1:01
6
$\begingroup$

Typically, you would add "guard cells", that is (for u) u(-1) and u(n+1) with your notation. Before each integration step:

u(n+1) = u(0)
u(-1) = u(n)

and similarly for the other variables. If you use higher order derivatives, you could also define u(-2) and u(n+2), etc.

$\endgroup$
1
  • 2
    $\begingroup$ +1 for using ghost nodes, incidentally my personal favorite $\endgroup$
    – nluigi
    Mar 10 '16 at 9:19

Not the answer you're looking for? Browse other questions tagged or ask your own question.