Undefined function 'minus' for input arguments of type 'function_handle'.

m_tot = 8; % in kg/sec
[M2,H] = B249VALVEcurve()
C_v249= [M2];
[M,H] = B251valveperfcurve()
C_v251 = [M];% a vector for opening from 1 - 100%
wInP = 400; % in kPa
k = 1;
PD_pipe2 = 0.96 % In kPa from toolbox online - assuming 3 inch dia, 20 ft length, flow rate for 5 kg/sec
PD_DIC1  = @(m_w2)(160*(m_w2 ^2) + (2.3*10^3)*m_w2 - 10^4)/1e3;
P_out_DIC1 = @(m_w2) wInP - PD_DIC1 - PD_pipe2;
rho_w =refpropm('D','T',345,'P',400,'water');
Q_w2 = @(m_w2)m_w2/rho_w;
Q_w2_hr =@(m_w2) 3600*Q_w2; % in m3/hr
m_w1 = @(m_w2) m_tot - m_w2;
PD_pipe1 = 0.96;
PD_oc = 5;
PD_IC1= @(m_w2) (590*(m_w1 ^2) + 0.42*m_w1+ 0.56)/1e3;
PD_IC2 =@(m_w2) (590*(m_w1 ^2) + 0.4*m_w1+ 0.51)/1e3;
PD_AC=  @(m_w2)(1400*(m_w1^2) + 0.37*m_w1 + 0.44)/1e3;
P_out_AC =@(m_w2) wInP - PD_IC1 - PD_pipe2 - PD_IC2 - PD_oc - PD_AC;
rho_w1 =refpropm('D','T',360,'P',280,'water'); % assuming 280, anyways rho shouldnt vary much
Q_w1 = m_tot/rho_w1;
Q_w1_hr= 3600*Q_w1; % in m3/hr
for i = 1:100
    Kv(i) = 0.865*C_v249(i); % gpm into m3/hr
    DP_B249(i) =  @(m_w2)((Q_w2_hr/Kv(i))^2)*(10^2); % Q is in m3/hr, and Press is in kPa
    DP_B249a = DP_B249(i);
    for j = 1:100
        Kv_251(j) = 0.865*C_v251(j);
        %         Kv_251(i,j) = Kv_251(j);
        %         DP_B249(i,j) = DP_B249(i) ;
        %         PD_DIC1(i,j) = PD_DIC1(i);
        DP_B251(j) = ((Q_w1_hr/Kv_251(j))^2)*(10^2); % Q is in m3/hr, and Press is in kPa
        DP_B251a = DP_B251(j);

          PD_1 = @(m_w2)PD_pipe1  +PD_oc  +PD_IC1  +PD_IC2 + PD_AC + DP_B251a ;
          PD_2 = @(m_w2)PD_pipe2 + PD_DIC1 + DP_B249a+ PD_AC+ DP_B251a ;
          fun = @(m_w2) (PD_2 - PD_1);
          m_w2 = fzero(fun,2);
          m_w2_f(k) = m_w2;
          m_w1_f(k) = m_tot - m_w2;
          B251valve_pos(k) = j; %valve position of 251 for this flow rate
          B249valve_pos(k) = i; % valve position of 249 for this flow rate
          k =  k+1;
          j = j+1;
      end
      i = i +1;
  end

[Walter Roberson]

            fun = @(m_w2) (PD_2(m_w2) - PD_1(m_w2));

[Roshni Khetan]

Hi Walter! I did this but it gives me an error:

Index exceeds matrix dimensions for this line: m_w2 = fzero(fun,2);

m_tot = 8; % in kg/sec RUN THE MODEL FOR DIFFERENT m_tot
[M2,H] = B249VALVEcurve()
C_v249= [M2];
[M,H] = B251valveperfcurve()
C_v251 = [M];% a vector for opening from 1 - 100%
wInP = 400; % in kPa
k = 1;
PD_pipe2 = 0.96 % In kPa from toolbox online - assuming 3 inch dia, 20 ft length, flow rate for 5 kg/sec
PD_DIC1  = @(m_w2)(160*(m_w2 ^2) + (2.3*10^3)*m_w2 - 10^4)/1e3;
rho_w =refpropm('D','T',345,'P',400,'water');
Q_w2 = @(m_w2)m_w2/rho_w;
Q_w2_hr =@(m_w2) 3600*Q_w2(m_w2); % in m3/hr
m_w1 = @(m_w2) m_tot - m_w2;
PD_pipe1 = 0.96;
PD_oc = 5;
PD_IC1= @(m_w2) (590*(m_w1(m_w2) ^2) + 0.42*m_w1(m_w2)+ 0.56)/1e3;
PD_IC2 =@(m_w2) (590*(m_w1(m_w2) ^2) + 0.4*m_w1(m_w2)+ 0.51)/1e3;
PD_AC=  @(m_w2)(1400*(m_w1(m_w2)^2) + 0.37*m_w1(m_w2) + 0.44)/1e3;
rho_w1 =refpropm('D','T',360,'P',280,'water'); % assuming 280, anyways rho shouldnt vary much
Q_w1 = m_tot/rho_w1;
Q_w1_hr= 3600*Q_w1; % in m3/hr
for i = 1:100
    Kv(i) = 0.865*C_v249(i); % gpm into m3/hr
    DP_B249(i) =  @(m_w2)((Q_w2_hr(m_w2)/Kv(i))^2)*(10^2); % Q is in m3/hr, and Press is in kPa
    DP_B249a = DP_B249(i);
    for j = 1:100
        Kv_251(j) = 0.865*C_v251(j);
        %         Kv_251(i,j) = Kv_251(j);
        %         DP_B249(i,j) = DP_B249(i) ;
        %         PD_DIC1(i,j) = PD_DIC1(i);
        DP_B251(j) = ((Q_w1_hr/Kv_251(j))^2)*(10^2); % Q is in m3/hr, and Press is in kPa
        DP_B251a = DP_B251(j);

          PD_1 = @(m_w2)PD_pipe1  +PD_oc  +PD_IC1(m_w2)  +PD_IC2(m_w2) + PD_AC + DP_B251a ;
          PD_2 = @(m_w2)PD_pipe2 + PD_DIC1(m_w2) + DP_B249a(m_w2)+ PD_AC+ DP_B251a ;
          fun = @(m_w2) (PD_2(m_w2)- PD_1(m_w2));
          m_w2 = fzero(@ (m_w2) fun(m_w2),4);
          m_w2_f(k) = m_w2;
          m_w1_f(k) = m_tot - m_w2;
          B251valve_pos(k) = j; %valve position of 251 for this flow rate
          B249valve_pos(k) = i; % valve position of 249 for this flow rate
          k =  k+1;
          j = j+1;
      end
      i = i +1;
  end

Answer 1

[Star Strider]

My guess is that this assignment (and others like it in your code) are throwing the error:

P_out_AC =@(m_w2) wInP - PD_IC1 - PD_pipe2 - PD_IC2 - PD_oc - PD_AC;

The solution is to provide arguments to your other functions in that assignment:

P_out_AC =@(m_w2) wInP - PD_IC1(m_w2) - PD_pipe2 - PD_IC2(m_w2) - PD_oc - PD_AC(m_w2);

Now they return values, and should work to calculate ‘P_out_AC’.

[Roshni Khetan]

Thank you for the response! I tried this but now it says:

Index exceeds matrix dimensions for this line: m_w2 = fzero(fun,2);

I dont understand why

%load flow characteristic curves
m_tot = 8; % in kg/sec RUN THE MODEL FOR DIFFERENT m_tot
[M2,H] = B249VALVEcurve()
C_v249= [M2];
[M,H] = B251valveperfcurve()
C_v251 = [M];% a vector for opening from 1 - 100%
wInP = 400; % in kPa
k = 1;
PD_pipe2 = 0.96 % In kPa from toolbox online - assuming 3 inch dia, 20 ft length, flow rate for 5 kg/sec
PD_DIC1  = @(m_w2)(160*(m_w2 ^2) + (2.3*10^3)*m_w2 - 10^4)/1e3;
rho_w =refpropm('D','T',345,'P',400,'water');
Q_w2 = @(m_w2)m_w2/rho_w;
Q_w2_hr =@(m_w2) 3600*Q_w2(m_w2); % in m3/hr
m_w1 = @(m_w2) m_tot - m_w2;
PD_pipe1 = 0.96;
PD_oc = 5;
PD_IC1= @(m_w2) (590*(m_w1(m_w2) ^2) + 0.42*m_w1(m_w2)+ 0.56)/1e3;
PD_IC2 =@(m_w2) (590*(m_w1(m_w2) ^2) + 0.4*m_w1(m_w2)+ 0.51)/1e3;
PD_AC=  @(m_w2)(1400*(m_w1(m_w2)^2) + 0.37*m_w1(m_w2) + 0.44)/1e3;
rho_w1 =refpropm('D','T',360,'P',280,'water'); % assuming 280, anyways rho shouldnt vary much
Q_w1 = m_tot/rho_w1;
Q_w1_hr= 3600*Q_w1; % in m3/hr
for i = 1:100
    Kv(i) = 0.865*C_v249(i); % gpm into m3/hr
    DP_B249(i) =  @(m_w2)((Q_w2_hr(m_w2)/Kv(i))^2)*(10^2); % Q is in m3/hr, and Press is in kPa
    DP_B249a = DP_B249(i);
    for j = 1:100
        Kv_251(j) = 0.865*C_v251(j);
        %         Kv_251(i,j) = Kv_251(j);
        %         DP_B249(i,j) = DP_B249(i) ;
        %         PD_DIC1(i,j) = PD_DIC1(i);
        DP_B251(j) = ((Q_w1_hr/Kv_251(j))^2)*(10^2); % Q is in m3/hr, and Press is in kPa
        DP_B251a = DP_B251(j);

          PD_1 = @(m_w2)PD_pipe1  +PD_oc  +PD_IC1(m_w2)  +PD_IC2(m_w2) + PD_AC + DP_B251a ;
          PD_2 = @(m_w2)PD_pipe2 + PD_DIC1(m_w2) + DP_B249a(m_w2)+ PD_AC+ DP_B251a ;
          fun = @(m_w2) (PD_2(m_w2)- PD_1(m_w2));
          m_w2 = fzero(fun,2);
          m_w2_f(k) = m_w2;
          m_w1_f(k) = m_tot - m_w2;
          B251valve_pos(k) = j; %valve position of 251 for this flow rate
          B249valve_pos(k) = i; % valve position of 249 for this flow rate
          k =  k+1;
          j = j+1;
      end
      i = i +1;
  end

[Star Strider]

You need to check the dimensions of your arrays. I cannot determine that from the information you posted.

[Walter Roberson]

At the command line give the command

dbstop if caught error

and then run your program. If it stops for any reason other than the "index exceeds" error, give the command "dbcont". When it stops because of the "index exceeds" error, examine the size of your matrices.