New structure

+model/slackmodel.m

@@ -4,14 +4,16 @@
    % Derives the following slack formulation
    %
    % minimize  f(x,s)
-   % subj to   c(x) - s = 0
+   % subj to   ceq(x) = ceq
+   %           c(x) - s = 0
    %           bL <= x <= bU
    %           cL <= s <= cU
    %
    % from the inequality-based formulation:
    %
    % minimize  f(x)
-   % subj to   cL <= c(x) <= cU
+   % subj to   ceq(x) = ceq
+   %           cL <= c(x) <= cU
    %           bL <=   x  <= bU
    %
    % Thus, the problem appears as such:
@@ -21,8 +23,7 @@
    properties
       nlp    % original inequality-based object
       islack % indictor of slack variables
-      n_I    % number of inequality constraints
-      c_I    % inequality constraints
+      nslack % number of slack variables
    end
 
    properties (SetAccess = private, Hidden = false)
@@ -34,23 +35,20 @@
 
       function self = slackmodel(nlp)
 
-         constraints_I = find(nlp.cL ~= nlp.cU);
-         nI = size(constraints_I,1);
-
          % Upper and lower bounds for the variables and slacks.
-         bL = [ nlp.bL; nlp.cL(constraints_I) ];
-         bU = [ nlp.bU; nlp.cU(constraints_I) ];
+         bL = [ nlp.bL; nlp.cL(~nlp.iFix) ];
+         bU = [ nlp.bU; nlp.cU(~nlp.iFix) ];
 
          % The linear and nonlinear constraints are equalities, ie,
          % 0 <= c(x) - s <= 0.
          cL = nlp.cL;
+         cL(~nlp.iFix) = zeros(sum(~nlp.iFix), 1);
          cU = nlp.cU;
-         cL(constraints_I) = zeros(nI, 1);
-         cU(constraints_I) = zeros(nI, 1);
+         cU(~nlp.iFix) = zeros(sum(~nlp.iFix), 1);
 
          % Initial point. Set slacks to be feasible.
          c = nlp.fcon(nlp.x0);
-         x0 = [ nlp.x0; c(constraints_I) ];
+         x0 = [ nlp.x0; c(~nlp.iFix) ];
 
          % Instantiate from the base class.
          self = [email protected](nlp.name, x0, cL, cU, bL, bU);
@@ -59,18 +57,19 @@
          self.linear = nlp.linear;
 
          % Create an indetifier for slack variables.
-         self.islack = [ false(nlp.n,1); true(nI,1) ];   
-
+         nS = sum(~nlp.iFix);
+         self.nslack = nS;
+         self.islack = [ false(nlp.n,1); true(nS,1) ];   
+
          % Jacobian sparsity pattern of the slack model.
-         Jx = nlp.gcon(nlp.x0);
-         J = [spones(Jx) sparse(nlp.m, nI)];
-         J(constraints_I , nlp.n +1: nlp.n + nI) = speye(nI);
-         self.Jpattern = J;
-
+         J = nlp.gcon(nlp.x0);
+         Js = sparse(nlp.m,nS);
+         Js(~nlp.iFix,:) = speye(nS);
+         self.Jpattern = [spones(J)  Js];
+         
          % Hessian sparsity pattern.
          y = ones(size(c));
          HL = nlp.hlag(nlp.x0, y);
-         nS = self.n - nlp.n;
          self.Hpattern = [ spones(HL)         sparse(nlp.n, nS)
                            sparse(nS, nlp.n)  sparse(nS   , nS) ];
 
@@ -81,10 +80,6 @@
          % Store the original NLP model.
          self.nlp = nlp;
 
-         % Store inequality stuff
-         self.n_I = nI;
-         self.c_I = constraints_I;
-
       end
 
       function f = fobj_local(self, xs)
@@ -95,48 +90,74 @@
       function g = gobj_local(self, xs)
          x = xs(~self.islack,:);
          gx = self.nlp.gobj(x);
-         g = [gx; sparse(self.n_I, 1)];
+         g = [gx; zeros(self.nslack, 1)];
       end
 
       function H = hobj_local(self, xs)
          x = xs(~self.islack,:);
          Hx = self.nlp.hobj(x);
-         nmZ = sparse(self.nlp.n, self.n_I);
-         mmZ = sparse(self.n_I, self.n_I);
+         nmZ = sparse(self.nlp.n, self.nslack);
+         mmZ = sparse(self.nslack, self.nslack);
          H = [ Hx     nmZ
                nmZ'   mmZ ];
       end
 
       function c = fcon_local(self, xs)
          x = xs(~self.islack,:);
          s = xs( self.islack,:);
-         c = self.nlp.fcon(x) ;
-         c(self.c_I) = c(self.c_I) - s;
+         c = self.nlp.fcon(x);
+         c(~self.nlp.iFix) = c(~self.nlp.iFix) - s;
       end
 
       function J = gcon_local(self, xs)
          x = xs(~self.islack,:);
          Jx = self.nlp.gcon(x);
-         J = [Jx sparse(self.nlp.m,self.n_I)];
-         J(self.c_I , self.nlp.n +1: self.nlp.n + self.n_I) = -speye(self.n_I);
+         Js = sparse(self.m, self.nslack);
+         Js(~self.nlp.iFix,:) = -speye(self.nslack);
+         J = [Jx Js];
+      end
+
+      function [Jprod, Jtprod] = gconprod_local(self, xs)
+         % J = self.gcon(x);
+         [Jxprod, Jxtprod] = self.nlp.gconprod(xs(~self.islack,:));
+         Js = sparse(self.m, self.nslack);
+         Js(~self.nlp.iFix,:) = -speye(self.nslack);
+         n = sum(~self.islack);
+
+         Jprod = @(v) Jxprod(v(1:n)) + Js*v(n+1:end);
+         Jtprod = @(v) [Jxtprod(v); Js'*v];
       end
 
       function HL = hlag_local(self, xs, y)
          x = xs(~self.islack,:);
          H = self.nlp.hlag(x, y);
-         nmZ = sparse(self.nlp.n, self.nI);
-         mmZ = sparse(self.nI,self.nI);
+         nmZ = zeros(self.nlp.n, self.nslack);
+         mmZ = zeros(self.nslack);
          HL = [ H     nmZ
                 nmZ'  mmZ ];
       end
+
+      function Hv = hconprod_local(self, xs, y, vv)
+         x = xs(~self.islack);
+         v = vv(~self.islack);
+         Hv = zeros(self.n, 1);
+         Hv(~self.islack) = self.nlp.hconprod(x, y, v);
+      end
 
       function Hv = hlagprod_local(self, xs, y, vv)
          x = xs(~self.islack);
          v = vv(~self.islack);
-         Hv = sparse(self.n, 1);
+         Hv = zeros(self.n, 1);
          Hv(~self.islack) = self.nlp.hlagprod(x, y, v);
       end
 
+      function z = ghivprod_local(self, xs, gxs, vxs)
+         x = xs (~self.islack);
+         g = gxs(~self.islack);
+         v = vxs(~self.islack);
+         z = self.nlp.ghivprod(x, g, v);
+      end
+
    end % methods
 
 end % classdef

Formulations/K1x.m

@@ -12,7 +12,7 @@
     end
 
     methods
-        function o = K1x(slack, options)
+        function o = K1x(options)
             o.diagHess = true;
         end
 

Formulations/K1y.m

@@ -1,4 +1,4 @@
-classdef K1y < pdcoO
+classdef K1y < handle
 
     properties
         M
@@ -11,8 +11,7 @@
     end
 
     methods
-        function o = K1y(slack, options)
-            o = o@pdcoO(slack, options);
+        function o = K1y(options)
             o.diagHess = true;
         end
 

Formulations/K1y_ls.m

@@ -1,4 +1,4 @@
-classdef K1y_ls < pdcoO
+classdef K1y_ls < handle
 
     properties
         M
@@ -12,8 +12,7 @@
     end
 
     methods
-        function o = K1y_ls(slack, options)
-            o = o@pdcoO(slack, options);
+        function o = K1y_ls(options)
             o.diagHess = true;
         end
 

Formulations/K2.m

@@ -1,4 +1,4 @@
-classdef K2 < pdcoO
+classdef K2 < handle
 
     properties
         M
@@ -12,8 +12,7 @@
     end
 
     methods
-        function o = K2(slack, options)
-            o = o@pdcoO(slack, options);
+        function o = K2(options)
             o.diagHess = false;
         end
 

Formulations/K25.m

@@ -1,4 +1,4 @@
-classdef K25 < pdcoO
+classdef K25 < handle
 
     properties
         M
@@ -11,8 +11,7 @@
     end
 
     methods
-        function o = K25(slack, options)
-            o = o@pdcoO(slack, options);
+        function o = K25(options)
             o.diagHess = false;
         end
 

Formulations/K2_CPCG.m

@@ -1,4 +1,4 @@
-classdef K2_CPCG < pdcoO
+classdef K2_CPCG < handle
 
     properties
         C
@@ -13,8 +13,7 @@
     end
 
     methods
-        function o = K2_CPCG(slack,options)
-            o = o@pdcoO(slack,options);
+        function o = K2_CPCG(options)
             o.diagHess = false;
         end
 

Formulations/K2_ls.m

@@ -1,4 +1,4 @@
-classdef K2_ls < pdcoO
+classdef K2_ls < handle
 
     properties
         M
@@ -12,8 +12,7 @@
     end
 
     methods
-        function o = K2_ls(slack, options)
-            o = o@pdcoO(slack, options);
+        function o = K2_ls(options)
             o.diagHess = true;
         end
 

Formulations/K35.m

@@ -1,4 +1,4 @@
-classdef K35 < pdcoO
+classdef K35 < handle
 
     properties
         M
@@ -12,8 +12,7 @@
     end
 
     methods
-        function o = K35(slack, options)
-            o = o@pdcoO(slack, options);
+        function o = K35(options)
             o.diagHess = false;
         end
 

Formulations/K35_CPCG.m

@@ -1,4 +1,4 @@
-classdef K35_CPCG < pdcoO
+classdef K35_CPCG < handle
 
     properties
         C
@@ -13,8 +13,7 @@
     end
 
     methods
-        function o = K35_CPCG(slack, options)
-            o = o@pdcoO(slack, options);
+        function o = K35_CPCG(options)
             o.diagHess = false;
         end
 

Solvers/CPCG.m

@@ -1,4 +1,4 @@
-classdef CPCG < pdcoO
+classdef CPCG < handle
 
     properties
         itnlim
@@ -11,23 +11,22 @@
     end
 
     methods
-        function o = CPCG(slack, options)
-            o = o@pdcoO(slack, options);
+        function o = CPCG(options)
 
-            if isfield(options, 'LSMRatol1')
-                o.atol1 = options.LSMRatol1;
+            if isfield(options, 'atol1')
+                o.atol1 = options.atol1;
             else
                 o.atol1 = 1e-10;
             end
 
-            if isfield(options, 'LSMRatol2')
-                o.atol2 = options.LSMRatol2;
+            if isfield(options, 'atol2')
+                o.atol2 = options.atol2;
             else
                 o.atol2 = 1e-15;
             end
 
-            if isfield(options, 'LSMRMaxIter')
-                o.itnlim = options.LSMRMaxIter * min(o.m, o.n);
+            if isfield(options, 'itnlim')
+                o.itnlim = options.itnlim * min(o.m, o.n);
             else
                 o.itnlim = 10;
             end

Solvers/Cholesky.m

@@ -6,7 +6,7 @@
     end
 
     methods
-        function o = Cholesky(slack, options)
+        function o = Cholesky(options)
             o.manage_op = false;
             o.need_precon = false;
             o.solver = '  Chol';  o.head3 = '     Chol';