closestProjectionOntoConfinementPartition

 find the closest point in some subdomain

 THIS IS NO USER FUNCTION

• projectVectorDec approximize hE-vectors to the domain covered by an improper decomposition

0001 function heMatrix = closestProjectionOntoConfinementPartition(...
0002     heMatrix, x0, A, projectionMatrix)
0003 % find the closest point in some subdomain
0004 %
0005 % THIS IS NO USER FUNCTION
0006 
0007 % The elk-library: convex geometry applied to crystallization modeling.
0008 %   Copyright (C) 2013 Alexander Reinhold
0009 %
0010 % This program is free software: you can redistribute it and/or modify it
0011 %   under the terms of the GNU General Public License as published by the
0012 %   Free Software Foundation, either version 3 of the License, or (at your
0013 %   option) any later version.
0014 %
0015 % This program is distributed in the hope that it will be useful, but
0016 %   WITHOUT ANY WARRANTY; without even the implied warranty of
0017 %   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
0018 %   General Public License for more details.
0019 %
0020 % You should have received a copy of the GNU General Public License along
0021 %   with this program.  If not, see <http://www.gnu.org/licenses/>.
0022 
0023 
0024 %% notes
0025 % Given two points x and xOpt, there distance (without sqrt) is:
0026 %     (xOpt - x)'*(xOpt - x)
0027 %   giving a necessarily quadratic problem. The quadprog minimizes:
0028 %     1/2*x'*H*x + f'*x
0029 %   where H is a positive definite matrix. Reformulating the distance,
0030 %   gives:
0031 %     xOpt'*xOpt - 2*xOpt'*x + x'*x
0032 %   for which we can ignore the last term and divide by two to obtain:
0033 %     H = 1
0034 %     f = -1*x
0035 
0036 nE = size(heMatrix, 1);
0037 nVec = size(heMatrix, 2);
0038 
0039 problem.solver = 'quadprog';
0040 problem.options = optimset('quadprog');
0041 problem.options.Algorithm = 'active-set';
0042 problem.options.Display = 'off';
0043 problem.H = eye(nE);
0044 problem.Aineq = A;
0045 problem.bineq = zeros(1, size(A, 1));
0046 problem.Aeq = eye(nE) - projectionMatrix;
0047 problem.beq = zeros(nE, 1);
0048 problem.x0 = x0;
0049 
0050 problem.f = -1*heMatrix(:, 1);
0051 objective = @(x)(1/2*x'*problem.H*x + problem.f'*x);
0052 objective(x0)
0053 objective(heMatrix(:,1))
0054 
0055 for iVec = 1:nVec
0056     iVec/nVec*100
0057     problem.f = -1*heMatrix(:, iVec);
0058     heMatrix(:, iVec) = quadprog(problem);
0059 %     [heMatrix(:, iVec),fval,exitflag,output,lambda] = quadprog(problem)
0060 %     lambda.lower,lambda.upper,lambda.eqlin,lambda.ineqlin
0061 end
0062