/* Fonction accel
   Syntaxe 1 : pair accel(path p, int t, int sign=0);
               If sign < 0, return the acceleration of the incoming path p at node t ; 
               if sign > 0, return the acceleration of the outgoing path. 
               if sign=0, the mean of these two accelerations is returned.
               
   Syntaxe 2 : pair accel(path p, real t);
               returns the acceleration of the path p at the point t.
*/
import geometry;
size(7.5cm,0);
path p=(-5,0)..(10,0)..(10,10)..(0,20)..(-10,20)..(-5,15)..cycle;
draw("p",p); dot(p,2bp+black);
pair po, pa;
for(real k=0; k<length(p); k+=.05)
{
  po=point(p,k);
  pa = po+(accel(p,k))/6;
  draw(po--pa,red);
}
void acceleration(path p, int k, real t){
     pair z0=point(p,k),
          c0=postcontrol(p,k),
          z1=point(p,k+1),
          c1=precontrol(p,k+1),
          m0=(1-t)*z0+t*c0,
          m1=(1-t)*c0+t*c1,
          m2=(1-t)*c1+t*z1,
          m3=(1-t)*m0+t*m1,
          m4=(1-t)*m1+t*m2,
          z2=(1-t)*m3+t*m4;
     path CourbeBezier = z0 .. controls c0 and c1 .. z1;
     draw(CourbeBezier,1bp+blue);
     path lignebrisee = z0 -- c0 -- c1 -- z1;
     draw(lignebrisee^^m0--m1--m2^^m3--m4,linetype("4 4"));
     dot(z0^^z1^^z2,blue); dot(c0^^c1,gray);
     dot(m0--m1--m2^^m3--m4,green);
     draw(z2--(z2+accel(CourbeBezier,t)/6),.8bp+.5green,Arrow());
     write((precontrol(p,k)+postcontrol(p,k))/2-point(p,k));
}
acceleration(p,1,.75);
acceleration(p,3,.75);
acceleration(p,5,.75);
addMargins(2mm,2mm);