# MATLAB: How to use fsolve within multiple for loops

for loopfsolveinitial guesses

I have a system of 4 non-linear equations with 4 unknowns x(1), x(2), x(3) and x(4). I am able to solve the equations using the following code, however, fsolve only gives me one possible solution. I determined analytically that there is one positive steady state but would like find the approximation to the SS numerically using loops for the initial guesses as well as for different parameter values and are struggling with the code for the for loops. Any help appreciated.
``% The function accepts as input the variables %(x, m0 ,m1, gamma, C0, r1, r2, Uc, p, r4, r5, r6, a, b, r7), % where x = [x(1) x(2) x(3) x(4)] = [C U E M], % and returns vector F(x) as output.%Parameters taken as inputm0=0.1; m1=1; gamma=1; C0=1; r1=0.01;r2=0.01; Uc=0.25; p=1.0; r4=1.0; r5=0.01; r6=1.0;a=0.1; b=0.1;r7=0.1;% Calling fsolvefun = @(x)SteadyStates(x,m0 ,m1, gamma, C0, r1, r2, Uc, p, r4, r5, r6, a, b, r7);x0 = [0.5;0.5;0.2;0.2]; % initial valuesx = fsolve(fun,x0)function F = SteadyStates(x, m0 ,m1, gamma, C0, r1, r2, Uc, p, r4, r5, r6, a, b, r7)C=x(1);U=x(2);E=x(3);M=x(4);F(1) = m0*(1+m1*gamma*x(1))/(1+gamma*x(1))*(C0-x(1))-r1*x(1);F(2) = r1*x(1)+r2*x(2)*(Uc-x(2))-x(2)^2/(p+x(2))+r4*x(3)/(1+x(2))-r5*x(2)-r6*x(2)*x(3);F(3) = a*r1*x(1)+x(2)^2/(p+x(2))-r4*x(3)/(1+x(2))-b*r6*x(2)*x(3);F(4) = b*r6*x(2)*x(3) - r7*x(4);end``

``(4*x2^4 + 4407*x2^3 + 43*x2^2 + 30*x2 - 40)/(4000*(x2^2 + x2 - 1))``