value of type 'sym' is not convertible to 'double'
489 views (last 30 days)
Show older comments
Hey,
I am trying to let the following code run but I always get the same error:
Unable to perform assignment because value of type 'sym' is not convertible to 'double'.
Error in spdiags>makeSparseGeneral (line 225)
V(offset) = B(i + d(k), k);
Error in spdiags (line 125)
res1 = makeSparseGeneral(B, d, A);
Error in Combined_Model>spdiag (line 157)
x = spdiags(v,0,rv,rv);
Error in Combined_Model (line 76)
Js = spdiag(f(kc, kp) - (delta+n)*kc)*d1Mx*theta - ...
Caused by: Error using symengine. Unable to convert expression containing symbolic variables into double array. Apply 'subs' function first to substitute values for variables.
Can anyone help me with this?
%% COMBINED MODEL
clear; close all; clc;
%% Parameters
rho = 0.04;
theta = 1.5;
nu = 0.35;
delta = 0.05;
A=1;
L=1;
alpha = 0.5;
gamma = 0.5;
n = 0.01;
%% Functions
syms kc
Ec = (A*kc.^gamma*L.^(1-gamma)) ;
syms kp
Ep = (A*kp.^gamma*L.^(1-gamma)) ;
E = Ec + Ep;
fs = A*(E.^gamma)*(L.^(1-gamma));
f_p1s_c = diff(fs, kc);
f_p1s_p = diff(fs, kp);
f = matlabFunction(fs,'vars',{kc kp});
f_p1_c = matlabFunction(f_p1s_c,'vars',{kc kp});
f_p1_p = matlabFunction(f_p1s_p,'vars',{kc kp});
u = @(c) (integral(exp(-(p-n))*c.^(1-theta)/(1-theta)));0;inf ;
% %% Steady state
% % SS for Ramsey, GR for Solow
% % Numerical solution using unidimensional Newton-Raphson
%
% metric=1; i=0; kSS=1;
% while (metric>=1e-10)
% kSS = ((delta +rho)/(A*alpha)).^(1/alpha - 1) ;
% metric = ((delta +rho)/(A*alpha)).^(1/alpha - 1); i=i+1;
% end
% cSS = f(kSS) - (delta+n)*kSS;
kSS_c = ((delta +rho)/(A*gamma)).^(1/(gamma - 1)) ;
kSS_p = ((delta +rho)/(A*gamma)).^(1/(gamma - 1)) ;
%kSS = [kSS_c kSS_p]
cSS = f(kSS_c, kSS_p) - (delta+n)*kSS_c ;
%% Grid, differentiation operators
Ik = 2000; kmin = 1e-5*kSS_c; kmax = 5*kSS_c;
k = linspace(kmin,kmax,Ik)';
dk = (kmax-kmin)/Ik;
d1Mx = spdiags(repmat([-1 sparse(1,1) 1],Ik,1),-1:1,Ik,Ik)/2;
d1Mx(1,1)=-1; d1Mx(1,2)=1; d1Mx(Ik,Ik-1)=-1; d1Mx(Ik,Ik)=1;
d1Mx = d1Mx/dk;
%% Main loop - Euler equation + Newton-Raphson (modified reflection method)
% Initialisation
metric=1; i=0;
c=k;
% Static (c-independent) component of the Jacobian
Js = spdiag(f(kc, kp) - (delta+n)*kc)*d1Mx*theta - ...
spdiag(f_p1_c(kc, kp)-delta - rho);
Jd = sparse(Ik,Ik);
tic
while metric > 1e-10
% Dynamic (c-dependent) component of the Jacobian
% Jd = spdiags([-0.5*c [0; 0.5*(c(1:N-2)-c(3:N)); 0] 0.5*c],-1:1,N,N)';
% Jd(1,1) = 2*c(1)-c(2); Jd(1,2) = -c(1);
% Jd(N,N-1) = c(N); Jd(N,N) = c(N-1) - 2*c(N);
Jd = -(spdiag(d1Mx*c)+spdiag(c)*d1Mx)*theta;
% Merit function/Error vector (from the Euler equation)
err = (d1Mx*c).*(f(k) - (delta+n)*k - c)*theta - ...
(f_p1(k)-delta - rho).*c;
c = c - (Js+Jd)\err;
metric = max(abs(err));
i=i+1;
end
toc
if metric<=1e-10
fprintf('Solver converged in %d iterations.\n',i);
else
fprintf('Solver failed to converge in %d iterations!\n',i);
return;
end
s = f(k) - (delta+n)*k - c;
%% Simple sparse diagonalisation
function[x] = spdiag(v)
[rv,cv] = size(v);
if (cv==1)
x = spdiags(v,0,rv,rv);
return;
elseif (rv==1)
x = spdiags(v',0,cv,cv);
return;
else x=v; return;
end
end
0 Comments
Answers (1)
Steven Lord
on 9 Aug 2022
You cannot create a sparse matrix that contains a symbolic expression. Either store the symbolic expression in a symbolic matrix or convert the symbolic expression into a numeric value and store that in the sparse matrix.
syms x
y = x.^2 % symbolic expression
S = sym(zeros(4)) % symbolic matrix
S(2, 3) = y % Assign a symbolic expression into a symbolic matrix
y2 = subs(y, x, 3) % numeric value
A = zeros(4) % numeric matrix
A(2, 3) = y2 % Assign a numeric value into a numeric matrix
You can store a numeric value in a symbolic array (which will automatically convert the numeric value into a symbolic expression) but storing a symbolic expression in a numeric array will not work unless the symbolic expression can be converted into a numeric value.
S(4, 1) = y2 % Storing a numeric value into a symbolic matrix works
A(4, 1) = y % Storing a symbolic expression in a numeric matrix does not work
See Also
Categories
Find more on Number Theory in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!