clc;
clear all;
R0 = [3102 5369 2625];          %[km]
V0 = [-6.426 0.7735 5.943];     %[km/s]
r = norm(R0);           % Initial radius  [km]
v = norm(V0);           % Initial speed   [km/s]
mu = 398600;            % Earth’s gravitational parameter [km^3/s^2]
a = mu/(2*mu/r - v^2);  % Semimajor Axis [km]
T = 2*pi*a^1.5/mu^0.5;  % Orbital period [s]
R_earth = 6378;         % Earth radius [km]
dt = T/1000;            % time step [s]
t = 0;                  % initial time
% Using fourth-order Runge–Kutta method to solve fundamental equation
% of relative two-body motion
F_r = @(R) -mu/(norm(R)^3)*R;
Rd = V0; R  = R0;
i = 1;
while (t <= T)
    Rv(i,:) = R;
    tv(i) =t;
    k_1 = dt*F_r(R);
    k_2 = dt*F_r(R+0.5*k_1);
    k_3 = dt*F_r(R+0.5*k_2);
    k_4 = dt*F_r(R+k_3);
    Rd  = Rd + (1/6)*(k_1+2*k_2+2*k_3+k_4);
    R   = R + Rd*dt;
    t   = t+dt;
    i = i+1;
end
rn = (Rv(:,1).^2+Rv(:,2).^2+Rv(:,3).^2).^0.5;
Alt = rn - R_earth; % Altitude Vector
Alt_max  = max(Alt); % [km]
fprintf('Maximum Altitude = %4.2f [km] \n',Alt_max);

Maximum Altitude = 6249.10 [km]

Ploting the altitude variation during one orbit

figure(1);
hold on;grid on;
plot(tv/3600,Alt);
ylabel('Altitude [km]');
xlabel('Time [hour]');
title('Satellite altitude variation during one orbit');