MATLAB: Index exceeds matrix dimensions, using fsolve

fsolvematrix manipulation

I don't see any reason why I should be exceeding matrix dimensions here. The constants in the function (z onwards till Dinf) are defined in a separate script,which I 'run' just before executing the function.
function [solution]=nano2(x,z,PHI,gamma,Kc,Dp,PHI_B,Mc,N,N_c,deltax_j,J_v,f,gc,T,Dinf)
A=x(1:3,1:N+2); %Number of rows to be modified with number of components
B=x(4:6,1:N+2);
Q=x(7:9,1:N+2);
R=x(10:12,1:N+2);
E=x(13:15,1:N+2);
F=x(16:18,1:N+3); %Might need to look at its dimensions.
PSI=x(19,1:N+3);
c=x(20:22,1:N+3);
zeta=x(23,1);
for c1=1:N_c
for c2=3:N+2
solution=[A(c1,c2)-((Dp(c1)/deltax_j)-(1/2)*Kc(c1)*J_v+(1/2)*z(c1)*Dp(c1)*(f/(gc*T))*((PSI(c2+1)-PSI(c2))/deltax_j));%Membrane active layer domain for ion i-linearization of Nernst-Planck equation--------(25-31)
B(c1,c2)+(Dp(c1)/deltax_j)-(1/2)*Kc(c1)*J_v+(1/2)*z(c1)*Dp(c1)*(f/(gc*T))*((PSI(c2+1)-PSI(c2))/deltax_j);
Q(c1,c2)-(1/2)*z(c1)*Dp(c1)*(f/(gc*T))*((c(c1,c2+1)+c(c1,c2))/deltax_j);
R(c1,c2)+(1/2)*z(c1)*Dp(c1)*((f/(gc*T)))*((c(c1,c2+1)+c(c1,c2))/deltax_j);
E(c1,c2)+J_v;
F(c1,c2)+(1/2)*z(c1)*Dp(c1)*(f/(gc*T))*(c(c1,c2)+c(c1,c2+1))*((PSI(c2+1)-PSI(c2))/deltax_j);
F(c1,c2)-A(c1,c2)*c(c1,c2)-B(c1,c2)*c(c1,c2+1)-Q(c1,c2)*PSI(c2)-R(c1,c2)*PSI(c2+1)-E(c1,c2)*c(c1,N+3);
% Feed/solution equilibrium
%Feed solution/membrane mass transfer coefficient for each component------(22)%
(Mc(c1)-J_v)*c(c1,2)+J_v*c(c1,N+3)+(z(c1)*c(c1,2)*Dinf(c1)*(f/(gc*T)))*zeta-Mc(c1)*c(c1,1);
%Thermodynamic equilibrium conditions at feed solution/membrane for each component--------(24)
(gamma(c1,3)*c(c1,3))+(gamma(c1,2)*c(c1,2)*z(c1)*(f/(gc*T))*PHI(c1)*PHI_B(c1)*exp((-z(c1)*f*PSI(3))/(gc*T))*PSI(3))-(gamma(c1,2)*c(c1,2)*PHI(c1)*PHI_B(c1)*exp((-z(c1)*f*PSI(3))/(gc*T))*(1+(((z(c1)*f*PSI(3))/(gc*T)))));
%Thermodynamic equilibrium conditions at membrane/permeate-solution interface for each component----(33)
(gamma(c1,N+2)*c(c1,N+2))-(PHI(c1)*exp(z(c1)*(f/(gc*T))*(PSI(N+3)-PSI(N+2))))*c(c1,N+3)+(gamma(c1,N+3)*c(c1,N+3)*z(c1)*(f/(gc*T))*PHI(c1)*exp(z(c1)*(f/(gc*T))*(PSI(N+3)-PSI(N+2))))*PSI(N+2)-(gamma(c1,N+3)*c(c1,N+3)*z(c1)*(f/(gc*T))*PHI(c1)*exp(z(c1)*(f/(gc*T))*(PSI(N+3)-PSI(N+2))))*PSI(N+3)-(gamma(c1,N+3)*c(c1,N+3)*PHI(c1)*exp((z(c1)*(f/(gc*T)))*(PSI(N+3)-PSI(N+2))))*(z(c1)*(f/(gc*T))*(PSI(N+3)-PSI(N+2)));
%Electroneutrality at each node--------(32)"
z(1)*c(1,c2)+z(2)*c(2,c2)+c_x;
%Electroneutrality in boundary layer on feed solution/membrane------(23)
z(1)*c(1,2)+z(2)*c(2,2)+z(3)*c(3,2);
%Electroneutrality in boundary layer on permeate solution/membrane------(34)
z(1)*c(1,N+3)+z(2)*c(2,N+3)+z(3)*c(3,N+3);
];
end
end
end

Best Answer

  • I suggest you not worry about fsolve just now. Call your nano2 function directly with your initial parameter estimates for it, then see if you can isolate the statement that is likely throwing the error. Also, look at the sizes of A through zeta and check c1 and c2.
    The debugger will likely be a significant help in code as complex as this.