// Copyright jr Licois et Fédération Denis Poisson
//  Modèle avec vitesse de type impulsion
// constante K  atténuation de la vitesse
// constante c c0 c1 a b  voir modèle vitesse
// constante  kb  disparition des enzymes
// constante kt  absorption des nutriments.
// constante C  correspond à 36\pi/\rho^2 ....
// x(1) = position
// x(2) = derivee position (ie vitesse)
// x(3) = nutriments
// x(4) = enzymes
// x(5) = digérés
function f=dig2(t,x)
  f(1)=x(2);
  f(2)=-K*x(2) + W(t-x(1)/c)*(c0+c1*x(3))/(a+b*x(1))*(1-x(2)/c) ;
  // pour eviter les problemes on prend la valeur absolue
  f(3)=-C*ph(x(1))*x(4)*abs(x(3))^(2/3);
  f(4)=-kb*x(4);
  f(5)=C*ph(x(1))*abs(x(3))^(2/3)*x(4)-kt*x(5);
endfunction

// activité fonction du ph
function kph=ph(x)
  kph = 2.0*abs(x)/3-5*abs(x-3/4)/3+67*abs(x-1)/68+abs(x-18)/68;
endfunction

// echelon de longueur alpha entre 0 et 1
function y=ech(x)
  y=(alpha+abs(x)-abs(x-alpha) )/2/alpha;
endfunction

// derivee de la precedente
function z=indic(x)
  if x  <0 then z=0; 
    elseif   x<alpha then z=1/alpha; else z= 0; end;
endfunction


// somme d'echelons espaces de tau
//function y=W(x)
//  if x < 0 then y=0; 
//  else y=indic(6*x-int(6*x));end;
//endfunction

function yw=W(x)
  //y= indic(2*x)+indic(2*x-1/2)+indic(2*x-1)+indic(2*x-3/2);
  yw = tau*indic(tau*x-fix(tau*x));
endfunction