function dy = euler_angles(t,y)
    dy = zeros(3,1);    % a column vector
    dy(1) = cos(y(3))/cos(y(2))*w1 - sin(y(3))/sin(y(2))*w2;
    dy(2) = sin(y(3))*w1+cos(y(3))*w2;
    dy(3) = -sin(y(2))*cos(y(3))/cos(y(2))*w1 + sin(y(2))*sin(y(3))/cos(y(2))*w2 + w3 ;