# -*- coding: utf-8 -*-
"""

Exercice 13-7b : résolution de l'équation des ondes avec 
contditions aux limites périodique (méthode des éléments finis)
Created on Tue Sep 12 13:29:49 2017

@author: ribot/grivet
"""

# Vitesse
c = 1
# Discrétisation en espace
xmin = 0; xmax = 1; h = 0.01
# la dernière valeur est exclue
xx = arange(xmin,xmax+h,h)
npt = size(xx); nx = npt-1
xxint = xx[0:nx]
# Matrice de discrétisation en espace 1D
K = 2*eye(nx,nx)-diag(ones(nx-1),1)-diag(ones(nx-1),-1); 
K[nx-1,0] = -1; K[0,nx-1] = -1
K1D = K/h**2;
# Discrétisation en temps
Tfin = 1; dt=0.95*h/c
ntps = ceil(Tfin/dt)
# Donnée initiale 
f = 5*sin(2*pi*xxint) + 4*sin(4*pi*xxint)
g = zeros(nx)
u0 = f; u1 = u0 + dt*g
# Résolution en temps - Euler Explicite
for t in range(int(ntps)):
   u = 2*u1-c**2*dt**2*dot(K1D,u1)-u0;
   u0 = u1
   u1=u; 
   if t % 10 == 0:
      uu = hstack((u,u[0]))
      plot(xx, uu);
uu = hstack((u,u[0]));
plot(xx, uu,'*')
uu = hstack((f,f[0]))
plot(xx, uu,'r')