Economic Operation of Fixed-bed Filter
Authors' Data: 
Housam BINOUS & Thomas F. Edgar
Department of Chemical Engineering
National Institute of Applied Sciences and Technology
Tunis, TUNISIA
University of Texas at Austin
Texas, USA
Emails: binoushousam@yahoo.com 
edgar@che.utexas.edu

Acknowledgement :
This problem was presented in Optimization of chemical processes, Edgar et al., second edition, McGrawHill, 2001. Results found by the authors are slightly different and are given in page 469 of their book.

GAMS: 
GAMS program gives the same solution. For GAMS code, please send request to binoushousam@yahoo.com.

Objective function

f=116 ((10^6 q/(1440 (Vf-Vb)/(Vf/Q+tb))))^0.86+4.73 10^3 (Vf/(Vf-Vb)-1) q

q=10

tb=10

Vb=200

Contour plot of objective function

cont=ContourPlot[f,{Vf,800,2000},{Q,1,20},ColorFunction ->Hue,PlotPoints->100,Contours->Table[i 5000,{i,1,30}],DisplayFunction->Identity]

Equality constraint

Vf=700 3 Log[10,10/0.081 1/3 1/Q]

Plot of Q versus Vf

plt=Plot[Vf,{Q,1,20},PlotStyle->Thickness[0.01],DisplayFunction->Identity]

plt2=Show[plt/.x_Line:>Map[Reverse,x,{2}],PlotRange->{{800,2000},{1,20}}]

Combining objective function and equality constraint in same figure

Show[cont,plt2,PlotRange->{{800,2000},{1,20}},DisplayFunction->$DisplayFunction]

Vf=.

Lagrangian function

L=f+ω (Vf-700 3 Log[10,10/0.081 1/3 1/Q])

First-order necessary condition for an extremum

Eq1=D[L,{Vf,1}]==0

Eq2=D[L,{Q,1}]==0

Eq3=D[L,{ω,1}]==0

Extremum

sol=Chop[FindRoot[{Eq1,Eq2,Eq3},{Vf,900},{Q,10},{ω,10}]]

Value of the objective function at the extremum

f/.sol

Checking second-order necessary and sufficient condition for optimality

Hessian and Gradient definitions

Off[General::spell1]

Grad[fon_,var_List]:=D[fon,#] & /@ var;
Hess[fon_,var_List]:=Grad[D[fon,#] &/@ var,var];

Hessian is positive definite at extremum

H=Hess[L,{Vf,Q}]/.sol

Hess[L,{Vf,Q}]/.sol//MatrixForm

Eigenvalues[H]

Gradient of the active constraints

G=Vf-700 3 Log[10,10/0.081 1/3 1/Q]

Grd=Grad[G,{Vf,Q}]/.sol

Solve[{x,y}.Grd==0,{x,y}]

(y^T)∇\(TraditionalForm\`\_x\^2\)L y  is positive for all vectors, y, such that Grd.y =0. This indicates that the extremum is a local minimum

{-64.4473 y,y}.H.{-64.4473 y,y}