kill(all)$
load(lrats)$

numerical methods :
=================

appendix : rk with resistance ' -k*m*v ',       version 1.0
peter.vlasschaert@gmail.com,30/07/2018.

info : ' rk = rk4 ', runge-kutta 4e th order
=================================
rk : command ' Maxima '

How explain : syntax 'rk'
====================
====================

general :
=======
--------------------------------------------------------------------------------------------------------
--------------------------------------------------------------------------------------------------------
domain : [indepent var,start_time,end_time,increment]
ODEi : rhs(ODEi)                                                   ,  i=1,2,..,n
--------------------------------------------------------------------------------------------------------
--------------------------------------------------------------------------------------------------------

one diff eq :
==========
ODE1                               =>

rk([ODE1],var,init,domain)  
--------------------------------------------------------------------------------------------------------
More diff eq's
============
ODE1,ODE2,.....,ODEn =>

rk([ODE1,ODE2,..,ODEn],[var1,var2,..,varn],[init1,init2,..,initn],domain)
-------------------------------------------------------------------------------------------------------

rem : notation same for ' ode2 '

example : 
--------------
=====================================
m*dv/dt = m*g-k*v => m*'diff(v(t),t) =m*g-m*k*v
=====================================
result    =>                                     ODE1 :   g-k*v
=====================================

p1:m*'diff(v,t) = m*g-m*k*v;

p1a:p1/m;

p2a:ratsimp(p1a);

rem : v=vx,v=vy ; gx=g=0,gy=-g

components :  differential equations {vx,vy}
=====================================

p3:subst([v=vx,g=0],p2a);
p3a:subst([v=vy,g=-g],p2a);

rhs : p3,p3a ?
============

ODE1:rhs(p3);
ODE2:rhs(p3a);

p3b:[ODE1,ODE2];

variables : ?
==========

g = 9.81 m/s^2   ( => MKS)

g:-9.81;
k: 0.6;

var:[vx,vy];

U =2 m/s            ( => MKS)
α = π/4 rad        ( => MKS)

U:2;
α: %pi/4;
vx0:U*cos(α);
vy0:U*sin(α);

init:[vx0,vy0];

start_time : 0;
end_time : 1.5;
increment : 0.1;

domain:[t,start_time,end_time,increment];

numerical solution : ODE1,ODE2.
============================
use rk ?

p4:rk(p3b,var,init,domain)$

result : table
==========

     ->  t,vx,vy
==========
----------------------------------------------------------------------------------------------
first row          ,  'p4a'  :  start_time                        ,  vx0     , vy0
second  row  , 'p4a'   :  start_time + increment ,new vx, new vy
---------------------------------------------------------------------------------------------

p4a:table_form(p4);

?  How extract data from  'p4a'
=========================
rem : table 

convert : to matrix ?
=================
------------------------------------------------------------------------------------------------------
rem : doesn't work ,only works ' two columns ' here 'three columns'
=>  matrix(p4a)  :  result : ' error'
------------------------------------------------------------------------------------------------------

other way ?

sol : p4;

 How to 'makelist ':  from table 'p4a'
===============================
  =>       first column ,second column   
get from : ' sol '

p5:makelist([z[1],z[2]],z,sol);

make graph :
===========

p5a:wxplot2d([discrete,p5],[style,points],[color,blue],[point_type,diamond],[xlabel,"t"],[ylabel,"vx(t)"],grid2d);

 How to 'makelist ':  from table 'p4a'
===============================
  =>       first column ,third column   
get from : ' sol '

p6:makelist([z[1],z[3]],z,sol);

make graph :
===========

p6a:wxplot2d([discrete,p6],[style,points],[color,red],[point_type,diamond],[xlabel,"t"],[ylabel,"vy(t)"],grid2d);

make graph : "t,vx(t),vy(t)". from ODE1,ODE2.
======================================

p7:[[discrete,p5],[discrete,p6]]$
p7a:[legend,"vx(t)","vy(t)"]$
p7b:wxplot2d(p7,[style,lines,lines],[color,blue,red],[xlabel,"t"],[ylabel,"x(t),y(t)"],p7a,grid2d);

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