#include <stdio.h>
#include <math.h>
#include <iostream>
#include <fstream>
using namespace std;
template<class Vector, class Real>
class test {
public:
  Vector compute(Vector& y, Real t) {
    Vector _result(1);
    _result [0]=-5.*y[0];
    return _result;
  };
  
};
template <class VectorFunction, class Vector, class Real>
void euler(Real& t, Real& tend, VectorFunction& f, Vector& y ){
  Real cou, h;
  long size = y.size();
  Vector z(size);
cout << z.size() << "  " << y.size() << endl;
  while ( t<tend ) {
    cou = f.cour( y, t );
    h=min(tend-t,cou);
    f.compute(y,t,z);
    for ( int i = 0; i < size; i++ ) {
      y[i]=y[i]+h*z[i];
    }
    t=t+h;
    f.isStepAccepted = true;
  };
  //delete z;
};
template <class VectorFunction, class Vector, class Real>
void midpoint(Real& t, Real& tend, VectorFunction& f, Vector& y ){
  Real cou, h;
  Vector z[1], z1[1];
  while ( t<tend ) {
    cou = cour( y, t, f );
    h=min(tend-t,cou);
    f.compute(y,t,z);
    z1=y+(0.5*h)*z;
    t=t+0.5*h
    f.compute(z1,t,z);
    y=y+h*z
    t=t+0.5*h;
  };
};
template <class VectorFunction, class Vector, class Real>
class extendVectorFunction {
  VectorFunction& m_f;
  public:
    extendVectorFunction(const VectorFunction& f) {
      m_f=f;
    };
    Vector operator*(const Vector& v){
      Vector _result=v+m_f*v;      
    };
};

