%**copyright Elijah Polak, EECS UCB **
%%**version:    0 revision:     date: 27-MAR-98 **
%Polak-He unified method of feasible directions with Armijo Step size
%for second reliability maximization problem
%function H = ph(x,alpha,beta,gamma,n)
%returns matrix H containing iterates, values of cost, values of
%the optimality function theta, and  values of the constraint violation
%function psi. 
% H = ph(x,alpha,beta, gamma, n)
% H(:j) =  H(:j) = [x(j)', cost(j), theta(i) psi(j)]'
%requires initial vector x, alpha, beta, gamma, n=# of iterations
%requires vector vector function f(x), matrix function gradf(x)
%f1 is cost, f2-fm, constraint functions
function H = ph(x,alpha,beta,gamma,n)
x = x(:);
k=0;
nu_mu = (1/7)*[1 1 1 1 1 1 1]';
u = [1 1 1 1]';
u = u/norm(u);
psimax = psi(x);
z = [nu_mu;u;u;u;u;psimax;x];
VLB = vlb(x,psimax);
VUB = vub(x,psimax);
for i=1:n,
   zhat = constr('fun',z,[],VLB,VUB,'gradfun');
   Z(:,i) = zhat;
save Z Z;
   theta = -fun(zhat);
   H(:,i) = [x' psimax theta]';
save H H
   nu = zhat(1:4);
   mu = zhat(5:7);
   U = [zhat(8:11) zhat(12:15) zhat(16:19) zhat(20:23)];
   G = [gradgx(x,U(:,1)) gradgx(x,U(:,2)) gradgx(x,U(:,3)),gradgx(x,U(:,4))];
   h = -G*nu + mu(1)*[x(2);x(1)] + mu(2)*[1;-2] + mu(3)*[-1;0.5];
   DIR(:,i) = h;
save DIR
   %unified step size calculation (U = unified function)
   Unew = new_max(x , x + beta^k*h);
       if Unew >  beta^k*alpha*theta
	   while Unew > beta^k*alpha*theta
	     k=k+1;
             Unew = new_max(x , x + beta^k*h);
	   	   if k > 100, break
		   end %if
           end %while
       else 
	   while Unew <=  beta^k*alpha*theta
             k=k-1;
             Unew = new_max(x , x + beta^k*h);
 		   if k < -100, break
		   end %if
	   end %while
           k=k+1;
        end %if
   x = x + beta^k*h;
   psimax = psi(x);
   z = [zhat(1:23);psimax;x];
   k
   end %for
   H(:,n+1) = [x' psimax theta]';
save H H;