matlab script runs but will not run as a function

Hi,

I have a simple function that takes a scalar input and gives a scalar output. The function is below. If I assign a value at the beginning of the script the file runs and gives the correct output but when I treat it as a function and assign the value as an input I get the error:

??? Error using ==> sym.cat>checkDimensions at 76
CAT arguments dimensions are not consistent.

Error in ==> sym.cat>catMany at 39
[resz, ranges] = checkDimensions(sz,dim);

Error in ==> sym.cat at 29
    y = catMany(dim, strs);

Error in ==> sym.vertcat at 26
    y = cat(1,args{:});

Error in ==> ramseysteadystate at 49
f2   = [f21; f22; f23; f24; f25; f26; f27; f28; f29; f210; f211; f212; f213; f214; f215; f216];

HERE IS THE CODE:

%%THIS FUNCTION COMPUTES THE SSR
function [ssr] = ramseysteadystate(theta_1)
%%FIX A VALUE FOR INFL
%theta_1 = 0;
inflrss  =  (1.1)/(1+exp(.2*(theta_1(1)+1)));
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%THIS PROGRAM FINDS THE RAMSEY STEADY STATE USING THE LINEAR PRJECTION (OLS) METHOD OF SCHMITT-GROHE AND URIBE
syms v a z tau x ma y h c l d n rd rl infl l1 l2 l3 q f w r vlag alag zlag taulag xlag malag ylag hlag clag llag dlag nlag rdlag rllag infllag l1lag l2lag l3lag qlag flag wlag rlag yfut inflfut l1fut l2fut l3fut qfut real;
syms beta_1 theta phi chi eta omega nu inflss real;
%%f1 IS THE CURRENT VALUE (GAMMA T) AND f2 IS THE LAGGED VALUE (GAMMA T-1) AND u IS THE PERIOD UTILITY FUNCTION;
u = a*(log(c) - eta*h - (eta/chi)*((v*c)/ma)^(chi));
%%NOW WE NEED TO WRITE OUT 16 OR (N-1) EQUATIONS THAT FORM THE CONSTRAINT VECTOR f1 OF THE RAMSEY PROBLEM;
f11  = (ma - ((nu^(1/omega))*(n^((omega-1)/omega))+((1-nu)^(1/omega))*(d^((omega-1)/omega)))^(omega/(omega-1)));
f12  = (y-z*h);
f13  = -1*(l - (1-tau)*d);
f14  = (rd - 1 - (rl-1)*(1-tau) + x);
f15  = -1*(w*h + rd*d + f - c - rl*l - d + l);
f16  = (f - y*(1 - (w/z) - (phi/2)*(((infl/inflss)-1)^2)));
f17  = ((1-theta)*l3*y + theta*l3*w*h - phi*l3*y*(((infl/inflss)^2)-(infl/inflss)) + beta_1*phi*l3fut*yfut*(((inflfut/inflss)^2)-(inflfut/inflss)));
f18  = (z*l1*inflfut - beta_1*r*l1fut);
f19  = -1*(l1*q-l3*f-beta_1*l1fut*qfut);
f110 = (l1 - l3 - l2*(ma^(1/omega))*(n^(-1/omega))*(nu^(1/omega)));
f111 = (l1 - rd*l3 - l2*(ma^(1/omega))*(d^(-1/omega))*((1-nu)^(1/omega)));
f112 = -1*(l1 - rl*l3);
f113 = (eta*a - l3*w);
f114 = -1*(a*eta*(1 - ((v*c/ma)^(chi))) - c*l3);
f115 = -1*(eta*a*((v*c/ma)^chi) - l2*ma);
f116 = (z*l3*inflfut - beta_1*l1fut);
% NOW WE NEED TO WRITE OUT 16 OR (N-1) EQUATIONS THAT FORM THE CONSTRAINT VECTOR f2 OF THE RAMSEY PROBLEM. NOTICE, THIS IS JUST f1 LAGGED.
f21  = (malag - ((nu^(1/omega))*(nlag^((omega-1)/omega))+((1-nu)^(1/omega))*(dlag^((omega-1)/omega)))^(omega/(omega-1)));
f22  = (ylag-zlag*hlag);
f23  = -1*(llag - (1-taulag)*dlag);
f24  = (rdlag - 1 - (rllag-1)*(1-taulag) + xlag);
f25  = -1*(wlag*hlag + rdlag*dlag + flag - clag - rllag*llag - dlag + llag);
f26  = (flag - ylag*(1 - (wlag/zlag) - (phi/2)*(((infllag/inflss)-1)^2)));
f27  = ((1-theta)*l3lag*ylag + theta*l3lag*wlag*hlag - phi*l3lag*ylag*(((infllag/inflss)^2)-(infllag/inflss)) + beta_1*phi*l3*y*(((infl/inflss)^2)-(infl/inflss)));
f28  = (zlag*l1lag*infl - beta_1*rlag*l1);
f29  = -1*(l1lag*qlag-l3lag*flag-beta_1*l1*q);
f210 = (l1lag - l3lag - l2lag*(malag^(1/omega))*(nlag^(-1/omega))*(nu^(1/omega)));
f211 = (l1lag - rdlag*l3lag - l2lag*(malag^(1/omega))*(dlag^(-1/omega))*((1-nu)^(1/omega)));
f212 = -1*(l1lag - rllag*l3lag);
f213 = (eta*alag - l3lag*wlag);
f214 = -1*(alag*eta*(1 - ((vlag*clag/malag)^(chi))) - clag*l3lag);
f215 = -1*(eta*alag*((vlag*clag/malag)^chi) - l2lag*malag);
f216 = (zlag*l3lag*infl - beta_1*l1);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% NOW COMPUTE THE FOC'S FOR THE RAMSEY PROBLEM WITH THE LAGRANGE MULTIPLIERS (GAMMA T AND GAMMA T-1) IMPLIED.
f1   = [f11; f12; f13; f14; f15; f16; f17; f18; f19; f110; f111; f112; f113; f114; f115; f116];
f2   = [f21; f22; f23; f24; f25; f26; f27; f28; f29; f210; f211; f212; f213; f214; f215; f216];
vars = [ma; y; h; c; l; d; n; rd; rl; infl; l1; l2; l3; q; f; w; r];
% FIND THE JACOBIAN/GRADIENT OF u, 1xN.
ju = jacobian(u, vars.');
%%FIND THE JACOBIAN OF f1 (N-1 x N)
jf1 = jacobian(f1.', vars.');
%%FIND THE JACOBIAN OF f2 (N-1 x N)
jf2 = jacobian(f2.', vars.');
%%NOW APPLY SCHMITT-GROHE AND URIBE'S LINEAR PROJECTION METHOD TO SOLVE FOR THE LAGRANGE MULTIPLIERS THAT SATISFY THE FOC'S.
%%THIS INSIGHT FOLLOWS FROM WRITING THE FOC'S IN MATRIX FORM AS: ju - gamma(jf1 + (1/beta)*jf2) = 0.
%%THE ABOVE IS SIMPLY THE FOCs FOR OLS IF WE MULTIPLY BY X = (jf1 + (1/beta)*jf2)' and let Y=ju' and gamma = Bols' then we have Y'X - Bols'*X'X = 0.
%%HENCE,  WE HAVE Bols = (X'X)^(-1) * X'Y.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%THE FIRST STEP IS TO DROP THE TIME SUBSCRIPTS FROM THE FOC's COMPUTED BY MATLAB
jf2ss = subs(jf2, {vlag, alag, zlag, taulag, xlag, malag, ylag, hlag, clag, llag, dlag, nlag, rdlag, rllag, infllag, l1lag, l2lag, l3lag, qlag, flag, wlag, rlag}, {v, a, z, tau, x, ma, y, h, c, l, d, n, rd, rl, infl, l1, l2, l3, q, f, w, r});
jf1ss = subs(jf1, {yfut, inflfut, l1fut, l2fut, l3fut, qfut}, {y, infl, l1, l2, l3, q});
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%NOW, WE WANT TO FIND THE STEADY STATE VALUES FOR THE N ENDOGENOUS VARIABLES USING THE N-1 CONSTRAINTS. 
%%TO IDENTIFY THIS UNDER-IDENTIFIED SYSTEM, FIX A VALUE FOR INFL AND SOLVE THE REMAINING N-1 EQUATIONS AS A FUNCTION OF THIS GIVEN VALUE FOR INFL.
%%I DENOTE SS VALUES OF x (AS OPPOSED TO THE SYMBOL X) BY xrss
%%VALUES FOR THE EXOGENOUS VARIABLES ARE GIVEN BY
vrss     = .4;    %.28;      %.4;
arss     = 1;
zrss     = 1.005; %1.0035;  %1.005;
taurss   = .02;
xrss     = .01;  %.005294; %.01; %.007;
betarss  = .99;   %.9979;   %.99;
thetarss = 6;
phirss   = 50;
chirss   = 2;
etarss   = 2.5;
omegarss = 1.5; 
nurss    = .225;  %.25;      %.225;
infl1ss  = 1.005; %1.039; %1.0059; %1.005; %THIS IS THE CALIBRATED VALUE FOR INFLATION
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
rrss     = (zrss*inflrss)/betarss;
wrss     = (zrss/thetarss)*((thetarss-1) + (1-betarss)*phirss*(((inflrss/infl1ss)^2) - (inflrss/infl1ss)));
l3rss    = arss*etarss/wrss;
l1rss    = zrss*l3rss*inflrss/betarss;
rlrss    = l1rss/l3rss;
rdrss    = 1 - xrss + (rlrss-1)*(1-taurss);
%%THE NEXT VARIABLES DO NOT HAVE CLOSED FORM SOLUTIONS AND I USE ANOTHER FUNCTION TO FIND THEIR VALUES
options=optimset('Display','iter','MaxIter',100000,'MaxFunEvals',100000,'TolX',1e-12,'TolFun',1e-12);
[cstar,fval]  = fminsearch(@(c1) rssnonlinearhelper(c1,etarss,arss,vrss,l1rss,l3rss,omegarss,nurss,rdrss),0,options);
crss  =    1/(1+exp(.2*(cstar+1)));
marss = (vrss*crss)*(((arss*etarss)/(arss-crss*l3rss))^.5);
l2rss = etarss*arss*((vrss*crss)^2)/((marss)^3);
nrss  = ((l2rss/(l1rss - l3rss))^omegarss)*marss*nurss;
drss  = ((l2rss/(l1rss - rdrss*l3rss))^omegarss)*marss*(1-nurss);
lrss     =     (1-taurss)*drss;
hrss     =     (crss + rlrss*lrss + drss - lrss -rdrss*drss)/(wrss + zrss*(1 - (wrss/zrss) - (phirss/2)*((inflrss/infl1ss) - 1)^2));
yrss     =     zrss*hrss;
frss     =     yrss*(1 - (wrss/zrss) - (phirss/2)*((inflrss/infl1ss) - 1)^2);
qrss     =     (l3rss*frss)/(l1rss*(1-betarss));
ystarrss   =  zrss/etarss;
gstarrss   =  yrss/ystarrss;
nvrss      =  taurss*drss;
brss       =  nrss+nvrss;
velssrss   =  inflrss;
veldivrss  =  inflrss;
gyrss      =  zrss*inflrss;
truprss    =  (nurss*(rlrss-1)^(1-omegarss) + (1-nurss)*(rlrss-rdrss)^(1-omegarss))^(1/(1-omegarss));
rarss     =  rlrss - truprss;
msrss      =  nrss+drss;
uarss      =  truprss/rlrss;
unrss      =  (rlrss-1)/rlrss;
udrss      =  (rlrss-rdrss)/rlrss;
erss       =  unrss*nrss+udrss*drss;
snrss      =  unrss*nrss/erss;
sdrss      =  udrss*drss/erss;
munrss     =  inflrss*zrss;
mudrss     =  inflrss*zrss;
muqrss     =  inflrss*zrss;
piprss     =  1;
gmarss     =  inflrss*zrss;
gmsrss     =  inflrss*zrss;
ngmdrss    =  inflrss*zrss;
gmdrss     =  zrss;
gtruprss   =  1;
hsrss      = (1/chirss)*(vrss*crss/marss)^2;

yaux = [ystarrss, gstarrss, nvrss, brss, velssrss, veldivrss, gyrss, truprss, rarss, msrss, uarss, unrss, udrss, erss, snrss, sdrss,...
        munrss, mudrss, muqrss, piprss, gmarss, gmsrss, ngmdrss, gmdrss, gtruprss];
y2rss = [marss, yrss, hrss, crss, lrss, drss, nrss, rdrss, rlrss, inflrss, l1rss, l2rss, l3rss, qrss, frss, wrss, rrss, vrss, arss, zrss, taurss, xrss];

    %%FINALLY, WRITE THE STEADY STATE VALUES INTO A VECTOR FOR LATER USE.
    y1rss = [betarss, thetarss, phirss, chirss, etarss, omegarss, nurss, infl1ss, vrss, arss, zrss, taurss, xrss, marss, yrss, hrss, crss, lrss, drss, nrss, rdrss, rlrss, inflrss, l1rss, l2rss, l3rss, qrss, frss, wrss, rrss];
    

    %%NOW EVALUATE THE FOCS AT THE STEADY STATE COMPUTED ABOVE
    numju    = double(subs(ju, {beta_1, theta, phi, chi, eta, omega, nu, inflss, v, a, z, tau, x, ma, y, h, c, l, d, n, rd, rl, infl, l1, l2, l3, q, f, w, r}, y1rss));
    numjf1   = double(subs(jf1ss, {beta_1, theta, phi, chi, eta, omega, nu, inflss, v, a, z, tau, x, ma, y, h, c, l, d, n, rd, rl, infl, l1, l2, l3, q, f, w, r}, y1rss));
    numjf2   = double(subs(jf2ss, {beta_1, theta, phi, chi, eta, omega, nu, inflss, v, a, z, tau, x, ma, y, h, c, l, d, n, rd, rl, infl, l1, l2, l3, q, f, w, r}, y1rss));
    betaramsey = .99;
    xols   = (numjf1 + (1/betaramsey)*numjf2)';
    yols   = numju';

    gamma = xols\yols;

    %%FIND THE SSR
    ssr = (yols - xols*gamma)'*(yols - xols*gamma);
yallrss = [crss, drss, rlrss, rdrss, marss, inflrss, zrss, taurss, arss, vrss, hrss, yrss, lrss, nrss, l1rss, l2rss, l3rss, qrss, frss, wrss, rrss, xrss...
      ystarrss, gstarrss, nvrss, gyrss, rarss, msrss, uarss, unrss, udrss, erss, snrss, sdrss, muqrss, munrss, mudrss, piprss, gmarss, gmsrss, ngmdrss, gmdrss, truprss, gtruprss...
      brss, velssrss, veldivrss];

If I comment the second line and un-comment fourth line the code runs fine as script.

Also, the code below is the function minimized in the middle of the function ramseysteadystate.m

%%THIS FUNCTION IS USED IN ramseysteadystate TO SOLVE FOR crss, marss, l2rss, nrss,drss
function [error1] = rssnonlinearhelper(c1,etarss,arss,vrss,l1rss,l3rss,omegarss,nurss,rdrss);

c0  =    1/(1+exp(.2*(c1(1)+1)));
ma0 = (vrss*c0)*(((arss*etarss)/(arss-c0*l3rss))^.5);
l20 = etarss*arss*((vrss*c0)^2)/((ma0)^3);
n0  = ((l20/(l1rss - l3rss))^omegarss)*ma0*nurss;
d0  = ((l20/(l1rss - rdrss*l3rss))^omegarss)*ma0*(1-nurss);
ma1 = ((nurss^(1/omegarss))*(n0^((omegarss-1)/omegarss)) + ((1-nurss)^(1/omegarss))*(d0^((omegarss-1)/omegarss)))^(omegarss/(omegarss-1));
error1 = abs(ma1-ma0);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

Any help is appreciated,

Lee

Answer 1

[Jan]

How do you define the input theta_1 in the caller? Is it a SYM variable, while it is a DOUBLE, when you define it by theta_1 = 0 inside the script version? While the large posted code is not useful to identify the problem, the detail of how the input is defined is not explained yet.

Try this:

theta_1 = double(0);
ssr = ramseysteadystate(theta_1)

[Andrew]

Hi,

Thanks for the assistance. Yes, theta_1 is a double at all times, it is never treated as a sym. I should have clarified that when I said theta_1 is a scalar.

I tried your fix but I get the same error.

I posted the whole code in case anyone wanted to verify for themselves the issue. I've never had this happen where code runs as a script but when I treat it as a function it fails to run. I thought others may want to see for themselves.

Is this a bug with with the symbolic toolbox? I am running R2011a, which had some issues with the symbolic toolbox.

Thanks for your time,

Lee