function fval = euler(x,c,fspace,s,e,w)

global beta r ro gamma

ns=size(s,1);
a=s(:,1);
y=exp(s(:,2));

aprime=(1+r)*a+y-x;



fval=x.^(-gamma);


for k=1:length(w)              % go over all possible values e (discretized distribution of disturbances) can take
      
        kk = k+zeros(ns,1);
              
        cprime=funeval(c,fspace,[aprime,ro*log(y)+e(kk)]);
        fval=fval-beta*w(k)*cprime.^(-gamma);
        
               
end