```julia using Revise using TaskBasedProduction using SpecialFunctions using QuadGK

θ = 1.0 κ = 0.5 z = 1.2 αVec = [0.1, 0.2, 0.3] labor_input=[0.5; 0.04; 0.19;;]

initialguess=findinitialguess(θ, κ, z, αVec; threshold=1e-2) q, xT, fval = prodfun(laborinput, θ, κ, z, αVec; initialguess=initialguess, xtol=1e-10)

println("Quantity Produced: ", q) println("Task Thresholds: ", xT) println("Approximation error: ", fval)

Call unitInputDemand and print the output

laborinput2 = q*unitInputDemand( xT, θ, κ, z, αVec) println("Labor Demand: ", laborinput2) println("Error", laborinput2-laborinput)

Call margProdLabor with labor demand and print the output

MPL= margProdLabor(labor_input, θ, κ, z, αVec) println("Marginal Products of Labor (with labor demand): ",MPL)

MPL= margProdLabor(labor_input, θ, κ, z, αVec, xT) println("Marginal Products of Labor (with labor demand): ",MPL)

Call elasticity_substitution with labor demand, MPL, xT and parameters of the gamma function

ϵsub, ϵcompl=elasticitysubcomp(laborinput, θ, κ, z, αVec, MPL, xT) println("Allen partial elasticity of substitution:", ϵsub) println("Hicks partial elasticity of substitution:", ϵ_compl)

GENERAL PARAMETERIZATION OF FUNCTIONS

Define the density function b_g(x)

bg(x) = (x^(κ-1) * exp(-x/θ)) / (θ^κ * gamma(κ)) eh1(x)=exp(0.1x) e_h2(x)=exp(0.2x) eh3(x)=exp(0.3*x) eh = [eh1, eh2, eh3] # Example eh functions

initialguessgen=findinitialguessgen(z, bg, eh; threshold=1e-2, verbose=false) qgen, xTgen,fval= prodfungeneral(laborinput,z,bg, eh; initialguess=initialguess_gen)

laborinputgeneral =qgen* unitInputDemandgeneral(xTgen, z, bg, eh) println("Labor Demand: ", laborinputgeneral) isapprox(laborinput, laborinputgeneral, atol=1e-6)

println("Quantity Produced: ", qgen) println("Task Thresholds: ", xTgen) MPLgen=margProdLaborgeneral(laborinputgeneral, z, bg, eh, xTgen, qgen) ϵsubgen, ϵcomplgen=elasticitysubcompgeneral(laborinputgeneral, z, bg, eh, MPLgen, xT_gen)