function u = f_RKG(z) global RKG xOld tau temp1 = pendulum(z(1: 4)+xOld); temp2 = pendulum(z(5: 8)+xOld); temp3 = pendulum(z(9:12)+xOld); u(1: 4) = z(1: 4) - tau*( RKG(1,1)*temp1 + RKG(1,2)*temp2 + RKG(1,3)*temp3 ); u(5: 8) = z(5: 8) - tau*( RKG(2,1)*temp1 + RKG(2,2)*temp2 + RKG(2,3)*temp3 ); u(9:12) = z(9:12) - tau*( RKG(3,1)*temp1 + RKG(3,2)*temp2 + RKG(3,3)*temp3 ); end