NAME PDL::Opt::QP - Quadratic programming solver for PDL SYNOPSIS use PDL::Opt::QP; ... DESCRIPTION ... This routine solves the quadratic programming optimization problem minimize f(x) = 0.5 x' D x - d' x x optionally constrained by: A' x = a B x >= b .... more docs to come .... }); pp_add_exported('', 'qp'); pp_addpm({At=>'Bot'},<<'EOD'); sub qp { my (\$Dmat, \$dvec, \$Amat, \$bvec, \$meq) = @_; my \$n = pdl \$Dmat->dim(1); # D is an [n x n] matrix my \$q = pdl \$Amat->dim(0); # A is an [n x q] matrix if( ! defined \$bvec ){ # if b is undef, create it \$bvec = zeros(\$q); } die("Dmat is not square!") if \$n != \$Dmat->dim(0); # Check D is [n x n] die("Dmat and dvec are incompatible!") if \$n != \$dvec->nelem; # Check d is [n] die("Amat and dvec are incompatible!") if \$n != \$Amat->dim(1); # Check A is [n x _] die("Amat and bvec are incompatible!") if \$q != \$bvec->nelem; # Check A is [_ x q] die("Value of meq is invalid!") if (\$meq > \$q) || (\$meq < 0 ); my \$iact = zeros(\$q); # Store which constraints are active my \$nact = pdl(0); # Store number of active constraints my \$r = \$n < \$q ? \$n : \$q; # Used to size work space my \$sol = zeros(\$n->sclr); # Store the solution [n] my \$lagr = zeros(\$q->sclr); # Store the Lagranges for the constraints my \$crval= pdl(0); # Value at min my \$work = zeros((2*\$n+\$r*(\$r+5)/2+2*\$q+1)->sclr); # Work space my \$iter = zeros(2); # Store info about interations my \$ierr = pdl(0); # Input: 1=Factorized; Output: error flag my \$res = qpgen2( \$Dmat->copy, \$dvec->copy, \$n, \$n, \$sol, \$lagr, \$crval, \$Amat->transpose->copy, \$bvec->copy, \$n, \$q, \$meq, \$iact, \$nact, \$iter, \$work, \$ierr ); croak "qp: constraints are inconsistent, no solution!" if \$ierr->sclr == 1; croak "qp: matrix D in quadratic function is not positive definite!" if \$ierr->sclr == 2; croak "qp: some problem with mininization" if \$ierr->sclr; return { x => \$sol, lagr => \$lagr, crval => \$crval, iact => \$iact, nact => \$nact, iter => \$iter, ierr => \$ierr, }; # TODO: process/return the results # # From R implementation: # # list(solution=res1\$sol, # value=res1\$crval, # unconstrained.solution=res1\$dvec, # iterations=res1\$iter, # Lagrangian = res1\$lagr, # iact=res1\$iact[1:res1\$nact]) } EOD pp_addpm({At=>'Bot'},<<'EOD'); SEE ALSO PDL, PDL::Opt::NonLinear BUGS Please report any bugs or suggestions at AUTHOR Mark Grimes, COPYRIGHT AND LICENSE This software is copyright (c) 2013 by Mark Grimes, . This is free software; you can redistribute it and/or modify it under the same terms as the Perl 5 programming language system itself.