/* Author: Klaus Reiner Schenk-Hoppé (klaus@iew.unizh.ch) http://www.iew.unizh.ch/home/klaus Last update: February 1, 2000 gcc -lm average_cons_prob.c */ #include #include #include #include /*** global variables ***/ FILE *fp_mat; double preiterations, iterations, steps; double stoch, fcase; double smin,smax; double alpha, A; double *stayprob; double Exi,xi,xi1,xi2; double Edelta, delta, delta1, delta2; double En, n, n1, n2; double s, *Ek, *Ec; double k_cd(double k, int i) { return ((1.0-delta) * k + s * xi * A * pow(k,alpha))/(1.0+n); } double c_cd(double k, int i) { return (1-s) * xi * A * pow(k,alpha); } double k_ces(double k, int i) { return ((1.0-delta)*k+s*xi*pow(1-A+A*pow(k,alpha),1.0/alpha))/(1.0+n); } double c_ces(double k, int i) { return (1-s)*xi*pow(1-A+A*pow(k,alpha),1.0/alpha); } void exp_kc_cd(int i) { unsigned long j; double k,c; double z; k = 100.0; Ek[i] = 0.0; Ec[i] = 0.0; if (stoch == 0){ xi = Exi; n = En; delta = Edelta; for(j=1;j stayprob[i]){ if (xi != xi1) xi = xi1; else xi = xi2; } k = k_cd(k,i); } for(j=1;j stayprob[i]){ if (xi != xi1) xi = xi1; else xi = xi2; } k = k_cd(k,i); } } if (stoch == 2){ xi = Exi; delta = Edelta; for(j=1;j stayprob[i]){ if (n != n1) n = n1; else n = n2; } k = k_cd(k,i); } for(j=1;j stayprob[i]){ if (n != n1) n = n1; else n = n2; } k = k_cd(k,i); } } if (stoch == 3){ xi = Exi; n = En; for(j=1;j stayprob[i]){ if (delta != delta1) delta = delta1; else delta = delta2; } k = k_cd(k,i); } for(j=1;j stayprob[i]){ if (delta != delta1) delta = delta1; else delta = delta2; } k = k_cd(k,i); } } } } void exp_kc_ces(int i) { unsigned long j; double k,c; double z; k = 2.0; Ek[i] = 0.0; Ec[i] = 0.0; if (stoch == 0){ xi = Exi; n = En; delta = Edelta; for(j=1;j stayprob[i]){ if (xi != xi1) xi = xi1; else xi = xi2; } k = k_ces(k,i); } for(j=1;j stayprob[i]){ if (xi != xi1) xi = xi1; else xi = xi2; } k = k_ces(k,i); } } if (stoch == 2){ xi = Exi; delta = Edelta; for(j=1;j stayprob[i]){ if (n != n1) n = n1; else n = n2; } k = k_ces(k,i); } for(j=1;j stayprob[i]){ if (n != n1) n = n1; else n = n2; } k = k_ces(k,i); } } if (stoch == 3){ xi = Exi; n = En; for(j=1;j stayprob[i]){ if (delta != delta1) delta = delta1; else delta = delta2; } k = k_ces(k,i); } for(j=1;j stayprob[i]){ if (delta != delta1) delta = delta1; else delta = delta2; } k = k_ces(k,i); } } } } typedef struct { long type; /* type */ long mrows; /* row dimension */ long ncols; /* column dimension */ long imagf; /* flag indicating imag part */ long namlen; /* name length (including NULL) */ } Fmatrix; /* * savemat - C language routine to save a matrix in a MAT-file. * * Here is an example that uses 'savemat' to save two matrices to disk, * the second of which is complex: * * FILE *fp; * double xyz[1000], ar[1000], ai[1000]; * fp = fopen("foo.mat","wb"); * savemat(fp, 2000, "xyz", 2, 3, 0, xyz, (double *)0); * savemat(fp, 2000, "a", 5, 5, 1, ar, ai); * fclose(fp); * * Author J.N. Little 11-3-86 */ void savemat(FILE *fp, int type, char *pname, int mrows, int ncols, int imagf, double *preal, double *pimag) /* FILE *fp; File pointer */ /* int type; Type flag: Normally 0 for PC, 1000 for Sun, Mac, and */ /* Apollo, 2000 for VAX D-float, 3000 for VAX G-float */ /* Add 1 for text variables. */ /* See LOAD in reference section of guide for more info. */ /* char *pname; pointer to matrix name */ /* int mrows; row dimension */ /* int ncols; column dimension */ /* int imagf; imaginary flag */ /* double *preal; pointer to real data */ /* double *pimag; pointer to imag data */ { Fmatrix x; int mn; x.type = type; x.mrows = mrows; x.ncols = ncols; x.imagf = imagf; x.namlen = strlen(pname) + 1; mn = x.mrows * x.ncols; fwrite(&x, sizeof(Fmatrix), 1, fp); fwrite(pname, sizeof(char), x.namlen, fp); fwrite(preal, sizeof(double), mn, fp); if (imagf) { fwrite(pimag, sizeof(double), mn, fp); } } int main(int argc, char *argv[]) { unsigned long i; double ds; if (argc < 9) { fprintf(stderr,"parameters: production-function stayprob_min stayprob_max steps stochastic-case pre-iterations iterations saving-rate\n production-function (1=Cobb-Douglas, 2=CES)\n stayprob_min stayprob_max steps (transition probability p_11 = p_22 = p in [stayprob_min, stayprob_max] on an equidistant grid with steps+1 vertices))\n stochastic-case (0=deterministic, 1=production shocks (xi), 2=stoch. population growth rate (n), 3=stoch. depreciation (delta))\n pre-iterations (recommended > 100)\n iterations (recommended > 500,000)\n saving-rate (fixed in (0,1))\n"); fflush(stderr); exit(1); } else{ sscanf(argv[1], "%lf", &fcase); sscanf(argv[2], "%lf", &smin); sscanf(argv[3], "%lf", &smax); sscanf(argv[4], "%lf", &steps); sscanf(argv[5], "%lf", &stoch); sscanf(argv[6], "%lf", &preiterations); sscanf(argv[7], "%lf", &iterations); sscanf(argv[8], "%lf", &s); } stayprob = malloc(sizeof(double)*(steps+1)); Ek = malloc(sizeof(double)*(steps+1)); Ec = malloc(sizeof(double)*(steps+1)); ds = (smax-smin)/((double)steps); /* xi1 = 0.75; xi2 = 1.25; */ Exi = 1.0; xi1 = 0.95; xi2 = 1.05; xi = xi2; En = 0.0; n1 = -0.07; n2 = 0.07; n = n2; Edelta = 0.5; delta1 = 1.0/3.0; delta2 = 2.0/3.0; delta = delta2; if (fcase<2) { /* Cobb-Douglas production function */ A = 10.0; alpha = 0.75; for (i=0; i < steps+1; i++) { stayprob[i] = smin + i * ds; exp_kc_cd(i); } } else{ /* CES production function */ A = 0.5; alpha = 0.5; for (i=0; i < steps+1; i++) { stayprob[i] = smin + i * ds; exp_kc_ces(i); } } if ((fp_mat = fopen("average_cons_prob.mat","wb")) == NULL) { fprintf(stderr,"Unable to open output file average_cons_prob.mat \n"); fflush(stderr); exit(1); } savemat(fp_mat, 1000, "productionfct", 1, 1, 0, &fcase, (double *)0); savemat(fp_mat, 1000, "alpha", 1, 1, 0, &alpha, (double *)0); savemat(fp_mat, 1000, "A", 1, 1, 0, &A, (double *)0); savemat(fp_mat, 1000, "s", 1, 1, 0, &s, (double *)0); savemat(fp_mat, 1000, "preiterations", 1, 1, 0, &preiterations, (double *)0); savemat(fp_mat, 1000, "iterations", 1, 1, 0, &iterations, (double *)0); savemat(fp_mat, 1000, "stayprob", 1, steps+1, 0, stayprob, (double *)0); savemat(fp_mat, 1000, "Ec", 1, steps+1, 0, Ec, (double *)0); savemat(fp_mat, 1000, "Ek", 1, steps+1, 0, Ek, (double *)0); fclose(fp_mat); return 0; }