title: MATLAB beginning having difficulty with function handles
body: <p>I'm working on a project for school that basically involves iteratively solving a cubic equation. I'm using MATLAB for it, but I've never really done much with MATLAB so I'm having some trouble with the logic of it.</p>
  
  <p>Here's my code:</p>
  
  <pre><code>% Redlich/Kwong EOS
  sigma = 1;
  epsilon = 0;
  omega = 0.08664;
  psi = 0.42748;
  
  beta = @(Psat_RK) omega*PsatRK/Pc/Tr(1); % Pc is in bar, may need a unit conversion later
  alpha = (Tr(1))^(-1/2);
  q = psi*alpha/omega/Tr(1);
  
  A = @(beta) (sigma + epsilon - 1)*beta - 1;
  B = @(beta) (sigma*epsilon - sigma - epsilon)*(beta^2) + (q - sigma - epsilon)*beta;
  C = @(beta) -(sigma*epsilon*(1+beta) + q)*(beta^2);
  
  Q = @(A,B)   ((A^2) - 3*B)/9;
  R = @(A,B,C) (2*(A^3) - 9*A*B + 27*C)/54;
  M = @(R,Q)   R^2 - Q^3;
  
  if M &gt; 0
      Z_single = ((-R+(M^0.5))^(1/3)) + ((-R-(M^0.5))^(1/3)) - (A/3);
  
      I = (1/(sigma-epsilon))*log((Z_single+sigma*beta)/(Z_single+epsilon*beta));
  end
  
  if M &lt; 0
     theta = acos(R/(Q^(1/3)));
     Z(1) = -2*(Q^0.5)*(cos(theta/3)) - (A/3);
     Z(2) = -2*(Q^0.5)*(cos((theta + 2*pi)/3)) - (A/3);
     Z(3) = -2*(Q^0.5)*(cos((theta - 2*pi)/3)) - (A/3);
  
     Z_liquid = min(Z)
     Z_vapor = max(Z)
  
     I_liquid = (1/(sigma-epsilon))*log((Z_liquid+sigma*beta)/(Z_liquid+epsilon*beta));
     I_vapor = (1/(sigma-epsilon))*log((Z_vapor+sigma*beta)/(Z_vapor+epsilon*beta));
  end
  
  ln_phi_liquid = Z_liquid - 1 - log(Z_liquid - beta) - q*I_liquid;
  ln_phi_vapor = Z_vapor - 1 - log(Z_vapor - beta) - q*I_vapor;
  
  objfun = (ln_phi_liquid - ln_phi_vapor);
  
  Psat_RK_solved = fsolve(objfun,10);
  </code></pre>
  
  <p>Basically, I'm trying to iterate on the value of <code>Psat_RK</code> until the <code>objfun</code> equals 0. I can post more details of the math if needed, but I figured this would be enough to get started. Thanks.</p>
  
  <p>Edit: Sorry, forgot to actually mention the problem.</p>
  
  <p>Here's the error I'm getting.</p>
  
  <blockquote>
    <p>Undefined operator '>' for input arguments of type 'function_handle'.</p>
    
    <p>Error in Proj2 (line 73)
    if M > 0</p>
  </blockquote>
  
  <p>I can't figure out how to establish in this line that M is being calculated off an anonymous function.</p>
  
  <hr>
  
  <p>EDIT: </p>
  
  <pre><code>sigma = 1;
  epsilon = 0;
  omega = 0.08664;
  psi = 0.42748;
  
  beta = @(Psat_RK) omega*PsatRK/Pc/Tr(1); % Pc is in bar, may need a unit conversion later
  alpha = (Tr(1))^(-1/2);
  q = psi*alpha/omega/Tr(1);
  
  A = @(Psat_RK) (sigma + epsilon - 1)*beta(Psat_RK) - 1;
  B = @(Psat_RK) (sigma*epsilon - sigma - epsilon)*(beta(Psat_RK)^2) + (q - sigma - epsilon)*beta(Psat_RK);
  C = @(Psat_RK) -(sigma*epsilon*(1+beta(Psat_RK)) + q)*(beta(Psat_RK)^2);
  
  Q = @(Psat_RK) ((A(Psat_RK)^2) - 3*B(Psat_RK))/9;
  R = @(Psat_RK) (2*(A(Psat_RK)^3) - 9*A(Psat_RK)*B(Psat_RK) + 27*C(Psat_RK))/54;
  M = @(Psat_RK) R(Psat_RK)^2 - Q(Psat_RK)^3;
  
  if M(Psat_RK) &gt; 0
      Z_single = ((-R+(M^0.5))^(1/3)) + ((-R-(M^0.5))^(1/3)) - (A/3);
  
      I = (1/(sigma-epsilon))*log((Z_single+sigma*beta)/(Z_single+epsilon*beta));
  end
  
  if M(Psat_RK) &lt; 0
     theta = acos(R/(Q^(1/3)));
     Z(1) = -2*(Q^0.5)*(cos(theta/3)) - (A/3);
     Z(2) = -2*(Q^0.5)*(cos((theta + 2*pi)/3)) - (A/3);
     Z(3) = -2*(Q^0.5)*(cos((theta - 2*pi)/3)) - (A/3);
  
     Z_liquid = min(Z)
     Z_vapor = max(Z)
  
     I_liquid = (1/(sigma-epsilon))*log((Z_liquid+sigma*beta)/(Z_liquid+epsilon*beta));
     I_vapor = (1/(sigma-epsilon))*log((Z_vapor+sigma*beta)/(Z_vapor+epsilon*beta));
  end
  
  ln_phi_liquid = Z_liquid - 1 - log(Z_liquid - beta) - q*I_liquid;
  ln_phi_vapor = Z_vapor - 1 - log(Z_vapor - beta) - q*I_vapor;
  
  objfun = (ln_phi_liquid - ln_phi_vapor);
  
  Psat_RK_solved = fsolve(objfun,10);
  </code></pre>
  
  <p>I know the code needs some work further down, but the code below shouldn't affect why it hangs at the first <code>if</code> statement, right? The error is:</p>
  
  <p>Undefined function or variable 'Psat_RK'.</p>
  
  <p>Error in Proj2 (line 122)
  if M(Psat_RK) > 0</p>