kill(all)$
load(lrats)$
load(vect)$

RE-ENTRY equations of motion 2D(g :neglected ,non-vertical)
peter.vlasschaert@gmail.com,28/08/2018

see also :RE-ENTRY equations of motion 2D( assume flat earth,vertical motion),website.

RE-ENTRY equations of motion 2D( assume flat earth).
==============================================
we start from : ODE1,ODE2
-----------------------------------------------------------------------------------

p1:'diff(Z,t,2)/g=(V^2*%e^(-Z/H)*sin(θ)*ρ0)/(2*B)-1;
p1a:'diff(X,t,2)/g=-(V^2*%e^(-Z/H)*cos(θ)*ρ0)/(2*B);

g : Z 'direction' ( 'Gravitation')
B = ballistic coef =  'B=(g*m)/(Cd*A[fr])', A[fr]=Af
V = velocity 
ρ = Density ( ρ0 @  Z = 0)
H = atmospheric scale height =' H=(R*T)/g ' = cte
θ = re-entry angle
X,Z = direction 

----------------------------------------------------------------------------------

p1b: K= (H*g*ρ0)/B;

assume :  g (neglected = 0, 'p2a')
===========================
--------------------------------------------------------
=>   good estimate.( = approximation)
--------------------------------------------------------

p2:V=%e^(-(K*%e^(-Z/H))/2)*sqrt(2*'integrate(%e^(K*%e^(-Z/H)),Z,Z,inf)*g+V0^2);

p2a:subst([g=0,abs(V0)= -V0],p2);

set : K1= K/2

p2a1:K/2 = K1;

p2a2:subst([p1b],p2a1);

depends([V0],t);

p2b:lratsubst([p2a1],p2a);

calculation : acceleration.
=====================

depends([Z],t);

p3:'diff(lhs(p2b),t)=diff(rhs(p2b),t);

assume : V0 = constant, 'diff(Z,t,1) = V , vertical  re-entry
---------------------------------------------------------------------------------

p3a:lhs(p3) =subst([('diff(V0,t,1)) = 0],rhs(p3));

p3b:subst([('diff(Z,t,1)) = rhs(p2b)],p3a);

rem :resubst : K1 values ,from 'p2a2'
-------------------------------------------------------

p3c:subst([rhs(p2a2)=lhs(p2a2)],p3b);

p3d:p3c/g;

find maximum deceleration. ( vertical entry)
====================================

rem : second part 'p3d' is a function ,  Z : F(Z)
--------------------------------------------------------------------

p4:F(Z) = rhs(p3d);

p4a:'diff(lhs(p4),Z) =diff(rhs(p4),Z);

   rem : max find 'p4a' : 'diff(F(Z),Z,1) = 0
==================================
=>   'diff(F(Z),Z,1) 'max' = 0  =>  Z = Zmax
==================================

p4a1:args(rhs(p4a));
p4a2:args(rhs(p4a))[1];
p4a3:args(rhs(p4a))[2];

rem : V0^2 >0, ρ0 > 0, %e^(-(H*%e^(-Z/H)*g*ρ0)/B-Z/H) ?

p4b:part(p4a2,3);

p4c:subst([Z=Zmax],p4b) = 0;

p4d:solve(p4c,Zmax);

insert , value 'p4d' into 'p3d'

p5:(lhs(p3d))[max] =subst([Z=Zmax],rhs(p3d));

p5a:lhs(p5)=subst(p4d,rhs(p5));

p5b:lhs(p5)=float(rhs(p5a));

rem : before : only vertical entry.

              Non -vertical Entry.
===============================
introduction : path length 's',Figure 1

scale height

rem : start from ,  'p2'

p6:subst([g=0],p2);

rem : exp(q) , q = number

change : q , for introduction of the ' scale height'

p6a1: Z/sin(θ[0])=s;
p6a2: H/sin(θ[0]);

p6b:lratsubst([p2a1,rhs(p2a2)=lhs(p2a2)],p6);

with use : 'p6a1,p6a2'

p6c:subst([Z=rhs(p6a1),H=p6a2],p6b);

rem : remove , s from the above : with 'p6a1'

p6:subst([rhs(p6a1)=lhs(p6a1),abs(V0)=V0],p6c);

find maximum deceleration. (non-vertical entry)
========================================

1e) acceleration 
==============

p7:'diff(lhs(p6),t)=diff(rhs(p6),t);

rem : V0 = constant ,('diff(V0,t,1)) = 0

p7a:lhs(p7) =subst([('diff(V0,t,1)) = 0],rhs(p7));

 2e)  ('diff(Z,t,1)) , rhs(p6)
====================

p7b:subst([('diff(Z,t,1))=rhs(p6)],p7a);

p7c:p7b/g;

3e) max deceleration
==================

p8:G(Z) = rhs(p7c);
p8a:'diff(lhs(p8),Z) =diff(rhs(p7c),Z);

when 'diff(G(Z),Z,1) = 0  ,then  Z=Zmax, see above. (why)
===============================================

p8b:((%e^(-Z/H)*g*ρ0)/(sin(θ[0])*B)-1/H) ;

from above : Z=Zmax , then  'p8b '= 0.
===============================

p8c:subst([Z=Zmax],p8b) =0;

p8d:solve(p8c,Zmax);

'diff(V,t,1)/g , 'max'  ,  Z =Zmax
=========================

p8e:(lhs(p7c))[max] =subst([Z=Zmax],rhs(p7c));

rem : introduction of the scale height , H => H/sin(θ[0])
=============================================

p8f:lratsubst([p8d[1],H=p6a2],p8e);

p8g:float(p8f);

--------------------------------------------------------------------------------------------------------------------------------------------------