Solve an system of nonlinear equations in for loop

I am trying to solve the system of the three non-linear equations for the different values of parameter A that varies in the iterval of (0; 10), other parameters w, e, d, n, t are fixed. And graphically represent all the solutions (at a certain interval) of each of the dependent variables as a function of parameter A, this is (A, Z) (A, F) and (A, Q). Here is my attempt to code. The problem that the figures are blank, please help me to rectify the failure in the code. Thanks.

      h=0.01;
      A_final=3;
      N=A_final/h;
      A(1)=0;
      syms Q F Z;
      w=pi;
      e=0.5;
      d=10;
      n=1;
      t=0;

      for i=1:N
          A(i+1)=A(i)+h;
          [sol_Q, sol_F, sol_Z]=vpasolve([Z+A(i)*w/(1-Q*(1+e))*cos(w*t+F)==0, cos(2*pi*n*Q+w*t+F)+(e-Q*(1+e)*Z)/(A(i)*w)==0, A(i)-1/(2*w*abs(sin(Q*n*pi)))*(((1-e)*Z)^2+(d*w+2*(1+e)*(Q-1)*Q*n*pi*Z)^2)^0.5==0], [Q, F, Z]);
      end
      figure(1); 
      plot(A,sol_Q)
      xlabel ('A')
      ylabel ('Q')
      figure(2); 
      plot(A,sol_F)
      xlabel ('A')
      ylabel ('F')
      figure(3); 
      plot(A,sol_Z)
      xlabel ('A')
      ylabel ('Z')

Answer 1

[Shantanu Gontia]

There are 3 problems in this piece of code, because of which the plot being generated is empty.

1. In each iteration of the loop, the variables sol_Q, sol_Z and sol_F are being overwritten. To plot the solution, you require the values of these solutions for all the values of A and so you should store the solutions as a vector as well.
2. The function vpasolve returns symbolic variables, while what you want to plot and store are double values. Therefore, an explicit conversion from the symbolic variables to the doubles is required. Additionally, it is possible that the function vpasolve does not return any solution, in which case you can store the solution as NaN if you do not want it plotted.
3. Although, not an error, since the number of values in the vector A is already known, it is more efficient to initialize the vector beforehand.
4. Also, you are calculating the number of points between the numbers 0 and A_final with a step-size of h. The required number of points is actually, N = (A_final-0)/h + 1;. The additional 1 is added because the point 0 is inclusive as well.

To this end I suggest the following changes in the code

% Initialize the vector A beforehand as well as the solution vectors
sol_Q = zeros(1,N);
sol_F = zeros(1,N);
sol_Z = zeros(1,N);
A = 0:h:A_final-h; 
for i=1:N
  [a,b,c]=vpasolve([Z+A(i)*w/(1-Q*(1+e))*cos(w*t+F)==0, cos(2*pi*n*Q+w*t+F)+(e-Q*(1+e)*Z)/(A(i)*w)==0, A(i)-1/(2*w*abs(sin(Q*n*pi)))*(((1-e)*Z)^2+(d*w+2*(1+e)*(Q-1)*Q*n*pi*Z)^2)^0.5==0], [Q, F, Z]);
  % Explicitly convert the symbolic variables to doubles and store them
  % in an vector. 
  if (~isempty(a))
    sol_Q(i) = double(a);
  else
    % If the the solution returned is empty, store NaN to avoid
    % plotting it.
    sol_Q(i) = NaN;  
  end
  
  if (~isempty(b))
    sol_F(i) = double(b);
  else
    sol_F(i) = NaN;
  end
  
  if (~isempty(b))
    sol_Z(i) = double(c);
  else
    sol_Z(i) = NaN;
  end
end