% Integrals_Prec_2d.m is a routine for constructing integrals which are used
% for evaluation of conditional expectation functions; it computes integrals 
% in the form I=E(exp(l*e)), l=0,1,...,5 for ordinary polynomials of degrees 
% from one to five and a normally distributed random shock to the productivity 
% level with mean 0 and variance sigma^2; 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:    "D" is the degree of polynomial for which the integrals
%                need to be constructed; (can be 1,2,3,4 or 5)
%          "vcv" is the variance of the shock
% Output:    "x" is the vector of precomputed integrals for complete 
%                polynomial of the given degree D
% -------------------------------------------------------------------------
% 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 x = Integrals_PREC_2d(D,vcv)

% Construct integrals for the coefficients that multiply powers of k;
% see equaton (22) of JMMT (2016)
% ------------------------------------------------------------------
p1 = exp(0^2*vcv/2); 
p2 = exp(0^2*vcv/2); 
p3 = exp(0^2*vcv/2); 
p4 = exp(0^2*vcv/2); 
p5 = exp(0^2*vcv/2);

% Construct integrals for the coefficients that multiply powers of z
% ------------------------------------------------------------------
q1 = exp(1^2*vcv/2); 
q2 = exp(2^2*vcv/2);
q3 = exp(3^2*vcv/2);
q4 = exp(4^2*vcv/2);
q5 = exp(5^2*vcv/2); 

% Construct the vector of precomputed integrals for complete polynomial of
% the degree D=1,...,5
%--------------------------------------------------------------------------
if D==1;
    x = [1 p1 q1];
elseif D==2;
    x = [1 p1 q1 p2 p1.*q1 q2];
elseif D==3;
    x = [1 p1 q1 p2 p1.*q1 q2 p3 p2.*q1 p1.*q2 q3];
elseif D==4;
    x = [1 p1 q1 p2 p1.*q1 q2 p3 p2.*q1 p1.*q2 q3 p4 p3.*q1 p2.*q2 p1.*q3 q4];
elseif D==5;
    x = [1 p1 q1 p2 p1.*q1 q2 p3 p2.*q1 p1.*q2 q3 p4 p3.*q1 p2.*q2 p1.*q3 q4 p5 p4.*q1 p3.*q2 p2.*q3 p1.*q4 q5];
end    

% Note that "Integrals_PREC_2d.m" is a function that is used together 
% with "Polynomial_2d.m" and "Polynomial_deriv_2d.m". The inputs of these 
% polynomial routines are in the form s = [k, z]. This ordering of state 
% variables ensures that precomputed integrals match the corresponding 
% polynomial terms in the approximating function of degree D
