clear

// ベクトルを図示
ome=0.5; // 角振動数
xmax = 8; // 計算横軸範囲
vmax = 4; // 計算軸範囲

x=-xmax:1:xmax;
v=-vmax:1:vmax;
dx=zeros(length(x),length(v));
dv=zeros(length(x),length(v));
for nv=1:length(v)
 for nx=1:length(x)
     dx(nx,nv)=v(nv);
     dv(nx,nv)=-ome^2 * x(nx) ;
 end
end

clf(); 
champ(x, v, dx, dv, 1)

// 軌道を上乗せ
t=0:0.01:13; // 計算時間範囲

A1=2; // 振幅その1
xe=A1*cos(ome*t)
ve=-ome*A1*sin(ome*t)
plot2d(xe,ve)

A1=4; // 振幅その2
xe=A1*cos(ome*t)
ve=-ome*A1*sin(ome*t)
plot2d(xe,ve)

A2=6; // 振幅その3
xe=A2*cos(ome*t)
ve=-ome*A2*sin(ome*t)
plot2d(xe,ve)

A3=8; // 振幅その4
xe=A3*cos(ome*t)
ve=-ome*A3*sin(ome*t)
plot2d(xe,ve)

zoom_rect([-xmax-0.5, -vmax-0.5, xmax+0.5, vmax+0.5]) // 図示範囲