% VF_Bellman_PREC.m is a routine for computing value function on the left 
% side of the Bellman equation; see "How to Solve Dynamic Stochastic Models
% Computing Expectations Just Once" by Kenneth L. Judd, Lilia Maliar, Serguei 
% Maliar and Inna Tsener, (2016), Quantitative Economics, forthcoming, 
% (henceforth, JMMT, 2016).
% -------------------------------------------------------------------------

% -------------------------------------------------------------------------
% Inputs:    "c0" is current consumption;
%            "k1" is next-period capital;
%            "z0" is current productivity;
%            "gam" is the utility-function parameter;
%            "beta" is the discount factor; 
%            "vcv" is the variance of the shock;
%            "rho" is the autocorrelation coefficient in the process for
%                shock;
%            "vf_coef" is the vector of polynomial coefficients in the 
%               approximated value function;  
%            "D" is the degree of polynomial approximation
%
% Output:    "V_new" is value function on the left side of the Bellman
%            equation 
% -------------------------------------------------------------------------
% Copyright © 2016 by Lilia Maliar, Serguei Maliar and Inna Tsener. All 
% rights reserved. The code may be used, modified and redistributed under 
% the terms provided in the file "License_Agreement.txt".
% -------------------------------------------------------------------------

function [V_new] = VF_Bellman_PREC(c0,k1,z0,gam,beta,vcv,rho,vf_coef,D)

% Precomputation of integrals in the Bellman equation
%-----------------------------------------------------
intV = Integrals_PREC_2d(D,vcv)'; % Precomputed integrals for constructing 
                                  % expectation functions
vf1_coef = vf_coef.*intV;         % A new vector of coefficients
                                  % with precomputed integrals

X1 = Polynomial_2d([k1 z0.^rho],D);% Construct polynomial   
EV = X1*vf1_coef;                 % Evaluate value function

% Compute value function using Bellman equation
%----------------------------------------------

if gam==1
    V_new = log(c0)+beta*EV; % Bellman equation if gam=1        
else
    V_new = (c0.^(1-gam)-1)/(1-gam)+beta*EV;
end                                       % Bellman equation otherwise