function sphere_rk4 ( )

%*****************************************************************************80
%
%% sphere_rk4 tests sphere_ode with an RK4 method.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    24 February 2020
%
%  Author:
%
%    John Burkardt
%
  fprintf ( 1, '\n' );
  fprintf ( 1, 'sphere_rk4:\n' );
  fprintf ( 1, '  MATLAB/Octave version %s\n', version ( ) );
  fprintf ( 1, '  Solve sphere_ode using Runge Kutta 4.\n' );
  fprintf ( 1, '  sphere_ode models motion on the surface of a sphere.\n' );

  x0 = cos ( 0.9 );
  y0 = 0.0;
  z0 = sin ( 0.9 );

  f = @sphere_deriv;
  tspan = [ 0.0, 45.0 ];
  xyz0 = [ x0; y0; z0 ];
  n = 200;

  [ t, xyz ] = rk4 ( f, tspan, xyz0, n );

  figure ( 1 );
  plot ( t, xyz, 'linewidth', 3 );
  grid ( 'on' );
  xlabel ( '<-- X -->', 'Fontsize', 16 );
  ylabel ( '<-- Y -->', 'Fontsize', 16 );
  title ( 'RK4: (X,Y,Z) motion on a sphere' );
  legend ( 'X(t)', 'Y(t)', 'Z(t)' );
  filename = 'sphere_rk4_plot.png';
  print ( '-dpng', filename );
  fprintf ( 1, '  Graphics saved as "%s"\n', filename );

  figure ( 2 );
  plot3 ( xyz(:,1), xyz(:,2), xyz(:,3), 'linewidth', 3 );
  grid ( 'on' );
  xlabel ( '<-- X -->', 'Fontsize', 16 );
  ylabel ( '<-- Y -->', 'Fontsize', 16 );
  zlabel ( '<-- Z -->', 'Fontsize', 16 );
  title ( 'RK4: Motion on a sphere', 'fontsize', 24 );
  view ( 3 );
  filename = 'sphere_rk4_plot3d.png';
  print ( '-dpng', filename );
  fprintf ( 1, '  Graphics saved as "%s"\n', filename );

  h = xyz(:,1).^2 + xyz(:,2).^2 + xyz(:,3).^2;

  figure ( 3 );
  clf ( );
  hold ( 'on' );
  plot ( t, h, 'r-', 'linewidth', 3 );
  plot ( [ t(1), t(end) ], [0.0, 0.0 ], 'k-', 'linewidth', 3 );
  hold ( 'off' );
  grid ( 'on' );
  xlabel ( '<-- T -->', 'Fontsize', 16 );
  ylabel ( '<-- H(T) -->', 'Fontsize', 16 );
  title ( 'RK4: H(T) conservation', 'fontsize', 24 );
  filename = 'sphere_rk4_conserved.png';
  print ( '-dpng', filename );
  fprintf ( 1, '  Graphics saved as "%s"\n', filename );

  return
end
