function delta=dUdt(~,U)
g=9.81;
delta=[U(2);-my(U(1))*g];