MATLAB: Fsolve for a generalizable number of equations.

fsolveoptimizationOptimization Toolbox

Hello,
I am attempting to set up an optimization code for solving an economics problem. The problem potentially has many equations and variables but follows a predictable pattern so I'm sure there is a much better way of solving the problem (maybe using a vector or array?).
This is the first problem of this magnitude that I have attempted. I have only used fsolve for single variable iterations prior to this and was unable to find an answer on this forum or Google. Likely due to the fact that I do not know the correct question to ask.
I have a toy example using 4 equations and 4 unknowns below:
% Price equations
p_1 = @(p1,p2,p3,p4,q1,q2,q3,q4) p1 - cost(1,1) + q1./sub(1,1) + q2./sub(2,1) ;
p_2 = @(p1,p2,p3,p4,q1,q2,q3,q4) p2 - cost(2,1) + q2./sub(2,2) + q3./sub(3,2) + q1./sub(1,2) ;
p_3 = @(p1,p2,p3,p4,q1,q2,q3,q4) p3 - cost(3,1) + q3./sub(3,3) + q4./sub(4,3) + q2./sub(2,3) ;
p_4 = @(p1,p2,p3,p4,q1,q2,q3,q4) p4 - cost(4,1) + q4./sub(4,4) + q3./sub(3,4) ;
% Quantity equations
q_1 = @(p1,p2,p3,p4,q1,q2,q3,q4) q1 - delta*((p2 - p1)/(x2 - x1) - (p1 - E)/(x1 - gamma)) ;
q_2 = @(p1,p2,p3,p4,q1,q2,q3,q4) q2 - delta*((p3 - p2)/(x3 - x2) - (p2 - p1)/(x2 - x1)) ;
q_3 = @(p1,p2,p3,p4,q1,q2,q3,q4) q3 - delta*((p4 - p3)/(x4 - x3) - (p3 - p2)/(x3 - x2)) ;
q_4 = @(p1,p2,p3,p4,q1,q2,q3,q4) q4 - delta*(v_max - (p4 - p3)/(x4 - x3)) ;
% Initial guess
w0 = [1000 1000 1000 1000 10000 10000 10000 10000] ;
% Bring price and quantity equations together
solution = @(w) [p_1(w(1),w(2),w(3),w(4),w(5),w(6),w(7),w(8)) ... ;
p_2(w(1),w(2),w(3),w(4),w(5),w(6),w(7),w(8)) ... ;
p_3(w(1),w(2),w(3),w(4),w(5),w(6),w(7),w(8)) ... ;
p_4(w(1),w(2),w(3),w(4),w(5),w(6),w(7),w(8)) ... ;
q_1(w(1),w(2),w(3),w(4),w(5),w(6),w(7),w(8)) ... ;
q_2(w(1),w(2),w(3),w(4),w(5),w(6),w(7),w(8)) ... ;
q_3(w(1),w(2),w(3),w(4),w(5),w(6),w(7),w(8)) ... ;
q_4(w(1),w(2),w(3),w(4),w(5),w(6),w(7),w(8))] ;
% Solve
fsolve(solution,w0)
cost, sub, gamma, E, v_max, x1, x2.x3. and x4 are all known constants. Solving this example works fine but this strategy is clearly infeasible for many products so I was wondering if there were a more efficient way to write this.
Thanks in advance!

Best Answer

A = [1, 0, 0, 0, 1/sub(1,1), 1/sub(2,1), 0, 0;...
0, 1, 0, 0, 1/sub(1,2), 1/sub(2,2), 1/sub(3,2), 0;...
0, 0, 1, 0, 0, 1/sub(2,3), 1/sub(3,3), 1/sub(4,3);...
0, 0, 0, 1, 0, 0, 1/sub(3,4), 1/sub(4,4);...
delta/(x2-x1)+delta/(x1-gamma),-delta/(x2-x1),0, 0, 1, 0, 0, 0;...
-delta/(x2-x1), delta/(x3-x2)+delta/(x2-x1),-delta/(x3-x2),0,0,1,0,0;...
0, -delta/x3-x2),delta/(x4-x3)+delta/(x3-x2),-delta/(x4-x3),0,0,1,0;...
0, 0,-delta/(x4-x3),delta/(x4-x3),0,0,0,0,1];
b = [cost(1,1);cost(2,1);cost(3,1);cost(4,1);delta*E/(x1-gamma);0;0;delta*vmax];
sol = A\b;
p = sol(1:4)
q = sol(5:8)
It's up to you to find the pattern how the matrix A and the vector b are built up.
Best wishes
Torsten.
Related Question