Save each loop values, without overwriting

Hi. I think that my program is not saving each loop value for the Matrix Qabal which comes from Q, and for QQ neither, but I just want to solve the Q or Qabal, which it should be a matrix (22x24)

MND is a matrix (14x25) and MTB(which is inside the function has 22 rows).

for j=2:25
  QQ=MND(:,j); %m3/s Vector containing fixed nodal discharges with each PF
  [Y,Df,Q,VK1,Vn1]=funF2(H,QQ,MTB,VK,Vn,Vii,g,vi);
  er=sum(abs(Y)); %Initial guess solution error
  nit=0; %Iteration counter initiation
  %Iterative solution search 
  while er>=eps && nit<=nma
    VK=VK1;
    Vn=Vn1;
    H(Vii,j)=H(Vii)-Df\Y; %New extended H vector (solution vector)
    [Y,Df,Q,VK1,Vn1,Qij]=funF2(H,QQ,MTB,VK,Vn,Vii,g,vi);
    Qabal(Hores,j)=Q(Hores);
    er=sum(abs(Y)); %New approximate solution error
    nit=nit+1; %Iteration counter update
    %Hores=Hores+1;
    HH=H(:,2:end);     %Only found H columns
  end
end

I don't know If i have explained it well, maybe is needed the other part of the program to find it out, but I think that with that should be enough.

Thank you.

Answer 1

[dpb]

for j=2:25
  QQ=MND(:,j-1); %m3/s Vector containing fixed nodal discharges with each PF
  [Y,Df,Q,VK1,Vn1]=funF2(H,QQ,MTB,VK,Vn,Vii,g,vi);

The above will overwrite Q each pass with the results for that output argument from funF2 given the input using the j column of MND and the other inputs. We don't know what the dimensions are, presuming it is also a vector and you wish to save the various columns into an array Q instead, you would need

[Y,Df,Q(:,j),VK1,Vn1]=funF2(H,QQ,MTB,VK,Vn,Vii,g,vi);

instead and if that is so, you should "preallocate" Q as

Q=zeros(size(MTB,1),size(MND,1)-1);  % preallocate

Similar holds for the other output array excepting you've written

   Qabal(Hores,j)=Q(Hores);

which will load individual elements but you commented out the defining line that incremented the row index Hores. But, you can only have a single element to assign to a specific array location and it appears your function returns a vector so you'll need to (it would seem) also use

   Qabal(:,j)=Q(Hores);

ERRATUM

Dang, missed an edit when pasted the line...intended to write

   Qabal(:,j)=Q;

Of course, this raises the question of just how the iteration is supposed to be working--looks like maybe it's trying to update individual locations somehow, not the full vector but we don't know the intent just from code (particularly if the code doesn't work as intended :) ).

[Jordi Escola]

I tried this before already but there is no way i can solve it. We are a group trying it but no solution...

Here is the "whole program", if you have time to check it, but don't worry if not, its quite long...

 %node/length/diameter/e
 MTB=[13,1,1000,0.380,0.0001,0,0;
    1,2,1200,0.305,0.0001,0,0;
    2,3,800,0.305,0.0001,0,0;
    3,4,1200,0.255,0.0001,0,0;
    2,5,1000,0.255,0.0001,0,0;
    1,6,1200,0.255,0.0001,0,0;
    1,7,1200,0.305,0.0001,0,0;
    4,5,1000,0.205,0.0001,0,0;
    5,6,1100,0.205,0.0001,0,0;
    4,8,1200,0.255,0.0001,0,0;
    5,9,1100,0.205,0.0001,0,0;
    6,9,2000,0.205,0.0001,0,0;
    6,10,1200,0.255,0.0001,0,0;
    6,11,2000,0.255,0.0001,0,0;
    7,11,2000,0.255,0.0001,0,0;
    11,10,800,0.150,0.0001,0,0;
    9,10,1000,0.205,0.0001,0,0;
    8,9,800,0.150,0.0001,0,0;
    8,12,1100,0.205,0.0001,0,0;
    9,12,1000,0.205,0.0001,0,0;
    10,12,2200,0.205,0.0001,0,0;
    12,14,800,0.305,0.0001,0,0];

 %height and Qout
 MND=[Matrix 14x25, where i wanna take just columns 2 to 25]

 Vii=[1,2,3,4,5,6,7,8,9,10,11,12]'; %Vector of nodes with unknown Hi

 g=9.81;                          %N/kg Acceleration of gravity
 vi=1.1e-6;                       %m2/s Fluid kinematic viscosity

 Hores=[Vector of 24];

 %Definition of computational parameters
 eps=1e-6; %Maximum allowed error in the Fi
 nma=30; %Maximum allowed iterative steps

 %Definition of intermediate variables

 nt=length(MTB(:,1)); %Number of pipes in the network
 nn=length(MND(:,1)); %Total number of nodes in the network
 ni=length(Vii); %Number of nodes with unknown Hi
 npf=length(PF);

 VK=zeros(nt,npf); %Vector containing the K value of each pipe
 Vn=zeros(nt,npf); %Vector containing the n value of each pipe

 %NEWTONS ITERATIVE METHOD
 %Definition of the starting point
 H=MND(:,1); %m Initial guess vector for H-eqs method (Hi=zi)

 %Initial approximation to K and n of each pipe
 f=0.02; %Initial guess for the friction factor
 n=1.95; %Initial guess for exponent n

 Q=zeros(size(MTB,1),size(MND,1)-1);  % preallocate

 for l=1:nt
     if MTB(l,7)~=0 %Pipe element with fixed K and n
     VK(l)=MTB(l,6);
     Vn(l)=MTB(l,7);
     else
     VK(l)=8*f*MTB(l,3)/g/pi^2/MTB(l,4)^5; %Initial guess for K (eq. 19)
     Vn(l)=n; %Initial guess for n
     end
 end
    for j=2:25
  QQ=MND(:,j-1); %m3/s Vector containing fixed nodal discharges with    

    [Y,Df,Q,VK1,Vn1]=funF2(H,QQ,MTB,VK,Vn,Vii,g,vi);
    er=sum(abs(Y)); %Initial guess solution error
    nit=0; %Iteration counter initiation

        %Iterative solution search 
        while er>=eps && nit<=nma

             VK=VK1;
             Vn=Vn1;
             H(Vii,j)=H(Vii)-Df\Y; %New extended H vector (solution
             [Y,Df,Q(:,j),VK1,Vn1,Qij]=funF2(H,QQ,MTB,VK,Vn,Vii,g,vi);
   Qabal(:,j)=Q(Hores);

             er=sum(abs(Y)); %New approximate solution error
             nit=nit+1; %Iteration counter update
             %Hores=Hores+1;

    HH=H(:,2:end);     %Only found H columns

        end
    end

And the Function funF2:

 %Auxiliar variables
 nn=length(H); %Number of nodes in the network
 nt=length(MTB(:,1)); %Number of pipes in the network
 Y=zeros(nn,24); %Column vector containing all Fi functions
 Df=zeros(nn); %Extended jacobian matrix initilization
 Q=zeros(nt,24); %Column vector to contain the flow of each pipe
 VK1=zeros(nt,24); %Column vector to contain new K values
 Vn1=zeros(nt,24); %Column vector to contain new n values
 %Iterative proces over all pipes
 for l=1:nt
     i=MTB(l,1);
     j=MTB(l,2);
     L=MTB(l,3); %m Pipe length
     D=MTB(l,4); %m Pipe inner diameter
     K=VK(l);
     n=Vn(l);
     Qij=sign(H(i)-H(j))*abs((H(i)-H(j))/K)^(1/n); %m3/s Flow from node i to node j
     Q(l)=Qij;
     %Computing the extended Fi terms
     Y(i)=Y(i)+Qij;
     Y(j)=Y(j)-Qij;
     %Computing the extended jacobian matrix terms
     Reij=4*abs(Qij)/pi/D/vi; %The actual Reynolds number
         if Reij<=2300
         dQdH=-1/K; %Laminar flow (eq.28)
         else
         dQdH=-abs((H(i)-H(j))/K)^(1/n-1)/n/K; %Turbulent flow (eq.26)
         end
 Df(i,j)=dQdH;
 Df(j,i)=dQdH;
 Df(i,i)=Df(i,i)-dQdH;
 Df(j,j)=Df(j,j)-dQdH;
%Update of K and n for use in the next iteration
   if MTB(l,7)~=0 %Pipe element with fixed K and n
   VK1(l)=MTB(l,6);
   Vn1(l)=MTB(l,7);
   elseif Reij<=2300 %LAMINAR FLOW
          VK1(l)=128*vi*L/pi/g/D^4; %Equation 17
     Vn1(l)=1; %Equation 18
     else %TURBULENT FLOW
     Re=[0.95,1.05]*Reij; %The [Re1, Re2] points
     f=1.325./log(MTB(l,5)/3.71/D+5.74./Re.^.9).^2; %Corresponding [f1, f2] points
     b=log(f(1)/f(2))/log(Re(2)/Re(1)); %Equation 11
     Q1=pi*D*vi*Re(1)/4;
     a=f(1)*Q1^b; %Equation 12
     VK1(l)=8*a*L/g/pi^2/D^5; %Equation 13
     Vn1(l)=2-b; %Equation 14
     end
 end
 Y=Y(Vii)+QQ(Vii); %m3/s Vector of Fi values for the ni nodes with unknown Hi
 Df=Df(Vii,Vii); %The jacobian matrix

Thanks anyway, for the time that it took.

[dpb]

Not sure what "this" refers to in the above but see the ERRATUM and comment following it to my previous Answer. Maybe the question raised there will lead to some enlightenment to resolve the issues.

As stated there, don't know any of us here can help without a more definitive statement of the problem trying to be solved, not just the code and a syntax issue with it.