% Parameters:
beta = .9;
R = 1.03;
theta = 2;
y = 1;
sigma = 0;

% Grid and initial guess:
X = [0:.1:10]';
copt = X+.1;

% Euler function iteration deterministic model:
stoch = 0;
figure(1); clf;
iter = 0; improv = 100;
while improv > .1
    cguess = copt;
    for i = 1:size(X)
        %euler(copt(i),1,X,cguess,stoch,beta,R,theta,y,sigma);
        [copt(i),accur] = fzero('euler',cguess(i),[],X(i),X,cguess,stoch,beta,R,theta,y,sigma);
        %accur
    end
    improv = sqrt( (copt-cguess)'*(copt-cguess) );
    iter = iter + 1;
    %disp(sprintf('Iteration: %g, improv = %g',iter,improv));
    %hold on; plot(X,copt);
end
plot(X,copt,'b:');

for sigma = .2:.2:1
sigma
copt = X+.1;

% Euler function iteration stochastic model:
stoch = 1;
iter = 0; improv = 100;
while improv > .1
    cguess = copt;
    for i = 1:size(X)
        %euler(copt(i),1,X,cguess,stoch,beta,R,theta,y,sigma);
        [copt(i),accur] = fzero('euler',cguess(i),[],X(i),X,cguess,stoch,beta,R,theta,y,sigma);
        %accur
    end
    improv = sqrt( (copt-cguess)'*(copt-cguess) );
    iter = iter + 1;
    %disp(sprintf('Iteration: %g, improv = %g',iter,improv));
    %hold on; plot(X,copt);
end
hold on; plot(X,copt,'r');

end