function [Jcol D B xD xB x] = col_obj(u1)

global d u NT
NT=41;

u = [u1 ; d(1:3)] ; %L, V, F, zF, qF

x0 = 0.5*ones(1,82)' ;
options = optimset('Display','none','TolFun',1e-10) ;
x = fsolve(@cola_lv,x0,options) ;
xB = 1 - x(1) ; 
xD = x(NT) ;

%Prices
pf = 1 ; pV = d(4) ; pB = 2 ; pD = 1 ;

F = d(1) ;
V = u(1) ;

KcB=10;  KcD=10;         % controller gains
MDs=0.5; MBs=0.5;        % Nominal holdups - these are rather small  
Ds=0.5; Bs=0.5;          % Nominal flows
MB=x(NT+1);  MD=x(2*NT); % Actual reboiler and condenser holdup
D=Ds+(MD-MDs)*KcD;       % Distillate flow
B=Bs+(MB-MBs)*KcB;       % Bottoms flow     

Jcol = pf*F + pV*V - pB*B - pD*D ; 

end

