#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Created on Tue Oct  3 18:56:42 2023

@author: pimazio
"""

import numpy as np
from numpy import linalg
from matplotlib import cm
from matplotlib import pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
from mpl_toolkits import mplot3d

x = np.linspace(-2,2,401) #spatial coordinate
Nx = np.size(x)

A_1 = 1
x_1 = 0.5
sigma_1 = 0.6 

A_2 = 1.2
x_2 = -0.5
sigma_2 = 0.3

dt = 0.01
Nt = 1001
tend = dt*(Nt-1)
t = np.linspace(0,tend,Nt)
omega_1 = 1.3
omega_2= 4.1

y_1 = A_1*np.exp(-((x-x_1)**2)/(2*sigma_1**2))
y_2 = A_2*np.exp(-((x-x_2)**2)/(2*sigma_2**2))

Y = np.zeros([Nx,Nt],dtype='d')

for tt in range(Nt):
    Y[:,tt] = y_1*np.sin(2*np.pi*omega_1*t[tt])+y_2*np.sin(2*np.pi*omega_2*t[tt])

plt.figure(1)    
plt.plot(x,y_1)
plt.plot(x,y_2)
plt.xlabel('x')
plt.ylabel('y_1,y_2')


T_grid,Y_grid = np.meshgrid(t,x)
plt.figure(2)
plt.contour(T_grid,Y_grid,np.abs(Y))
plt.xlabel('t')
plt.ylabel('x')
plt.xlim(0,2);

U,S,VT = linalg.svd(Y,full_matrices=False)   #Este modo es el economico. Si pongo True, me tira las matrices sin recortar
print(np.shape(U))
print(np.shape(S))
print(np.shape(VT))

#Buscamos los valores singulares
plt.figure(3)
plt.plot(S,'o-')
plt.xlabel('Indice')
plt.ylabel('Valores Singulares')
plt.xlim(-1,10)

plt.figure(4)
plt.plot(x,U[:,0])
plt.plot(x,U[:,1])
plt.xlabel('x')
plt.ylabel('u1(x),u2(x)')
plt.title('Modos POD')

plt.figure(5)
plt.plot(t,VT[0,:])
plt.plot(t,VT[1,:])
plt.xlabel('t')
plt.xlim(0,5)
plt.ylabel('a1(t),a2(t)')
plt.title('Modos POD')
plt.show()