I feel difficulty to find the roots of peng robinson eqution given below ..can any body tell me how to find the roots with respect to volume.
    5 views (last 30 days)
  
       Show older comments
    
P *V^3 + (b*P - R*T) *V ^2 - (3*b^2P + 2*b*R*T - a)*V + (b^3P + R*T*b^2 - a*b) = 0
Answers (1)
  John D'Errico
      
      
 on 5 Sep 2016
        If the constants are known, then just use roots! WTP?
If the constants are not known, then it will get messy, since you have a non-constant coefficient cubic polynomial. Using solve will do it though, since the degree is low enough. It will look a bit messy though.
syms P V b R T a
vpa(solve(P *V^3 + (b*P - R*T) *V ^2 - (3*b^2*P + 2*b*R*T - a)*V + (b^3*P + R*T*b^2 - a*b) == 0,V),16)
ans =
                                                             ((0.1111111111111111*(R*T - 1.0*P*b)^2)/P^2 + (0.3333333333333333*(3.0*P*b^2 + 2.0*R*T*b - 1.0*a))/P)/((0.03703703703703704*(R*T - 1.0*P*b)^3)/P^3 + (((0.03703703703703704*(R*T - 1.0*P*b)^3)/P^3 - (0.5*P*b^3 + 0.5*R*T*b^2 - 0.5*a*b)/P + (0.1666666666666667*(R*T - 1.0*P*b)*(3.0*P*b^2 + 2.0*R*T*b - 1.0*a))/P^2)^2 - 1.0*((0.1111111111111111*(R*T - 1.0*P*b)^2)/P^2 + (1.0*P*b^2 + 0.6666666666666667*R*T*b - 0.3333333333333333*a)/P)^3)^(1/2) - (0.5*P*b^3 + 0.5*R*T*b^2 - 0.5*a*b)/P + (0.1666666666666667*(R*T - 1.0*P*b)*(3.0*P*b^2 + 2.0*R*T*b - 1.0*a))/P^2)^(1/3) + ((0.03703703703703704*(R*T - 1.0*P*b)^3)/P^3 - (0.5*(P*b^3 + R*T*b^2 - 1.0*a*b))/P + (((0.03703703703703704*(R*T - 1.0*P*b)^3)/P^3 - (0.5*(P*b^3 + R*T*b^2 - 1.0*a*b))/P + (0.1666666666666667*(R*T - 1.0*P*b)*(3.0*P*b^2 + 2.0*R*T*b - 1.0*a))/P^2)^2 - 1.0*((0.1111111111111111*(R*T - 1.0*P*b)^2)/P^2 + (0.3333333333333333*(3.0*P*b^2 + 2.0*R*T*b - 1.0*a))/P)^3)^(1/2) + (0.1666666666666667*(R*T - 1.0*P*b)*(3.0*P*b^2 + 2.0*R*T*b - 1.0*a))/P^2)^(1/3) + (0.3333333333333333*(R*T - 1.0*P*b))/P
 (((0.1111111111111111*(R*T - 1.0*P*b)^2)/P^2 + (0.3333333333333333*(3.0*P*b^2 + 2.0*R*T*b - 1.0*a))/P)*(- 0.5 - 0.8660254037844386i))/((0.03703703703703704*(R*T - 1.0*P*b)^3)/P^3 + (((0.03703703703703704*(R*T - 1.0*P*b)^3)/P^3 - (0.5*P*b^3 + 0.5*R*T*b^2 - 0.5*a*b)/P + (0.1666666666666667*(R*T - 1.0*P*b)*(3.0*P*b^2 + 2.0*R*T*b - 1.0*a))/P^2)^2 - 1.0*((0.1111111111111111*(R*T - 1.0*P*b)^2)/P^2 + (1.0*P*b^2 + 0.6666666666666667*R*T*b - 0.3333333333333333*a)/P)^3)^(1/2) - (0.5*P*b^3 + 0.5*R*T*b^2 - 0.5*a*b)/P + (0.1666666666666667*(R*T - 1.0*P*b)*(3.0*P*b^2 + 2.0*R*T*b - 1.0*a))/P^2)^(1/3) - ((0.03703703703703704*(R*T - 1.0*P*b)^3)/P^3 - (0.5*(P*b^3 + R*T*b^2 - 1.0*a*b))/P + (((0.03703703703703704*(R*T - 1.0*P*b)^3)/P^3 - (0.5*(P*b^3 + R*T*b^2 - 1.0*a*b))/P + (0.1666666666666667*(R*T - 1.0*P*b)*(3.0*P*b^2 + 2.0*R*T*b - 1.0*a))/P^2)^2 - 1.0*((0.1111111111111111*(R*T - 1.0*P*b)^2)/P^2 + (0.3333333333333333*(3.0*P*b^2 + 2.0*R*T*b - 1.0*a))/P)^3)^(1/2) + (0.1666666666666667*(R*T - 1.0*P*b)*(3.0*P*b^2 + 2.0*R*T*b - 1.0*a))/P^2)^(1/3)*(0.5 - 0.8660254037844386i) + (0.3333333333333333*(R*T - 1.0*P*b))/P
 (((0.1111111111111111*(R*T - 1.0*P*b)^2)/P^2 + (0.3333333333333333*(3.0*P*b^2 + 2.0*R*T*b - 1.0*a))/P)*(- 0.5 + 0.8660254037844386i))/((0.03703703703703704*(R*T - 1.0*P*b)^3)/P^3 + (((0.03703703703703704*(R*T - 1.0*P*b)^3)/P^3 - (0.5*P*b^3 + 0.5*R*T*b^2 - 0.5*a*b)/P + (0.1666666666666667*(R*T - 1.0*P*b)*(3.0*P*b^2 + 2.0*R*T*b - 1.0*a))/P^2)^2 - 1.0*((0.1111111111111111*(R*T - 1.0*P*b)^2)/P^2 + (1.0*P*b^2 + 0.6666666666666667*R*T*b - 0.3333333333333333*a)/P)^3)^(1/2) - (0.5*P*b^3 + 0.5*R*T*b^2 - 0.5*a*b)/P + (0.1666666666666667*(R*T - 1.0*P*b)*(3.0*P*b^2 + 2.0*R*T*b - 1.0*a))/P^2)^(1/3) - ((0.03703703703703704*(R*T - 1.0*P*b)^3)/P^3 - (0.5*(P*b^3 + R*T*b^2 - 1.0*a*b))/P + (((0.03703703703703704*(R*T - 1.0*P*b)^3)/P^3 - (0.5*(P*b^3 + R*T*b^2 - 1.0*a*b))/P + (0.1666666666666667*(R*T - 1.0*P*b)*(3.0*P*b^2 + 2.0*R*T*b - 1.0*a))/P^2)^2 - 1.0*((0.1111111111111111*(R*T - 1.0*P*b)^2)/P^2 + (0.3333333333333333*(3.0*P*b^2 + 2.0*R*T*b - 1.0*a))/P)^3)^(1/2) + (0.1666666666666667*(R*T - 1.0*P*b)*(3.0*P*b^2 + 2.0*R*T*b - 1.0*a))/P^2)^(1/3)*(0.5 + 0.8660254037844386i) + (0.3333333333333333*(R*T - 1.0*P*b))/P
2 Comments
  John D'Errico
      
      
 on 5 Sep 2016
				So what is the problem? I showed you how to do it. I won't do your homework for you.
See Also
Categories
				Find more on Mathematics and Optimization in Help Center and File Exchange
			
	Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!

