x=0:2*pi/1000:2*pi;
k1=40;k2=41;k3=42;k4=43;k5=44;
v=10;
alpha=0.15;
om1=sqrt(v^2*k1^2+alpha*k1^4);
om2=sqrt(v^2*k2^2+alpha*k2^4);
om3=sqrt(v^2*k3^2+alpha*k3^4);
om4=sqrt(v^2*k4^2+alpha*k4^4);
om5=sqrt(v^2*k5^2+alpha*k5^4);

k=0:1:70;
omeganond=v*k;
omega=omeganond;
for j=1:71 
    omega(j)=sqrt(v^2*k(j)^2+alpha*k(j)^4); 
end
hold on
plot(k,omega,'r');
plot(k,omeganond);
hold off

%%
ntim=10000;Tmax=10*2*pi;
t=1:Tmax/ntim:Tmax;

% G = zeros(ntim);

for j=1:ntim-1
    
    y1=sin(k1*x-om1*t(j));
    y2=0*sin(k2*x-om2*t(j));
    y3=0*sin(k3*x-om3*t(j));
    y4=0*sin(k4*x-om4*t(j));
    y5=sin(k5*x-om5*t(j));
    
    y=y1+y5+y3+y4+y5;
    
    plot(x,y,x,y1+10,x,y2+13,x,y3+16,x,y4+19,x,y5+22);
    GG(j)=getframe;
end