#!/bin/sh
# This is a shell archive (produced by GNU sharutils 4.2).
# To extract the files from this archive, save it to some FILE, remove
# everything before the `!/bin/sh' line above, then type `sh FILE'.
#
# Made on 1998-07-23 15:11 CST by .
# Source directory was `/usr/local/imports/algona1/staff/bturlach/lang/splus/Develop/QuadProg'.
#
# Existing files will *not* be overwritten unless `-c' is specified.
#
# This shar contains:
# length mode name
# ------ ---------- ------------------------------------------
# 1430 -r--r--r-- QuadProg/achck.f
# 1468 -r--r--r-- QuadProg/aind.f
# 18522 -r--r--r-- QuadProg/solve.QP.compact.f
# 17890 -r--r--r-- QuadProg/solve.QP.f
# 3581 -r--r--r-- QuadProg/util.f
# 4034 -r--r--r-- QuadProg/solve.QP.compact.d
# 3255 -r--r--r-- QuadProg/solve.QP.d
# 7036 -rw-r--r-- QuadProg/QuadProg.q
# 6834 -rw-r--r-- QuadProg/Makefile
# 25266 -rw-r--r-- QuadProg/COPYING.LIB
# 231 -rw-r--r-- QuadProg/README
# 3833 -rw-r--r-- QuadProg/README.1st
# 2553 -rw-r--r-- QuadProg/README.INSTALL
# 104 -rw-r--r-- QuadProg/README.pd
# 3234 -rw-r--r-- QuadProg/ChangeLog
#
save_IFS="${IFS}"
IFS="${IFS}:"
gettext_dir=FAILED
locale_dir=FAILED
first_param="$1"
for dir in $PATH
do
if test "$gettext_dir" = FAILED && test -f $dir/gettext \
&& ($dir/gettext --version >/dev/null 2>&1)
then
set `$dir/gettext --version 2>&1`
if test "$3" = GNU
then
gettext_dir=$dir
fi
fi
if test "$locale_dir" = FAILED && test -f $dir/shar \
&& ($dir/shar --print-text-domain-dir >/dev/null 2>&1)
then
locale_dir=`$dir/shar --print-text-domain-dir`
fi
done
IFS="$save_IFS"
if test "$locale_dir" = FAILED || test "$gettext_dir" = FAILED
then
echo=echo
else
TEXTDOMAINDIR=$locale_dir
export TEXTDOMAINDIR
TEXTDOMAIN=sharutils
export TEXTDOMAIN
echo="$gettext_dir/gettext -s"
fi
touch -am 1231235999 $$.touch >/dev/null 2>&1
if test ! -f 1231235999 && test -f $$.touch; then
shar_touch=touch
else
shar_touch=:
echo
$echo 'WARNING: not restoring timestamps. Consider getting and'
$echo "installing GNU \`touch', distributed in GNU File Utilities..."
echo
fi
rm -f 1231235999 $$.touch
#
if mkdir _sh14571; then
$echo 'x -' 'creating lock directory'
else
$echo 'failed to create lock directory'
exit 1
fi
# ============= QuadProg/achck.f ==============
if test ! -d 'QuadProg'; then
$echo 'x -' 'creating directory' 'QuadProg'
mkdir 'QuadProg'
fi
if test -f 'QuadProg/achck.f' && test "$first_param" != -c; then
$echo 'x -' SKIPPING 'QuadProg/achck.f' '(file already exists)'
else
$echo 'x -' extracting 'QuadProg/achck.f' '(text)'
sed 's/^X//' << 'SHAR_EOF' > 'QuadProg/achck.f' &&
c Copyright (C) 1998
c Berwin A. Turlach
c $Id: achck.f,v 1.2 1998/07/23 05:05:58 bturlach Exp $
c
c This library is free software; you can redistribute it and/or
c modify it under the terms of the GNU Library General Public
c License as published by the Free Software Foundation; either
c version 2 of the License, or (at your option) any later version.
c
c This library is distributed in the hope that it will be useful,
c but WITHOUT ANY WARRANTY; without even the implied warranty of
c MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
c Library General Public License for more details.
c
c You should have received a copy of the GNU Library General Public
c License along with this library; if not, write to the Free Software
c Foundation, Inc., 59 Temple Place, Suite 330, Boston,
c MA 02111-1307 USA
c
c this routine checks whether all constraints are fulfilled.
c
X subroutine achck(sol,n,amat,aind,bvec,m,q,meq,prec,ok)
X implicit none
X integer m, aind(m+1,*), q, meq, n, i, j
X double precision sol(*), amat(m,*), bvec(*), sum, prec
X logical ok
X ok = .FALSE.
X do i=1,q
X sum = -bvec(i)
X do j=1,aind(1,i)
X sum = sum + amat(j,i)*sol(aind(j+1,i))
X enddo
X if( i .LE. meq ) sum = -abs(sum)
X if( sum .LT. -prec ) return
X enddo
X ok = .TRUE.
X return
X end
SHAR_EOF
$shar_touch -am 0723151198 'QuadProg/achck.f' &&
chmod 0444 'QuadProg/achck.f' ||
$echo 'restore of' 'QuadProg/achck.f' 'failed'
if ( md5sum --help 2>&1 | grep 'sage: md5sum \[' ) >/dev/null 2>&1 \
&& ( md5sum --version 2>&1 | grep -v 'textutils 1.12' ) >/dev/null; then
md5sum -c << SHAR_EOF >/dev/null 2>&1 \
|| $echo 'QuadProg/achck.f:' 'MD5 check failed'
ce936bd0ac4ff21abd68464c640796bd QuadProg/achck.f
SHAR_EOF
else
shar_count="`LC_ALL= LC_CTYPE= LANG= wc -c < 'QuadProg/achck.f'`"
test 1430 -eq "$shar_count" ||
$echo 'QuadProg/achck.f:' 'original size' '1430,' 'current size' "$shar_count!"
fi
fi
# ============= QuadProg/aind.f ==============
if test -f 'QuadProg/aind.f' && test "$first_param" != -c; then
$echo 'x -' SKIPPING 'QuadProg/aind.f' '(file already exists)'
else
$echo 'x -' extracting 'QuadProg/aind.f' '(text)'
sed 's/^X//' << 'SHAR_EOF' > 'QuadProg/aind.f' &&
c Copyright (C) 1997, 1998
c Berwin A. Turlach
c $Id: aind.f,v 1.3 1998/07/23 05:06:32 bturlach Exp $
c
c This library is free software; you can redistribute it and/or
c modify it under the terms of the GNU Library General Public
c License as published by the Free Software Foundation; either
c version 2 of the License, or (at your option) any later version.
c
c This library is distributed in the hope that it will be useful,
c but WITHOUT ANY WARRANTY; without even the implied warranty of
c MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
c Library General Public License for more details.
c
c You should have received a copy of the GNU Library General Public
c License along with this library; if not, write to the Free Software
c Foundation, Inc., 59 Temple Place, Suite 330, Boston,
c MA 02111-1307 USA
c
c this routine checks whether Aind has valid entries, i.e.,
c 1) 1<= Aind(1,i) <= n for i=1,...,q (number of constraints)
c 2) 1<= Aind(j,i) <= n for j=2,...,Aind(1,i)+1, i=1,...,q
c
c Aind is a m times q matrix
c
X subroutine aind(ind,m,q,n,ok)
X implicit none
X integer m, ind(m,*), q, n, i, j
X logical ok
X ok = .FALSE.
X do i=1,q
X if( ind(1,i) .LT. 1 .OR. ind(1,i) .GT. n ) return
X do j=2,ind(1,i)+1
X if( ind(j,i) .LT. 1 .OR. ind(j,i) .GT. n ) return
X enddo
X enddo
X ok = .TRUE.
X return
X end
SHAR_EOF
$shar_touch -am 0723151198 'QuadProg/aind.f' &&
chmod 0444 'QuadProg/aind.f' ||
$echo 'restore of' 'QuadProg/aind.f' 'failed'
if ( md5sum --help 2>&1 | grep 'sage: md5sum \[' ) >/dev/null 2>&1 \
&& ( md5sum --version 2>&1 | grep -v 'textutils 1.12' ) >/dev/null; then
md5sum -c << SHAR_EOF >/dev/null 2>&1 \
|| $echo 'QuadProg/aind.f:' 'MD5 check failed'
a32e1fb3b21b90a0906d10ec0dcda561 QuadProg/aind.f
SHAR_EOF
else
shar_count="`LC_ALL= LC_CTYPE= LANG= wc -c < 'QuadProg/aind.f'`"
test 1468 -eq "$shar_count" ||
$echo 'QuadProg/aind.f:' 'original size' '1468,' 'current size' "$shar_count!"
fi
fi
# ============= QuadProg/solve.QP.compact.f ==============
if test -f 'QuadProg/solve.QP.compact.f' && test "$first_param" != -c; then
$echo 'x -' SKIPPING 'QuadProg/solve.QP.compact.f' '(file already exists)'
else
$echo 'x -' extracting 'QuadProg/solve.QP.compact.f' '(text)'
sed 's/^X//' << 'SHAR_EOF' > 'QuadProg/solve.QP.compact.f' &&
c Copyright (C) 1995, 1996, 1997, 1998
c Berwin A. Turlach
c $Id: solve.QP.compact.f,v 1.13 1998/07/23 05:08:23 bturlach Exp $
c
c This library is free software; you can redistribute it and/or
c modify it under the terms of the GNU Library General Public
c License as published by the Free Software Foundation; either
c version 2 of the License, or (at your option) any later version.
c
c This library is distributed in the hope that it will be useful,
c but WITHOUT ANY WARRANTY; without even the implied warranty of
c MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
c Library General Public License for more details.
c
c You should have received a copy of the GNU Library General Public
c License along with this library; if not, write to the Free Software
c Foundation, Inc., 59 Temple Place, Suite 330, Boston,
c MA 02111-1307 USA
c
c This routine uses the Goldfarb/Idnani algorithm to solve the
c following minimization problem:
c
c minimize -d^T x + 1/2 * x^T D x
c where A1^T x = b1
c A2^T x >= b2
c
c the matrix D is assumed to be positive definite. Especially,
c w.l.o.g. D is assumed to be symmetric.
c
c Input parameter:
c dmat nxn matrix, the matrix D from above (dp)
c *** WILL BE DESTROYED ON EXIT ***
c The user has two possibilities:
c a) Give D (ierr=0), in this case we use routines from LINPACK
c to decompose D.
c b) To get the algorithm started we need R^-1, where D=R^TR.
c So if it is cheaper to calculate R^-1 in another way (D may
c be a band matrix) then with the general routine, the user
c may pass R^{-1}. Indicated by ierr not equal to zero.
c dvec nx1 vector, the vector d from above (dp)
c *** WILL BE DESTROYED ON EXIT ***
c contains on exit the solution to the initial, i.e.,
c unconstrained problem
c fddmat scalar, the leading dimension of the matrix dmat
c n the dimension of dmat and dvec (int)
c amat lxq matrix (dp)
c *** ENTRIES CORRESPONDING TO EQUALITY CONSTRAINTS MAY HAVE
c CHANGED SIGNES ON EXIT ***
c iamat (l+1)xq matrix (int)
c these two matrices store the matrix A in compact form. the format
c is: [ A=(A1 A2) ]
c iamat(1,i) is the number of non-zero elements in column i of A
c iamat(k,i) for k>=2, is equal to j if the (k-1)-th non-zero
c element in column i of A is A(i,j)
c amat(k,i) for k>=1, is equal to the k-th non-zero element
c in column i of A.
c
c bvec qx1 vector, the vector of constants b in the constraints (dp)
c [ b = (b1^T b2^T)^T ]
c *** ENTRIES CORRESPONDING TO EQUALITY CONSTRAINTS MAY HAVE
c CHANGED SIGNES ON EXIT ***
c fdamat the first dimension of amat as declared in the calling program.
c fdamat >= n (and iamat must have fdamat+1 as first dimension)
c q integer, the number of constraints.
c meq integer, the number of equality constraints, 0 <= meq <= q.
c ierr integer, code for the status of the matrix D:
c ierr = 0, we have to decompose D
c ierr != 0, D is already decomposed into D=R^TR and we were
c given R^{-1}.
c
c Output parameter:
c sol nx1 the final solution (x in the notation above)
c crval scalar, the value of the criterion at the minimum
c iact qx1 vector, the constraints which are active in the final
c fit (int)
c nact scalar, the number of constraints active in the final fit (int)
c iter 2x1 vector, first component gives the number of "main"
c iterations, the second one says how many constraints were
c deleted after they became active
c ierr integer, error code on exit, if
c ierr = 0, no problems
c ierr = 1, the minimization problem has no solution
c ierr = 2, problems with decomposing D, in this case sol
c contains garbage!!
c
c Working space:
c work vector with length at least 2*n+r*(r+5)/2 + 2*q +1
c where r=min(n,q)
c
X subroutine qpgen1(dmat, dvec, fddmat, n, sol, crval, amat, iamat,
X * bvec, fdamat, q, meq, iact, nact, iter, work, ierr)
X implicit none
X integer n, i, j, l, l1, fdamat, fddmat,
X * info, q, iamat(fdamat+1,*), iact(*), iter(*), it1,
X * ierr, nact, iwzv, iwrv, iwrm, iwsv, iwuv, nvl,
X * r, iwnbv, meq
X double precision dmat(fddmat,*), dvec(*),sol(*), bvec(*),
X * work(*), temp, sum, t1, tt, gc, gs, crval,
X * nu, amat(fdamat,*)
X logical t1inf, t2min
X r = min(n,q)
X l = 2*n + (r*(r+5))/2 + 2*q + 1
c
c store the initial dvec to calculate below the unconstrained minima of
c the critical value.
c
X do 10 i=1,n
X work(i) = dvec(i)
X 10 continue
X do 11 i=n+1,l
X work(i) = 0.d0
X 11 continue
X do 12 i=1,q
X iact(i)=0
X 12 continue
c
c get the initial solution
c
X if( ierr .EQ. 0 )then
X call dpofa(dmat,fddmat,n,info)
X if( info .NE. 0 )then
X ierr = 2
X goto 999
X endif
X call dposl(dmat,fddmat,n,dvec)
X call dpori(dmat,fddmat,n)
X else
c
c Matrix D is already factorized, so we have to multiply d first with
c R^-T and then with R^-1. R^-1 is stored in the upper half of the
c array dmat.
c
X do 20 j=1,n
X sol(j) = 0.d0
X do 21 i=1,j
X sol(j) = sol(j) + dmat(i,j)*dvec(i)
X 21 continue
X 20 continue
X do 22 j=1,n
X dvec(j) = 0.d0
X do 23 i=j,n
X dvec(j) = dvec(j) + dmat(j,i)*sol(i)
X 23 continue
X 22 continue
X endif
c
c set lower triangular of dmat to zero, store dvec in sol and
c calculate value of the criterion at unconstrained minima
c
X crval = 0.d0
X do 30 j=1,n
X sol(j) = dvec(j)
X crval = crval + work(j)*sol(j)
X work(j) = 0.d0
X do 32 i=j+1,n
X dmat(i,j) = 0.d0
X 32 continue
X 30 continue
X crval = -crval/2.d0
X ierr = 0
c
c calculate some constants, i.e., from which index on the different
c quantities are stored in the work matrix
c
X iwzv = n
X iwrv = iwzv + n
X iwuv = iwrv + r
X iwrm = iwuv + r+1
X iwsv = iwrm + (r*(r+1))/2
X iwnbv = iwsv + q
c
c calculate the norm of each column of the A matrix
c
X do 51 i=1,q
X sum = 0.d0
X do 52 j=1,iamat(1,i)
X sum = sum + amat(j,i)*amat(j,i)
X 52 continue
X work(iwnbv+i) = sqrt(sum)
X 51 continue
X nact = 0
X iter(1) = 0
X iter(2) = 0
X 50 continue
c
c start a new iteration
c
X iter(1) = iter(1)+1
c
c calculate all constraints and check which are still violated
c for the equality constraints we have to check whether the normal
c vector has to be negated (as well as bvec in that case)
c
X l = iwsv
X do 60 i=1,q
X l = l+1
X sum = -bvec(i)
X do 61 j = 1,iamat(1,i)
X sum = sum + amat(j,i)*sol(iamat(j+1,i))
X 61 continue
X if (i .GT. meq) then
X work(l) = sum
X else
X work(l) = -abs(sum)
X if (sum .GT. 0.d0) then
X do 62 j=1,iamat(1,i)
X amat(j,i) = -amat(j,i)
X 62 continue
X bvec(i) = -bvec(i)
X endif
X endif
X 60 continue
c
c as safeguard against rounding errors set already active constraints
c explicitly to zero
c
X do 70 i=1,nact
X work(iwsv+iact(i)) = 0.d0
X 70 continue
c
c we weight each violation by the number of non-zero elements in the
c corresponding row of A. then we choose the violated constraint which
c has maximal absolute value, i.e., the minimum.
c by obvious commenting and uncommenting we can choose the strategy to
c take always the first constraint which is violated. ;-)
c
X nvl = 0
X temp = 0.d0
X do 71 i=1,q
X if (work(iwsv+i) .LT. temp*work(iwnbv+i)) then
X nvl = i
X temp = work(iwsv+i)/work(iwnbv+i)
X endif
c if (work(iwsv+i) .LT. 0.d0) then
c nvl = i
c goto 72
c endif
X 71 continue
X 72 if (nvl .EQ. 0) goto 999
c
c calculate d=J^Tn^+ where n^+ is the normal vector of the violated
c constraint. J is stored in dmat in this implementation!!
c if we drop a constraint, we have to jump back here.
c
X 55 continue
X do 80 i=1,n
X sum = 0.d0
X do 81 j=1,iamat(1,nvl)
X sum = sum + dmat(iamat(j+1,nvl),i)*amat(j,nvl)
X 81 continue
X work(i) = sum
X 80 continue
c
c Now calculate z = J_2 d_2
c
X l1 = iwzv
X do 90 i=1,n
X work(l1+i) =0.d0
X 90 continue
X do 92 j=nact+1,n
X do 93 i=1,n
X work(l1+i) = work(l1+i) + dmat(i,j)*work(j)
X 93 continue
X 92 continue
c
c and r = R^{-1} d_1, check also if r has positive elements (among the
c entries corresponding to inequalities constraints).
c
X t1inf = .TRUE.
X do 95 i=nact,1,-1
X sum = work(i)
X l = iwrm+(i*(i+3))/2
X l1 = l-i
X do 96 j=i+1,nact
X sum = sum - work(l)*work(iwrv+j)
X l = l+j
X 96 continue
X sum = sum / work(l1)
X work(iwrv+i) = sum
X if (iact(i) .LE. meq) goto 95
X if (sum .LE. 0.d0) goto 95
X 7 t1inf = .FALSE.
X it1 = i
X 95 continue
c
c if r has positive elements, find the partial step length t1, which is
c the maximum step in dual space without violating dual feasibility.
c it1 stores in which component t1, the min of u/r, occurs.
c
X if ( .NOT. t1inf) then
X t1 = work(iwuv+it1)/work(iwrv+it1)
X do 100 i=1,nact
X if (iact(i) .LE. meq) goto 100
X if (work(iwrv+i) .LE. 0.d0) goto 100
X temp = work(iwuv+i)/work(iwrv+i)
X if (temp .LT. t1) then
X t1 = temp
X it1 = i
X endif
X 100 continue
X endif
c
c test if the z vector is equal to zero
c
X sum = 0.d0
X do 110 i=iwzv+1,iwzv+n
X sum = sum + work(i)*work(i)
X 110 continue
X temp = 1000.d0
X sum = sum+temp
X if (temp .EQ. sum) then
c
c No step in pmrimal space such that the new constraint becomes
c feasible. Take step in dual space and drop a constant.
c
X if (t1inf) then
c
c No step in dual space possible either, problem is not solvable
c
X ierr = 1
X goto 999
X else
c
c we take a partial step in dual space and drop constraint it1,
c that is, we drop the it1-th active constraint.
c then we continue at step 2(a) (marked by label 55)
c
X do 111 i=1,nact
X work(iwuv+i) = work(iwuv+i) - t1*work(iwrv+i)
X 111 continue
X work(iwuv+nact+1) = work(iwuv+nact+1) + t1
X goto 700
X endif
X else
c
c compute full step length t2, minimum step in primal space such that
c the constraint becomes feasible.
c keep sum (which is z^Tn^+) to update crval below!
c
X sum = 0.d0
X do 120 i = 1,iamat(1,nvl)
X sum = sum + work(iwzv+iamat(i+1,nvl))*amat(i,nvl)
X 120 continue
X tt = -work(iwsv+nvl)/sum
X t2min = .TRUE.
X if (.NOT. t1inf) then
X if (t1 .LT. tt) then
X tt = t1
X t2min = .FALSE.
X endif
X endif
c
c take step in primal and dual space
c
X do 130 i=1,n
X sol(i) = sol(i) + tt*work(iwzv+i)
X 130 continue
X crval = crval + tt*sum*(tt/2.d0 + work(iwuv+nact+1))
X do 131 i=1,nact
X work(iwuv+i) = work(iwuv+i) - tt*work(iwrv+i)
X 131 continue
X work(iwuv+nact+1) = work(iwuv+nact+1) + tt
c
c if it was a full step, then we check wheter further constraints are
c violated otherwise we can drop the current constraint and iterate once
c more
X if(t2min) then
c
c we took a full step. Thus add constraint nvl to the list of active
c constraints and update J and R
c
X nact = nact + 1
X iact(nact) = nvl
c
c to update R we have to put the first nact-1 components of the d vector
c into column (nact) of R
c
X l = iwrm + ((nact-1)*nact)/2 + 1
X do 150 i=1,nact-1
X work(l) = work(i)
X l = l+1
X 150 continue
c
c if now nact=n, then we just have to add the last element to the new
c row of R.
c Otherwise we use Givens transformations to turn the vector d(nact:n)
c into a multiple of the first unit vector. That multiple goes into the
c last element of the new row of R and J is accordingly updated by the
c Givens transformations.
c
X if (nact .EQ. n) then
X work(l) = work(n)
X else
X do 160 i=n,nact+1,-1
c
c we have to find the Givens rotation which will reduce the element
c (l1) of d to zero.
c if it is already zero we don't have to do anything, except of
c decreasing l1
c
X if (work(i) .EQ. 0.d0) goto 160
X gc = max(abs(work(i-1)),abs(work(i)))
X gs = min(abs(work(i-1)),abs(work(i)))
X temp = sign(gc*sqrt(1+gs*gs/(gc*gc)), work(i-1))
X gc = work(i-1)/temp
X gs = work(i)/temp
c
c The Givens rotation is done with the matrix (gc gs, gs -gc).
c If gc is one, then element (i) of d is zero compared with element
c (l1-1). Hence we don't have to do anything.
c If gc is zero, then we just have to switch column (i) and column (i-1)
c of J. Since we only switch columns in J, we have to be careful how we
c update d depending on the sign of gs.
c Otherwise we have to apply the Givens rotation to these columns.
c The i-1 element of d has to be updated to temp.
c
X if (gc .EQ. 1.d0) goto 160
X if (gc .EQ. 0.d0) then
X work(i-1) = gs * temp
X do 170 j=1,n
X temp = dmat(j,i-1)
X dmat(j,i-1) = dmat(j,i)
X dmat(j,i) = temp
X 170 continue
X else
X work(i-1) = temp
X nu = gs/(1.d0+gc)
X do 180 j=1,n
X temp = gc*dmat(j,i-1) + gs*dmat(j,i)
X dmat(j,i) = nu*(dmat(j,i-1)+temp) - dmat(j,i)
X dmat(j,i-1) = temp
X 180 continue
X endif
X 160 continue
c
c l is still pointing to element (nact,nact) of the matrix R.
c So store d(nact) in R(nact,nact)
X work(l) = work(nact)
X endif
X else
c
c we took a partial step in dual space. Thus drop constraint it1,
c that is, we drop the it1-th active constraint.
c then we continue at step 2(a) (marked by label 55)
c but since the fit changed, we have to recalculate now "how much"
c the fit violates the chosen constraint now.
c
X sum = -bvec(nvl)
X do 190 j = 1,iamat(1,nvl)
X sum = sum + sol(iamat(j+1,nvl))*amat(j,nvl)
X 190 continue
X if( nvl .GT. meq ) then
X work(iwsv+nvl) = sum
X else
X work(iwsv+nvl) = -abs(sum)
X if( sum .GT. 0.d0) then
X do 191 j=1,iamat(1,nvl)
X amat(j,nvl) = -amat(j,nvl)
X 191 continue
X bvec(i) = -bvec(i)
X endif
X endif
X goto 700
X endif
X endif
X goto 50
c
c Drop constraint it1
c
X 700 continue
c
c if it1 = nact it is only necessary to update the vector u and nact
c
X if (it1 .EQ. nact) goto 799
c
c After updating one row of R (column of J) we will also come back here
c
X 797 continue
c
c we have to find the Givens rotation which will reduce the element
c (it1+1,it1+1) of R to zero.
c if it is already zero we don't have to do anything except of updating
c u, iact, and shifting column (it1+1) of R to column (it1)
c l will point to element (1,it1+1) of R
c l1 will point to element (it1+1,it1+1) of R
c
X l = iwrm + (it1*(it1+1))/2 + 1
X l1 = l+it1
X if (work(l1) .EQ. 0.d0) goto 798
X gc = max(abs(work(l1-1)),abs(work(l1)))
X gs = min(abs(work(l1-1)),abs(work(l1)))
X temp = sign(gc*sqrt(1+gs*gs/(gc*gc)), work(l1-1))
X gc = work(l1-1)/temp
X gs = work(l1)/temp
c
c The Givens rotatin is done with the matrix (gc gs, gs -gc).
c If gc is one, then element (it1+1,it1+1) of R is zero compared with
c element (it1,it1+1). Hence we don't have to do anything.
c if gc is zero, then we just have to switch row (it1) and row (it1+1)
c of R and column (it1) and column (it1+1) of J. Since we swithc rows in
c R and columns in J, we can ignore the sign of gs.
c Otherwise we have to apply the Givens rotation to these rows/columns.
c
X if (gc .EQ. 1.d0) goto 798
X if (gc .EQ. 0.d0) then
X do 710 i=it1+1,nact
X temp = work(l1-1)
X work(l1-1) = work(l1)
X work(l1) = temp
X l1 = l1+i
X 710 continue
X do 711 i=1,n
X temp = dmat(i,it1)
X dmat(i,it1) = dmat(i,it1+1)
X dmat(i,it1+1) = temp
X 711 continue
X else
X nu = gs/(1.d0+gc)
X do 720 i=it1+1,nact
X temp = gc*work(l1-1) + gs*work(l1)
X work(l1) = nu*(work(l1-1)+temp) - work(l1)
X work(l1-1) = temp
X l1 = l1+i
X 720 continue
X do 721 i=1,n
X temp = gc*dmat(i,it1) + gs*dmat(i,it1+1)
X dmat(i,it1+1) = nu*(dmat(i,it1)+temp) - dmat(i,it1+1)
X dmat(i,it1) = temp
X 721 continue
X endif
c
c shift column (it1+1) of R to column (it1) (that is, the first it1
c elements). The posit1on of element (1,it1+1) of R was calculated above
c and stored in l.
c
X 798 continue
X l1 = l-it1
X do 730 i=1,it1
X work(l1)=work(l)
X l = l+1
X l1 = l1+1
X 730 continue
c
c update vector u and iact as necessary
c Continue with updating the matrices J and R
c
X work(iwuv+it1) = work(iwuv+it1+1)
X iact(it1) = iact(it1+1)
X it1 = it1+1
X if (it1 .LT. nact) goto 797
X 799 work(iwuv+nact) = work(iwuv+nact+1)
X work(iwuv+nact+1) = 0.d0
X iact(nact) = 0
X nact = nact-1
X iter(2) = iter(2)+1
X goto 55
X 999 continue
X return
X end
SHAR_EOF
$shar_touch -am 0723151198 'QuadProg/solve.QP.compact.f' &&
chmod 0444 'QuadProg/solve.QP.compact.f' ||
$echo 'restore of' 'QuadProg/solve.QP.compact.f' 'failed'
if ( md5sum --help 2>&1 | grep 'sage: md5sum \[' ) >/dev/null 2>&1 \
&& ( md5sum --version 2>&1 | grep -v 'textutils 1.12' ) >/dev/null; then
md5sum -c << SHAR_EOF >/dev/null 2>&1 \
|| $echo 'QuadProg/solve.QP.compact.f:' 'MD5 check failed'
ca0858f32328ecd11d533de9a925dbe6 QuadProg/solve.QP.compact.f
SHAR_EOF
else
shar_count="`LC_ALL= LC_CTYPE= LANG= wc -c < 'QuadProg/solve.QP.compact.f'`"
test 18522 -eq "$shar_count" ||
$echo 'QuadProg/solve.QP.compact.f:' 'original size' '18522,' 'current size' "$shar_count!"
fi
fi
# ============= QuadProg/solve.QP.f ==============
if test -f 'QuadProg/solve.QP.f' && test "$first_param" != -c; then
$echo 'x -' SKIPPING 'QuadProg/solve.QP.f' '(file already exists)'
else
$echo 'x -' extracting 'QuadProg/solve.QP.f' '(text)'
sed 's/^X//' << 'SHAR_EOF' > 'QuadProg/solve.QP.f' &&
c Copyright (C) 1995, 1996, 1997, 1998
c Berwin A. Turlach
c $Id: solve.QP.f,v 1.15 1998/07/23 05:23:26 bturlach Exp $
c
c This library is free software; you can redistribute it and/or
c modify it under the terms of the GNU Library General Public
c License as published by the Free Software Foundation; either
c version 2 of the License, or (at your option) any later version.
c
c This library is distributed in the hope that it will be useful,
c but WITHOUT ANY WARRANTY; without even the implied warranty of
c MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
c Library General Public License for more details.
c
c You should have received a copy of the GNU Library General Public
c License along with this library; if not, write to the Free Software
c Foundation, Inc., 59 Temple Place, Suite 330, Boston,
c MA 02111-1307 USA
c
c This routine uses the Goldfarb/Idnani algorithm to solve the
c following minimization problem:
c
c minimize -d^T x + 1/2 * x^T D x
c where A1^T x = b1
c A2^T x >= b2
c
c the matrix D is assumed to be positive definite. Especially,
c w.l.o.g. D is assumed to be symmetric.
c
c Input parameter:
c dmat nxn matrix, the matrix D from above (dp)
c *** WILL BE DESTROYED ON EXIT ***
c The user has two possibilities:
c a) Give D (ierr=0), in this case we use routines from LINPACK
c to decompose D.
c b) To get the algorithm started we need R^-1, where D=R^TR.
c So if it is cheaper to calculate R^-1 in another way (D may
c be a band matrix) then with the general routine, the user
c may pass R^{-1}. Indicated by ierr not equal to zero.
c dvec nx1 vector, the vector d from above (dp)
c *** WILL BE DESTROYED ON EXIT ***
c contains on exit the solution to the initial, i.e.,
c unconstrained problem
c fddmat scalar, the leading dimension of the matrix dmat
c n the dimension of dmat and dvec (int)
c amat nxq matrix, the matrix A from above (dp) [ A=(A1 A2) ]
c *** ENTRIES CORRESPONDING TO EQUALITY CONSTRAINTS MAY HAVE
c CHANGED SIGNES ON EXIT ***
c bvec qx1 vector, the vector of constants b in the constraints (dp)
c [ b = (b1^T b2^T)^T ]
c *** ENTRIES CORRESPONDING TO EQUALITY CONSTRAINTS MAY HAVE
c CHANGED SIGNES ON EXIT ***
c fdamat the first dimension of amat as declared in the calling program.
c fdamat >= n !!
c q integer, the number of constraints.
c meq integer, the number of equality constraints, 0 <= meq <= q.
c ierr integer, code for the status of the matrix D:
c ierr = 0, we have to decompose D
c ierr != 0, D is already decomposed into D=R^TR and we were
c given R^{-1}.
c
c Output parameter:
c sol nx1 the final solution (x in the notation above)
c crval scalar, the value of the criterion at the minimum
c iact qx1 vector, the constraints which are active in the final
c fit (int)
c nact scalar, the number of constraints active in the final fit (int)
c iter 2x1 vector, first component gives the number of "main"
c iterations, the second one says how many constraints were
c deleted after they became active
c ierr integer, error code on exit, if
c ierr = 0, no problems
c ierr = 1, the minimization problem has no solution
c ierr = 2, problems with decomposing D, in this case sol
c contains garbage!!
c
c Working space:
c work vector with length at least 2*n+r*(r+5)/2 + 2*q +1
c where r=min(n,q)
c
X subroutine qpgen2(dmat, dvec, fddmat, n, sol, crval, amat,
X * bvec, fdamat, q, meq, iact, nact, iter, work, ierr)
X implicit none
X integer n, i, j, l, l1, fdamat, fddmat,
X * info, q, iact(*), iter(*), it1,
X * ierr, nact, iwzv, iwrv, iwrm, iwsv, iwuv, nvl,
X * r, iwnbv, meq
X double precision dmat(fddmat,*), dvec(*),sol(*), bvec(*),
X * work(*), temp, sum, t1, tt, gc, gs, crval,
X * nu, amat(fdamat,*)
X logical t1inf, t2min
X r = min(n,q)
X l = 2*n + (r*(r+5))/2 + 2*q + 1
c
c store the initial dvec to calculate below the unconstrained minima of
c the critical value.
c
X do 10 i=1,n
X work(i) = dvec(i)
X 10 continue
X do 11 i=n+1,l
X work(i) = 0.d0
X 11 continue
X do 12 i=1,q
X iact(i)=0
X 12 continue
c
c get the initial solution
c
X if( ierr .EQ. 0 )then
X call dpofa(dmat,fddmat,n,info)
X if( info .NE. 0 )then
X ierr = 2
X goto 999
X endif
X call dposl(dmat,fddmat,n,dvec)
X call dpori(dmat,fddmat,n)
X else
c
c Matrix D is already factorized, so we have to multiply d first with
c R^-T and then with R^-1. R^-1 is stored in the upper half of the
c array dmat.
c
X do 20 j=1,n
X sol(j) = 0.d0
X do 21 i=1,j
X sol(j) = sol(j) + dmat(i,j)*dvec(i)
X 21 continue
X 20 continue
X do 22 j=1,n
X dvec(j) = 0.d0
X do 23 i=j,n
X dvec(j) = dvec(j) + dmat(j,i)*sol(i)
X 23 continue
X 22 continue
X endif
c
c set lower triangular of dmat to zero, store dvec in sol and
c calculate value of the criterion at unconstrained minima
c
X crval = 0.d0
X do 30 j=1,n
X sol(j) = dvec(j)
X crval = crval + work(j)*sol(j)
X work(j) = 0.d0
X do 32 i=j+1,n
X dmat(i,j) = 0.d0
X 32 continue
X 30 continue
X crval = -crval/2.d0
X ierr = 0
c
c calculate some constants, i.e., from which index on the different
c quantities are stored in the work matrix
c
X iwzv = n
X iwrv = iwzv + n
X iwuv = iwrv + r
X iwrm = iwuv + r+1
X iwsv = iwrm + (r*(r+1))/2
X iwnbv = iwsv + q
c
c calculate the norm of each column of the A matrix
c
X do 51 i=1,q
X sum = 0.d0
X do 52 j=1,n
X sum = sum + amat(j,i)*amat(j,i)
X 52 continue
X work(iwnbv+i) = sqrt(sum)
X 51 continue
X nact = 0
X iter(1) = 0
X iter(2) = 0
X 50 continue
c
c start a new iteration
c
X iter(1) = iter(1)+1
c
c calculate all constraints and check which are still violated
c for the equality constraints we have to check whether the normal
c vector has to be negated (as well as bvec in that case)
c
X l = iwsv
X do 60 i=1,q
X l = l+1
X sum = -bvec(i)
X do 61 j = 1,n
X sum = sum + amat(j,i)*sol(j)
X 61 continue
X if (i .GT. meq) then
X work(l) = sum
X else
X work(l) = -abs(sum)
X if (sum .GT. 0.d0) then
X do 62 j=1,n
X amat(j,i) = -amat(j,i)
X 62 continue
X bvec(i) = -bvec(i)
X endif
X endif
X 60 continue
c
c as safeguard against rounding errors set already active constraints
c explicitly to zero
c
X do 70 i=1,nact
X work(iwsv+iact(i)) = 0.d0
X 70 continue
c
c we weight each violation by the number of non-zero elements in the
c corresponding row of A. then we choose the violated constraint which
c has maximal absolute value, i.e., the minimum.
c by obvious commenting and uncommenting we can choose the strategy to
c take always the first constraint which is violated. ;-)
c
X nvl = 0
X temp = 0.d0
X do 71 i=1,q
X if (work(iwsv+i) .LT. temp*work(iwnbv+i)) then
X nvl = i
X temp = work(iwsv+i)/work(iwnbv+i)
X endif
c if (work(iwsv+i) .LT. 0.d0) then
c nvl = i
c goto 72
c endif
X 71 continue
X 72 if (nvl .EQ. 0) goto 999
c
c calculate d=J^Tn^+ where n^+ is the normal vector of the violated
c constraint. J is stored in dmat in this implementation!!
c if we drop a constraint, we have to jump back here.
c
X 55 continue
X do 80 i=1,n
X sum = 0.d0
X do 81 j=1,n
X sum = sum + dmat(j,i)*amat(j,nvl)
X 81 continue
X work(i) = sum
X 80 continue
c
c Now calculate z = J_2 d_2
c
X l1 = iwzv
X do 90 i=1,n
X work(l1+i) =0.d0
X 90 continue
X do 92 j=nact+1,n
X do 93 i=1,n
X work(l1+i) = work(l1+i) + dmat(i,j)*work(j)
X 93 continue
X 92 continue
c
c and r = R^{-1} d_1, check also if r has positive elements (among the
c entries corresponding to inequalities constraints).
c
X t1inf = .TRUE.
X do 95 i=nact,1,-1
X sum = work(i)
X l = iwrm+(i*(i+3))/2
X l1 = l-i
X do 96 j=i+1,nact
X sum = sum - work(l)*work(iwrv+j)
X l = l+j
X 96 continue
X sum = sum / work(l1)
X work(iwrv+i) = sum
X if (iact(i) .LE. meq) goto 95
X if (sum .LE. 0.d0) goto 95
X 7 t1inf = .FALSE.
X it1 = i
X 95 continue
c
c if r has positive elements, find the partial step length t1, which is
c the maximum step in dual space without violating dual feasibility.
c it1 stores in which component t1, the min of u/r, occurs.
c
X if ( .NOT. t1inf) then
X t1 = work(iwuv+it1)/work(iwrv+it1)
X do 100 i=1,nact
X if (iact(i) .LE. meq) goto 100
X if (work(iwrv+i) .LE. 0.d0) goto 100
X temp = work(iwuv+i)/work(iwrv+i)
X if (temp .LT. t1) then
X t1 = temp
X it1 = i
X endif
X 100 continue
X endif
c
c test if the z vector is equal to zero
c
X sum = 0.d0
X do 110 i=iwzv+1,iwzv+n
X sum = sum + work(i)*work(i)
X 110 continue
X temp = 1000.d0
X sum = sum+temp
X if (temp .EQ. sum) then
c
c No step in pmrimal space such that the new constraint becomes
c feasible. Take step in dual space and drop a constant.
c
X if (t1inf) then
c
c No step in dual space possible either, problem is not solvable
c
X ierr = 1
X goto 999
X else
c
c we take a partial step in dual space and drop constraint it1,
c that is, we drop the it1-th active constraint.
c then we continue at step 2(a) (marked by label 55)
c
X do 111 i=1,nact
X work(iwuv+i) = work(iwuv+i) - t1*work(iwrv+i)
X 111 continue
X work(iwuv+nact+1) = work(iwuv+nact+1) + t1
X goto 700
X endif
X else
c
c compute full step length t2, minimum step in primal space such that
c the constraint becomes feasible.
c keep sum (which is z^Tn^+) to update crval below!
c
X sum = 0.d0
X do 120 i = 1,n
X sum = sum + work(iwzv+i)*amat(i,nvl)
X 120 continue
X tt = -work(iwsv+nvl)/sum
X t2min = .TRUE.
X if (.NOT. t1inf) then
X if (t1 .LT. tt) then
X tt = t1
X t2min = .FALSE.
X endif
X endif
c
c take step in primal and dual space
c
X do 130 i=1,n
X sol(i) = sol(i) + tt*work(iwzv+i)
X 130 continue
X crval = crval + tt*sum*(tt/2.d0 + work(iwuv+nact+1))
X do 131 i=1,nact
X work(iwuv+i) = work(iwuv+i) - tt*work(iwrv+i)
X 131 continue
X work(iwuv+nact+1) = work(iwuv+nact+1) + tt
c
c if it was a full step, then we check wheter further constraints are
c violated otherwise we can drop the current constraint and iterate once
c more
X if(t2min) then
c
c we took a full step. Thus add constraint nvl to the list of active
c constraints and update J and R
c
X nact = nact + 1
X iact(nact) = nvl
c
c to update R we have to put the first nact-1 components of the d vector
c into column (nact) of R
c
X l = iwrm + ((nact-1)*nact)/2 + 1
X do 150 i=1,nact-1
X work(l) = work(i)
X l = l+1
X 150 continue
c
c if now nact=n, then we just have to add the last element to the new
c row of R.
c Otherwise we use Givens transformations to turn the vector d(nact:n)
c into a multiple of the first unit vector. That multiple goes into the
c last element of the new row of R and J is accordingly updated by the
c Givens transformations.
c
X if (nact .EQ. n) then
X work(l) = work(n)
X else
X do 160 i=n,nact+1,-1
c
c we have to find the Givens rotation which will reduce the element
c (l1) of d to zero.
c if it is already zero we don't have to do anything, except of
c decreasing l1
c
X if (work(i) .EQ. 0.d0) goto 160
X gc = max(abs(work(i-1)),abs(work(i)))
X gs = min(abs(work(i-1)),abs(work(i)))
X temp = sign(gc*sqrt(1+gs*gs/(gc*gc)), work(i-1))
X gc = work(i-1)/temp
X gs = work(i)/temp
c
c The Givens rotation is done with the matrix (gc gs, gs -gc).
c If gc is one, then element (i) of d is zero compared with element
c (l1-1). Hence we don't have to do anything.
c If gc is zero, then we just have to switch column (i) and column (i-1)
c of J. Since we only switch columns in J, we have to be careful how we
c update d depending on the sign of gs.
c Otherwise we have to apply the Givens rotation to these columns.
c The i-1 element of d has to be updated to temp.
c
X if (gc .EQ. 1.d0) goto 160
X if (gc .EQ. 0.d0) then
X work(i-1) = gs * temp
X do 170 j=1,n
X temp = dmat(j,i-1)
X dmat(j,i-1) = dmat(j,i)
X dmat(j,i) = temp
X 170 continue
X else
X work(i-1) = temp
X nu = gs/(1.d0+gc)
X do 180 j=1,n
X temp = gc*dmat(j,i-1) + gs*dmat(j,i)
X dmat(j,i) = nu*(dmat(j,i-1)+temp) - dmat(j,i)
X dmat(j,i-1) = temp
X 180 continue
X endif
X 160 continue
c
c l is still pointing to element (nact,nact) of the matrix R.
c So store d(nact) in R(nact,nact)
X work(l) = work(nact)
X endif
X else
c
c we took a partial step in dual space. Thus drop constraint it1,
c that is, we drop the it1-th active constraint.
c then we continue at step 2(a) (marked by label 55)
c but since the fit changed, we have to recalculate now "how much"
c the fit violates the chosen constraint now.
c
X sum = -bvec(nvl)
X do 190 j = 1,n
X sum = sum + sol(j)*amat(j,nvl)
X 190 continue
X if( nvl .GT. meq ) then
X work(iwsv+nvl) = sum
X else
X work(iwsv+nvl) = -abs(sum)
X if( sum .GT. 0.d0) then
X do 191 j=1,n
X amat(j,nvl) = -amat(j,nvl)
X 191 continue
X bvec(i) = -bvec(i)
X endif
X endif
X goto 700
X endif
X endif
X goto 50
c
c Drop constraint it1
c
X 700 continue
c
c if it1 = nact it is only necessary to update the vector u and nact
c
X if (it1 .EQ. nact) goto 799
c
c After updating one row of R (column of J) we will also come back here
c
X 797 continue
c
c we have to find the Givens rotation which will reduce the element
c (it1+1,it1+1) of R to zero.
c if it is already zero we don't have to do anything except of updating
c u, iact, and shifting column (it1+1) of R to column (it1)
c l will point to element (1,it1+1) of R
c l1 will point to element (it1+1,it1+1) of R
c
X l = iwrm + (it1*(it1+1))/2 + 1
X l1 = l+it1
X if (work(l1) .EQ. 0.d0) goto 798
X gc = max(abs(work(l1-1)),abs(work(l1)))
X gs = min(abs(work(l1-1)),abs(work(l1)))
X temp = sign(gc*sqrt(1+gs*gs/(gc*gc)), work(l1-1))
X gc = work(l1-1)/temp
X gs = work(l1)/temp
c
c The Givens rotatin is done with the matrix (gc gs, gs -gc).
c If gc is one, then element (it1+1,it1+1) of R is zero compared with
c element (it1,it1+1). Hence we don't have to do anything.
c if gc is zero, then we just have to switch row (it1) and row (it1+1)
c of R and column (it1) and column (it1+1) of J. Since we swithc rows in
c R and columns in J, we can ignore the sign of gs.
c Otherwise we have to apply the Givens rotation to these rows/columns.
c
X if (gc .EQ. 1.d0) goto 798
X if (gc .EQ. 0.d0) then
X do 710 i=it1+1,nact
X temp = work(l1-1)
X work(l1-1) = work(l1)
X work(l1) = temp
X l1 = l1+i
X 710 continue
X do 711 i=1,n
X temp = dmat(i,it1)
X dmat(i,it1) = dmat(i,it1+1)
X dmat(i,it1+1) = temp
X 711 continue
X else
X nu = gs/(1.d0+gc)
X do 720 i=it1+1,nact
X temp = gc*work(l1-1) + gs*work(l1)
X work(l1) = nu*(work(l1-1)+temp) - work(l1)
X work(l1-1) = temp
X l1 = l1+i
X 720 continue
X do 721 i=1,n
X temp = gc*dmat(i,it1) + gs*dmat(i,it1+1)
X dmat(i,it1+1) = nu*(dmat(i,it1)+temp) - dmat(i,it1+1)
X dmat(i,it1) = temp
X 721 continue
X endif
c
c shift column (it1+1) of R to column (it1) (that is, the first it1
c elements). The posit1on of element (1,it1+1) of R was calculated above
c and stored in l.
c
X 798 continue
X l1 = l-it1
X do 730 i=1,it1
X work(l1)=work(l)
X l = l+1
X l1 = l1+1
X 730 continue
c
c update vector u and iact as necessary
c Continue with updating the matrices J and R
c
X work(iwuv+it1) = work(iwuv+it1+1)
X iact(it1) = iact(it1+1)
X it1 = it1+1
X if (it1 .LT. nact) goto 797
X 799 work(iwuv+nact) = work(iwuv+nact+1)
X work(iwuv+nact+1) = 0.d0
X iact(nact) = 0
X nact = nact-1
X iter(2) = iter(2)+1
X goto 55
X 999 continue
X return
X end
SHAR_EOF
$shar_touch -am 0723151198 'QuadProg/solve.QP.f' &&
chmod 0444 'QuadProg/solve.QP.f' ||
$echo 'restore of' 'QuadProg/solve.QP.f' 'failed'
if ( md5sum --help 2>&1 | grep 'sage: md5sum \[' ) >/dev/null 2>&1 \
&& ( md5sum --version 2>&1 | grep -v 'textutils 1.12' ) >/dev/null; then
md5sum -c << SHAR_EOF >/dev/null 2>&1 \
|| $echo 'QuadProg/solve.QP.f:' 'MD5 check failed'
f6cea73038b863d6dfcc23d817ff0edd QuadProg/solve.QP.f
SHAR_EOF
else
shar_count="`LC_ALL= LC_CTYPE= LANG= wc -c < 'QuadProg/solve.QP.f'`"
test 17890 -eq "$shar_count" ||
$echo 'QuadProg/solve.QP.f:' 'original size' '17890,' 'current size' "$shar_count!"
fi
fi
# ============= QuadProg/util.f ==============
if test -f 'QuadProg/util.f' && test "$first_param" != -c; then
$echo 'x -' SKIPPING 'QuadProg/util.f' '(file already exists)'
else
$echo 'x -' extracting 'QuadProg/util.f' '(text)'
sed 's/^X//' << 'SHAR_EOF' > 'QuadProg/util.f' &&
X subroutine dpori(a,lda,n)
X integer lda,n
X double precision a(lda,1)
c
c dpori computes the inverse of the factor of a
c double precision symmetric positive definite matrix
c using the factors computed by dpofa.
c
c modification of dpodi by BaT 05/11/95
c
c on entry
c
c a double precision(lda, n)
c the output a from dpofa
c
c lda integer
c the leading dimension of the array a .
c
c n integer
c the order of the matrix a .
c
c on return
c
c a if dpofa was used to factor a then
c dpodi produces the upper half of inverse(a) .
c elements of a below the diagonal are unchanged.
c
c error condition
c
c a division by zero will occur if the input factor contains
c a zero on the diagonal and the inverse is requested.
c it will not occur if the subroutines are called correctly
c and if dpoco or dpofa has set info .eq. 0 .
c
c linpack. this version dated 08/14/78 .
c cleve moler, university of new mexico, argonne national lab.
c modified by Berwin A. Turlach 05/11/95
c
c subroutines and functions
c
c blas daxpy,dscal
c fortran mod
c
c internal variables
c
X double precision t
X integer j,k,kp1
c
c compute inverse(r)
c
X do 100 k = 1, n
X a(k,k) = 1.0d0/a(k,k)
X t = -a(k,k)
X call dscal(k-1,t,a(1,k),1)
X kp1 = k + 1
X if (n .lt. kp1) go to 90
X do 80 j = kp1, n
X t = a(k,j)
X a(k,j) = 0.0d0
X call daxpy(k,t,a(1,k),1,a(1,j),1)
X 80 continue
X 90 continue
X 100 continue
X return
X end
X
X subroutine dposl(a,lda,n,b)
X integer lda,n
X double precision a(lda,1),b(1)
c
c dposl solves the double precision symmetric positive definite
c system a * x = b
c using the factors computed by dpoco or dpofa.
c
c on entry
c
c a double precision(lda, n)
c the output from dpoco or dpofa.
c
c lda integer
c the leading dimension of the array a .
c
c n integer
c the order of the matrix a .
c
c b double precision(n)
c the right hand side vector.
c
c on return
c
c b the solution vector x .
c
c error condition
c
c a division by zero will occur if the input factor contains
c a zero on the diagonal. technically this indicates
c singularity but it is usually caused by improper subroutine
c arguments. it will not occur if the subroutines are called
c correctly and info .eq. 0 .
c
c to compute inverse(a) * c where c is a matrix
c with p columns
c call dpoco(a,lda,n,rcond,z,info)
c if (rcond is too small .or. info .ne. 0) go to ...
c do 10 j = 1, p
c call dposl(a,lda,n,c(1,j))
c 10 continue
c
c linpack. this version dated 08/14/78 .
c cleve moler, university of new mexico, argonne national lab.
c
c subroutines and functions
c
c blas daxpy,ddot
c
c internal variables
c
X double precision ddot,t
X integer k,kb
c
c solve trans(r)*y = b
c
X do 10 k = 1, n
X t = ddot(k-1,a(1,k),1,b(1),1)
X b(k) = (b(k) - t)/a(k,k)
X 10 continue
c
c solve r*x = y
c
X do 20 kb = 1, n
X k = n + 1 - kb
X b(k) = b(k)/a(k,k)
X t = -b(k)
X call daxpy(k-1,t,a(1,k),1,b(1),1)
X 20 continue
X return
X end
SHAR_EOF
$shar_touch -am 0723151198 'QuadProg/util.f' &&
chmod 0444 'QuadProg/util.f' ||
$echo 'restore of' 'QuadProg/util.f' 'failed'
if ( md5sum --help 2>&1 | grep 'sage: md5sum \[' ) >/dev/null 2>&1 \
&& ( md5sum --version 2>&1 | grep -v 'textutils 1.12' ) >/dev/null; then
md5sum -c << SHAR_EOF >/dev/null 2>&1 \
|| $echo 'QuadProg/util.f:' 'MD5 check failed'
bbf1c42d679725a411a775e4f97074e0 QuadProg/util.f
SHAR_EOF
else
shar_count="`LC_ALL= LC_CTYPE= LANG= wc -c < 'QuadProg/util.f'`"
test 3581 -eq "$shar_count" ||
$echo 'QuadProg/util.f:' 'original size' '3581,' 'current size' "$shar_count!"
fi
fi
# ============= QuadProg/solve.QP.compact.d ==============
if test -f 'QuadProg/solve.QP.compact.d' && test "$first_param" != -c; then
$echo 'x -' SKIPPING 'QuadProg/solve.QP.compact.d' '(file already exists)'
else
$echo 'x -' extracting 'QuadProg/solve.QP.compact.d' '(text)'
sed 's/^X//' << 'SHAR_EOF' > 'QuadProg/solve.QP.compact.d' &&
X.\" Copyright (C) 1995, 1996, 1997, 1998
X.\" Berwin A. Turlach
X.\" $Id: solve.QP.compact.d,v 1.7 1998/07/23 05:11:59 bturlach Exp $
X.\"
X.\" This library is free software; you can redistribute it and/or
X.\" modify it under the terms of the GNU Library General Public
X.\" License as published by the Free Software Foundation; either
X.\" version 2 of the License, or (at your option) any later version.
X.\"
X.\" This library is distributed in the hope that it will be useful,
X.\" but WITHOUT ANY WARRANTY; without even the implied warranty of
X.\" MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
X.\" Library General Public License for more details.
X.\"
X.\" You should have received a copy of the GNU Library General Public
X.\" License along with this library; if not, write to the Free Software
X.\" Foundation, Inc., 59 Temple Place, Suite 330, Boston,
X.\" MA 02111-1307 USA
X.\"
X.BG
X.FN solve.QP.compact
X.TL
Solve a Quadratic Programming Problem
X.DN
Routine to solve the following Quadratic Programming problem:
X minimize -d^T b + 1/2 b^T D b
X where A^T b >= b0
X
Unlike solve.QP, this routine expects the matrix A in a compact
storage form.
X.CS
solve.QP.compact(Dmat, dvec, Amat, Aind, bvec, meq=0, factorized=F)
X.RA
X.AG Dmat
matrix appearing in the quadratic function to be minimized.
X.AG dvec
vector appearing in the quadratic function to be minimized.
X.AG Amat
matrix containing the non-zero elements of the matrix A that defines
the constraints. If mi denotes the number of non-zero elements in the
i-th column of A then the first mi entries of the i-th column of
`Amat' hold these non-zero elements.
(If maxmi denoted the maximum of all mi, then each column of `Amat'
may have arbitrary elements from row mi+1 to row maxmi in the i-th
column)
X.AG Aind
matrix of integers. The first element of each column gives the number of
non-zero elements in the corresponding column of the matrix A. The
following entries in each column contain the indexes of the rows in
which these non-zero elements are.
X.OA
X.AG bvec
vector holding the values of b0 (defaults to zero).
X.AG meq
the first `meq' constraints are treated as equality constraints,
all further as inequality constraints (defaults to 0).
X.AG factorized
logical flag: if `TRUE', then we are passing
R^-1 (where D = R^T R) instead of the matrix D in the argument `Dmat'.
X.RT
a list with the following components:
X.AG solution
vector containing the solution of the quadratic programming problem.
X.AG value
scalar, the value of the quadratic function at the solution
X.AG unconstrained.solution
vector containing the unconstrained minimizer of the quadratic function.
X.AG iterations
vector of length 2, the first component contains the number of
iterations the algorithm needed, the second indicates how often
constraints became inactive after becoming active first.
X.AG iact
vector with the indices of the active constraints at the solution.
X.DT
This routine implements the dual method of Goldfarb and Idnani (1982,
1983) for solving quadratic programming problems.
X.SH REFERENCES
Goldfarb, D. and Idnani, A. (1982).
Dual and Primal-Dual Methods for Solving Strictly Convex
Quadratic Programs. In
X.ul
Numerical Analysis
J.P. Hennart, ed. Springer-Verlag, Berlin. pp. 226-239.
X.PP
Goldfarb, D. and Idnani, A. (1983).
A numerically stable dual method for solving strictly convex quadratic
programs.
X.ul
Mathematical Programming
\fB27\fR, 1-33.
X.SA
`solve.QP'
X.EX
#
# Assume we want to minimize: -(0 5 0) %*% b + 1/2 b^T b
# under the constraints: A^T b >= b0
# with b0 = (-8,2,0)^T
# and (-4 2 0)
# A = (-3 1 -2)
# ( 0 0 1)
# we can use solve.QP.compact as follows:
#
Dmat <- matrix(0,3,3)
diag(Dmat) <- 1
dvec <- c(0,5,0)
Aind <- c(2,2,2)
Aind <- rbind(Aind,c(1,1,2))
Aind <- rbind(Aind,c(2,2,3))
Amat <- c(-4,2,-2)
Amat <- rbind(Amat,c(-3,1,1))
bvec <- c(-8,2,0)
solve.QP.compact(Dmat,dvec,Amat,Aind,bvec=bvec)
X
X.KW optimize
X.WR
SHAR_EOF
$shar_touch -am 0723151198 'QuadProg/solve.QP.compact.d' &&
chmod 0444 'QuadProg/solve.QP.compact.d' ||
$echo 'restore of' 'QuadProg/solve.QP.compact.d' 'failed'
if ( md5sum --help 2>&1 | grep 'sage: md5sum \[' ) >/dev/null 2>&1 \
&& ( md5sum --version 2>&1 | grep -v 'textutils 1.12' ) >/dev/null; then
md5sum -c << SHAR_EOF >/dev/null 2>&1 \
|| $echo 'QuadProg/solve.QP.compact.d:' 'MD5 check failed'
ab25f6dfaaa8cb3e48e9c439bcc08a77 QuadProg/solve.QP.compact.d
SHAR_EOF
else
shar_count="`LC_ALL= LC_CTYPE= LANG= wc -c < 'QuadProg/solve.QP.compact.d'`"
test 4034 -eq "$shar_count" ||
$echo 'QuadProg/solve.QP.compact.d:' 'original size' '4034,' 'current size' "$shar_count!"
fi
fi
# ============= QuadProg/solve.QP.d ==============
if test -f 'QuadProg/solve.QP.d' && test "$first_param" != -c; then
$echo 'x -' SKIPPING 'QuadProg/solve.QP.d' '(file already exists)'
else
$echo 'x -' extracting 'QuadProg/solve.QP.d' '(text)'
sed 's/^X//' << 'SHAR_EOF' > 'QuadProg/solve.QP.d' &&
X.\" Copyright (C) 1995, 1996, 1997, 1998
X.\" Berwin A. Turlach
X.\" $Id: solve.QP.d,v 1.7 1998/07/23 05:12:20 bturlach Exp $
X.\"
X.\" This library is free software; you can redistribute it and/or
X.\" modify it under the terms of the GNU Library General Public
X.\" License as published by the Free Software Foundation; either
X.\" version 2 of the License, or (at your option) any later version.
X.\"
X.\" This library is distributed in the hope that it will be useful,
X.\" but WITHOUT ANY WARRANTY; without even the implied warranty of
X.\" MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
X.\" Library General Public License for more details.
X.\"
X.\" You should have received a copy of the GNU Library General Public
X.\" License along with this library; if not, write to the Free Software
X.\" Foundation, Inc., 59 Temple Place, Suite 330, Boston,
X.\" MA 02111-1307 USA
X.\"
X.BG
X.FN solve.QP
X.TL
Solve a Quadratic Programming Problem
X.DN
Routine to solve the following Quadratic Programming problem:
X minimize -d^T b + 1/2 b^T D b
X where A^T b >= b0
X.CS
solve.QP(Dmat, dvec, Amat, bvec, meq=0, factorized=F)
X.RA
X.AG Dmat
matrix appearing in the quadratic function to be minimized.
X.AG dvec
vector appearing in the quadratic function to be minimized.
X.AG Amat
matrix defining the constraints under which we want to
minimize the quadratic function.
X.OA
X.AG bvec
vector holding the values of b0 (defaults to zero).
X.AG meq
the first `meq' constraints are treated as equality constraints,
all further as inequality constraints (defaults to 0).
X.AG factorized
logical flag: if `TRUE', then we are passing
R^-1 (where D = R^T R) instead of the matrix D in the argument `Dmat'.
X.RT
a list with the following components:
X.AG solution
vector containing the solution of the quadratic programming problem.
X.AG value
scalar, the value of the quadratic function at the solution
X.AG unconstrained.solution
vector containing the unconstrained minimizer of the quadratic function.
X.AG iterations
vector of length 2, the first component contains the number of
iterations the algorithm needed, the second indicates how often
constraints became inactive after becoming active first.
X.AG iact
vector with the indices of the active constraints at the solution.
X.DT
This routine implements the dual method of Goldfarb and Idnani (1982,
1983) for solving quadratic programming problems.
X.SH REFERENCES
Goldfarb, D. and Idnani, A. (1982).
Dual and Primal-Dual Methods for Solving Strictly Convex
Quadratic Programs. In
X.ul
Numerical Analysis
J.P. Hennart, ed. Springer-Verlag, Berlin. pp. 226-239.
X.PP
Goldfarb, D. and Idnani, A. (1983).
A numerically stable dual method for solving strictly convex quadratic
programs.
X.ul
Mathematical Programming
\fB27\fR, 1-33.
X.SA
`solve.QP.compact'
X.EX
#
# Assume we want to minimize: -(0 5 0) %*% b + 1/2 b^T b
# under the constraints: A^T b >= b0
# with b0 = (-8,2,0)^T
# and (-4 2 0)
# A = (-3 1 -2)
# ( 0 0 1)
# we can use solve.QP as follows:
#
Dmat <- matrix(0,3,3)
diag(Dmat) <- 1
dvec <- c(0,5,0)
Amat <- matrix(c(-4,-3,0,2,1,0,0,-2,1),3,3)
bvec <- c(-8,2,0)
solve.QP(Dmat,dvec,Amat,bvec=bvec)
X
X.KW optimize
X.WR
SHAR_EOF
$shar_touch -am 0723151198 'QuadProg/solve.QP.d' &&
chmod 0444 'QuadProg/solve.QP.d' ||
$echo 'restore of' 'QuadProg/solve.QP.d' 'failed'
if ( md5sum --help 2>&1 | grep 'sage: md5sum \[' ) >/dev/null 2>&1 \
&& ( md5sum --version 2>&1 | grep -v 'textutils 1.12' ) >/dev/null; then
md5sum -c << SHAR_EOF >/dev/null 2>&1 \
|| $echo 'QuadProg/solve.QP.d:' 'MD5 check failed'
e7070ef9905a3e426130a67b8b4b2447 QuadProg/solve.QP.d
SHAR_EOF
else
shar_count="`LC_ALL= LC_CTYPE= LANG= wc -c < 'QuadProg/solve.QP.d'`"
test 3255 -eq "$shar_count" ||
$echo 'QuadProg/solve.QP.d:' 'original size' '3255,' 'current size' "$shar_count!"
fi
fi
# ============= QuadProg/QuadProg.q ==============
if test -f 'QuadProg/QuadProg.q' && test "$first_param" != -c; then
$echo 'x -' SKIPPING 'QuadProg/QuadProg.q' '(file already exists)'
else
$echo 'x -' extracting 'QuadProg/QuadProg.q' '(text)'
sed 's/^X//' << 'SHAR_EOF' > 'QuadProg/QuadProg.q' &&
### Copyright (C) 1995, 1996, 1997, 1998
### Berwin A. Turlach
### $Id: First.lib.template,v 1.11 1998/07/23 05:13:14 bturlach Exp $
###
### This library is free software; you can redistribute it and/or
### modify it under the terms of the GNU Library General Public
### License as published by the Free Software Foundation; either
### version 2 of the License, or (at your option) any later version.
###
### This library is distributed in the hope that it will be useful,
### but WITHOUT ANY WARRANTY; without even the implied warranty of
### MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
### Library General Public License for more details.
###
### You should have received a copy of the GNU Library General Public
### License along with this library; if not, write to the Free Software
### Foundation, Inc., 59 Temple Place, Suite 330, Boston,
### MA 02111-1307 USA
###
".First.lib"<-
X function(a, b)
{
X X <- paste(a, "/", b, "/", "QuadProg_l.o", sep = "")
X dyn.load(X)
X cat("Library of quadratic programming routines. Release 1.4\n")
X cat("Copyright (C) 1995, 1996, 1997, 1998\n")
X cat("Berwin A. Turlach \n\n")
}
X
### Copyright (C) 1995, 1996, 1997, 1998
### Berwin A. Turlach
### $Id: solve.QP.compact.sf,v 1.12 1998/07/23 05:10:29 bturlach Exp $
###
### This library is free software; you can redistribute it and/or
### modify it under the terms of the GNU Library General Public
### License as published by the Free Software Foundation; either
### version 2 of the License, or (at your option) any later version.
###
### This library is distributed in the hope that it will be useful,
### but WITHOUT ANY WARRANTY; without even the implied warranty of
### MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
### Library General Public License for more details.
###
### You should have received a copy of the GNU Library General Public
### License along with this library; if not, write to the Free Software
### Foundation, Inc., 59 Temple Place, Suite 330, Boston,
### MA 02111-1307 USA
###
solve.QP.compact <- function(Dmat, dvec, Amat, Aind, bvec, meq=0,
X factorized=F){
X n <- nrow(Dmat)
X q <- ncol(Amat)
X anrow <- nrow(Amat)
X if( missing(bvec) )
X bvec <- rep(0,q)
X
X if( n != ncol(Dmat) )
X stop("Dmat is not symmetric!")
X if( n != length(dvec) )
X stop("Dmat and dvec are incompatible!")
X if( (anrow+1 != nrow(Aind)) || (q != ncol(Aind)) || (q != length(bvec)) )
X stop("Amat, Aind and bvec are incompatible!")
X Aindok <- .Fortran("aind",
X as.integer(Aind), as.integer(anrow+1),
X as.integer(q), as.integer(n),
X ok=as.logical(T))$ok
X if( !Aindok )
X stop("Aind contains illegal indexes!")
X if( (meq > q) || (meq < 0 ) )
X stop("Value of meq is invalid!")
X
X iact <- rep(0,q)
X nact <- 0
X r <- min(n,q)
X sol <- rep(0,n)
X crval <- 0
X work <- rep(0,2*n+r*(r+5)/2+2*q+1)
X iter <- rep(0,2)
X
X res1 <- .Fortran("qpgen1",
X as.double(Dmat), dvec=as.double(dvec),
X as.integer(n), as.integer(n),
X sol=as.double(sol), crval=as.double(crval),
X as.double(Amat), as.integer(Aind), as.double(bvec),
X as.integer(anrow), as.integer(q), as.integer(meq),
X iact=as.integer(iact), nact=as.integer(nact),
X iter=as.integer(iter),
X work=as.double(work), ierr=as.integer(factorized))
X
X if( res1$ierr == 1)
X stop("constraints are inconsistent, no solution!")
X else if( res1$ierr == 2)
X stop("matrix D in quadratic function is not positive definite!")
X
X prec <- .Machine$double.eps^0.25 ## this should be achievable
X check <- .Fortran("achck",
X as.double(res1$sol), as.integer(n),
X as.double(Amat), as.integer(Aind), as.double(bvec),
X as.integer(anrow), as.integer(q), as.integer(meq),
X as.double(prec), ok=as.logical(T))$ok
X if( !check )
X stop("problem is ill-conditioned, cannot find solution!")
X
X
X list(solution=res1$sol,
X value=res1$crval,
X unconstrained.solution=res1$dvec,
X iterations=res1$iter,
X iact=res1$iact[1:res1$nact])
}
### Copyright (C) 1995, 1996, 1997, 1998
### Berwin A. Turlach
### $Id: solve.QP.sf,v 1.13 1998/07/23 05:10:03 bturlach Exp $
###
### This library is free software; you can redistribute it and/or
### modify it under the terms of the GNU Library General Public
### License as published by the Free Software Foundation; either
### version 2 of the License, or (at your option) any later version.
###
### This library is distributed in the hope that it will be useful,
### but WITHOUT ANY WARRANTY; without even the implied warranty of
### MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
### Library General Public License for more details.
###
### You should have received a copy of the GNU Library General Public
### License along with this library; if not, write to the Free Software
### Foundation, Inc., 59 Temple Place, Suite 330, Boston,
### MA 02111-1307 USA
###
solve.QP <- function(Dmat, dvec, Amat, bvec, meq=0, factorized=F){
X
X n <- nrow(Dmat)
X q <- ncol(Amat)
X if( missing(bvec) )
X bvec <- rep(0,q)
X
X if( n != ncol(Dmat) )
X stop("Dmat is not symmetric!")
X if( n != length(dvec) )
X stop("Dmat and dvec are incompatible!")
X if( n != nrow(Amat) )
X stop("Amat and dvec are incompatible!")
X if( q != length(bvec) )
X stop("Amat and bvec are incompatible!")
X if( (meq > q) || (meq < 0 ) )
X stop("Value of meq is invalid!")
X
X iact <- rep(0,q)
X nact <- 0
X r <- min(n,q)
X sol <- rep(0,n)
X crval <- 0
X work <- rep(0,2*n+r*(r+5)/2+2*q+1)
X iter <- rep(0,2)
X
X res1 <- .Fortran("qpgen2",
X as.double(Dmat), dvec=as.double(dvec),
X as.integer(n), as.integer(n),
X sol=as.double(sol), crval=as.double(crval),
X as.double(Amat), as.double(bvec), as.integer(n),
X as.integer(q), as.integer(meq),
X iact=as.integer(iact), nact=as.integer(nact),
X iter=as.integer(iter), work=as.double(work),
X ierr=as.integer(factorized))
X
X if( res1$ierr == 1)
X stop("constraints are inconsistent, no solution!")
X else if( res1$ierr == 2)
X stop("matrix D in quadratic function is not positive definite!")
X
X check <- t(Amat)%*%res1$sol-bvec
X if( meq > 0 )
X check[1:meq] <- -abs(check[1:meq])
X prec <- .Machine$double.eps^0.25 ## this should be achievable
X if( any(check < -prec) )
X stop("problem is ill-conditioned, cannot find solution!")
X
X
X list(solution=res1$sol,
X value=res1$crval,
X unconstrainted.solution=res1$dvec,
X iterations=res1$iter,
X iact=res1$iact[1:res1$nact])
}
SHAR_EOF
$shar_touch -am 0723151198 'QuadProg/QuadProg.q' &&
chmod 0644 'QuadProg/QuadProg.q' ||
$echo 'restore of' 'QuadProg/QuadProg.q' 'failed'
if ( md5sum --help 2>&1 | grep 'sage: md5sum \[' ) >/dev/null 2>&1 \
&& ( md5sum --version 2>&1 | grep -v 'textutils 1.12' ) >/dev/null; then
md5sum -c << SHAR_EOF >/dev/null 2>&1 \
|| $echo 'QuadProg/QuadProg.q:' 'MD5 check failed'
b3b872a38c6030c543629d3a6220f941 QuadProg/QuadProg.q
SHAR_EOF
else
shar_count="`LC_ALL= LC_CTYPE= LANG= wc -c < 'QuadProg/QuadProg.q'`"
test 7036 -eq "$shar_count" ||
$echo 'QuadProg/QuadProg.q:' 'original size' '7036,' 'current size' "$shar_count!"
fi
fi
# ============= QuadProg/Makefile ==============
if test -f 'QuadProg/Makefile' && test "$first_param" != -c; then
$echo 'x -' SKIPPING 'QuadProg/Makefile' '(file already exists)'
else
$echo 'x -' extracting 'QuadProg/Makefile' '(text)'
sed 's/^X//' << 'SHAR_EOF' > 'QuadProg/Makefile' &&
# Makefile created by S-PLUS utility CHAPTER, version 3.25
# S-PLUS Version 3.4 Release 1 for Sun SPARC, SunOS 4.1.3_U1 : 1996
X
SHOME=/usr/local/splus3.4
X
#WHICH_LOAD=static.load
WHICH_LOAD=dyn.load
#WHICH_LOAD=dyn.load.shared
X
EXTRA_OBJ_FILES=
EXTRA_C_NAMES=
EXTRA_F_NAMES=
X
# This is a Makefile produced by the S-PLUS utility CHAPTER. It guides
# compilation of C, Ratfor and Fortran source code, loading of the resulting
# object code, and installation of functions and helpfiles.
#
# Overview:
# ========
#
# You will need to attend to a few 'make' macro settings at the top of this
# file first. See the Detailed Instructions below for these. Then, if you
# want to compile and load the default set of C and Fortran routines, type
# (1) 'make install' to install functions on .Data and thereby enable
# building of the default list of compiled routine names, then (2) 'make
# load' to compile the object code, build the set of names to load, and load
# this object code. Other, secondary, 'make' targets are 'install.funs' and
# 'install.help' ('install' makes both of those), 'clean' (remove only
# intermediate object code) and 'virgin' (remove object code, the local
# standalone S-PLUS executable if you built it, and the directory
# .Data).
#
# Detailed Instructions:
# =====================
#
# 1. On the line near the top that starts with "SHOME=", change the
# directory pathname on the right of the equals sign to the pathname of the
# top S-PLUS directory on your machine; use Splus SHOME to get this
# pathname.
#
# 2. Near the top, the macro WHICH_LOAD specifies whether 'make' targets
# 'load' and 'all' should make
# (a) a standalone S-PLUS binary executable named
# "local.Sqpe" which includes the local C and Fortran code, or
# (b) a dyn.loadable file "QuadProg_l.o" containing only this local code, or
# (c) a shared library "QuadProg.so" containing this local code
# and loadable with the dyn.load.shared function (this is the only form
# of dynamic loading available on Splus 3.3 for the Iris and Dec Alpha).
# To make a complete standalone executable, WHICH_LOAD should be set to
# "static.load"; the other possible values are "dyn.load" and
# "dyn.load.shared". Make sure a # comment symbol (#) appears in front
# of the setting you DON'T want.
# The default has been set to "static.load".
#
# Incidentally, you can always make targets 'static.load', 'dyn.load',
# and 'dyn.load.shared' regardless of the value of WHICH_LOAD.
# The macro just specifies what the generic targets 'load' and 'all'
# should make.
#
# 3. Near the top, the macro EXTRA_OBJ_FILES names object files (ending in
# ".o") which should be loaded in addition to those which are part of this
# chapter. Add any additional object file names you want to load, separated
# by spaces, on the right hand side of the equals sign. There MUST be
# corresponding C (".c"), Ratfor (".r") or Fortran (".f") files present in
# this directory.
#
# If you add file names to EXTRA_OBJ_FILES, you must also add the names of
# the corresponding C functions or Fortran subroutines to the macros
# EXTRA_C_NAMES and EXTRA_F_NAMES respectively, underneath EXTRA_OBJ_FILES.
# The names added should be as they appear at the source-code level; do not
# add leading or trailing underscores unless these are present in the
# source-level name.
#
# You won't generally need to change any more 'make' macros, or indeed
# anything else in this Makefile.
#
# 4. The file QuadProg_i.c is a C language source file which will ensure
# that the make target 'load' will load C and Fortran routines referenced
# only through .C or .Fortran calls in S functions. If QuadProg_i.c does
# not exist, it will be made during 'make load' based on the objects in
# .Data. If QuadProg_i.c exists already and you have since added new C
# or Fortran routine names to EXTRA_C_NAMES and EXTRA_F_NAMES, you must
# remove QuadProg_i.c so that it will get remade using this new
# information. Do this now.
#
# 5. Type 'make install.funs' to create the S functions on .Data. This
# must precede 'make load' in Step 6. (You can install the helpfiles at
# this point as well by typing making 'install' instead of 'install.funs'.)
#
# 6. Type 'make load'. This will build the list of routines to load, then
# compile the C, Ratfor and Fortran source files specified by this Makefile,
# and lastly load the routines named on the list into either a standalone
# S-PLUS executable or a dyn.loadable file or a shared library.
#
# 7. To remove unnecessary object files, type 'make clean'. This will not
# remove the dyn.loadable file QuadProg_l.o, the shared library
# QuadProg.so, or the standalone local.Sqpe. To clean out everything
# including .Data and start over, use 'make virgin'.
#
# 8. See the helpfile for the "library" function in S-PLUS for hints on
# how to install the new code as a library section.
X
X
# ========================= End instructions. =============================
X
default : all
include $(SHOME)/newfun/lib/S_makefile
chapters=../QuadProg
SRC= achck.f aind.f solve.QP.compact.f solve.QP.f util.f
OBJ= achck.o aind.o solve.QP.compact.o solve.QP.o util.o $(EXTRA_OBJ_FILES)
CFLAGS=-Xs -O2
FFLAGS=-O2
FUNS=force_ld.q QuadProg.q
HELPS= solve.QP.compact.d solve.QP.d
FLAGS=
RM=-rm
X
all load: $(WHICH_LOAD)
X
static.load: QuadProg.a
X Splus LOAD $(FLAGS) CHAPTERS='"$(chapters)"'
X
dyn.load QuadProg_l.o: QuadProg.a
X ld -r -o QuadProg_l.o QuadProg_i.o QuadProg.a $(FLAGS)
X @echo dynamically loadable file in QuadProg_l.o
X
dyn.load.shared QuadProg.so: $(SRC)
X Splus SHLIB -o QuadProg.so $(SRC)
X
QuadProg.a: QuadProg_i.o $(OBJ)
X Splus LIBRARY QuadProg.a QuadProg_i.o $(OBJ)
X
QuadProg_i.c:
X -mkdir .Data
X Splus make.init QuadProg .Data
X
install : install.funs install.help
X
funs install.funs : $(FUNS) .Data
X Splus QINSTALL .Data $(FUNS)
X
force_ld.q :
X ( echo "force.loading <- function(){" ;\
X set $(EXTRA_C_NAMES) terminator ;\
X while test $$1 != terminator ;\
X do \
X echo $$1 | sed 's:.*:.C("&"):' ;\
X shift ;\
X done ;\
X set $(EXTRA_F_NAMES) terminator ;\
X while test $$1 != terminator ;\
X do \
X echo $$1 | sed 's:.*:.Fortran("&"):' ;\
X shift ;\
X done ;\
X echo } \
X ) > $@
X
X
# force_ld.q :
# echo "force.loading <- function(){ X\
# stop(\"should not be executed\") X\
# # Add .C and .Fortran calls here to force routines to X\
# # be loaded. For example, .C(\"fun1\"), .Fortran(\"fun2\"). X\
# }" | tr X \\012 > $@
X
help install.help : $(HELPS) .Data/.Help
X Splus HINSTALL .Data/.Help $(HELPS)
X Splus help.findsum .Data
X
X.Data :
X -mkdir $@
X.Data/.Help : .Data
X -mkdir $@
X
virgin : clean virgin.std
# add additional cleanup rules/targets above here for target 'virgin'
clean :
X $(RM) -f $(OBJ) QuadProg_i.o S_load_time.[oc] core
virgin.std :
X $(RM) -rf .Data
X $(RM) -f QuadProg.a local.Sqpe QuadProg_l.o QuadProg.so QuadProg_i.c force_ld.q
SHAR_EOF
$shar_touch -am 0723151198 'QuadProg/Makefile' &&
chmod 0644 'QuadProg/Makefile' ||
$echo 'restore of' 'QuadProg/Makefile' 'failed'
if ( md5sum --help 2>&1 | grep 'sage: md5sum \[' ) >/dev/null 2>&1 \
&& ( md5sum --version 2>&1 | grep -v 'textutils 1.12' ) >/dev/null; then
md5sum -c << SHAR_EOF >/dev/null 2>&1 \
|| $echo 'QuadProg/Makefile:' 'MD5 check failed'
f1a1a30636e4fe41eda63a4eeecbae32 QuadProg/Makefile
SHAR_EOF
else
shar_count="`LC_ALL= LC_CTYPE= LANG= wc -c < 'QuadProg/Makefile'`"
test 6834 -eq "$shar_count" ||
$echo 'QuadProg/Makefile:' 'original size' '6834,' 'current size' "$shar_count!"
fi
fi
# ============= QuadProg/COPYING.LIB ==============
if test -f 'QuadProg/COPYING.LIB' && test "$first_param" != -c; then
$echo 'x -' SKIPPING 'QuadProg/COPYING.LIB' '(file already exists)'
else
$echo 'x -' extracting 'QuadProg/COPYING.LIB' '(text)'
sed 's/^X//' << 'SHAR_EOF' > 'QuadProg/COPYING.LIB' &&
X GNU LIBRARY GENERAL PUBLIC LICENSE
X Version 2, June 1991
X
X Copyright (C) 1991 Free Software Foundation, Inc.
X 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
X Everyone is permitted to copy and distribute verbatim copies
X of this license document, but changing it is not allowed.
X
[This is the first released version of the library GPL. It is
X numbered 2 because it goes with version 2 of the ordinary GPL.]
X
X Preamble
X
X The licenses for most software are designed to take away your
freedom to share and change it. By contrast, the GNU General Public
Licenses are intended to guarantee your freedom to share and change
free software--to make sure the software is free for all its users.
X
X This license, the Library General Public License, applies to some
specially designated Free Software Foundation software, and to any
other libraries whose authors decide to use it. You can use it for
your libraries, too.
X
X When we speak of free software, we are referring to freedom, not
price. Our General Public Licenses are designed to make sure that you
have the freedom to distribute copies of free software (and charge for
this service if you wish), that you receive source code or can get it
if you want it, that you can change the software or use pieces of it
in new free programs; and that you know you can do these things.
X
X To protect your rights, we need to make restrictions that forbid
anyone to deny you these rights or to ask you to surrender the rights.
These restrictions translate to certain responsibilities for you if
you distribute copies of the library, or if you modify it.
X
X For example, if you distribute copies of the library, whether gratis
or for a fee, you must give the recipients all the rights that we gave
you. You must make sure that they, too, receive or can get the source
code. If you link a program with the library, you must provide
complete object files to the recipients so that they can relink them
with the library, after making changes to the library and recompiling
it. And you must show them these terms so they know their rights.
X
X Our method of protecting your rights has two steps: (1) copyright
the library, and (2) offer you this license which gives you legal
permission to copy, distribute and/or modify the library.
X
X Also, for each distributor's protection, we want to make certain
that everyone understands that there is no warranty for this free
library. If the library is modified by someone else and passed on, we
want its recipients to know that what they have is not the original
version, so that any problems introduced by others will not reflect on
the original authors' reputations.
X
X Finally, any free program is threatened constantly by software
patents. We wish to avoid the danger that companies distributing free
software will individually obtain patent licenses, thus in effect
transforming the program into proprietary software. To prevent this,
we have made it clear that any patent must be licensed for everyone's
free use or not licensed at all.
X
X Most GNU software, including some libraries, is covered by the ordinary
GNU General Public License, which was designed for utility programs. This
license, the GNU Library General Public License, applies to certain
designated libraries. This license is quite different from the ordinary
one; be sure to read it in full, and don't assume that anything in it is
the same as in the ordinary license.
X
X The reason we have a separate public license for some libraries is that
they blur the distinction we usually make between modifying or adding to a
program and simply using it. Linking a program with a library, without
changing the library, is in some sense simply using the library, and is
analogous to running a utility program or application program. However, in
a textual and legal sense, the linked executable is a combined work, a
derivative of the original library, and the ordinary General Public License
treats it as such.
X
X Because of this blurred distinction, using the ordinary General
Public License for libraries did not effectively promote software
sharing, because most developers did not use the libraries. We
concluded that weaker conditions might promote sharing better.
X
X However, unrestricted linking of non-free programs would deprive the
users of those programs of all benefit from the free status of the
libraries themselves. This Library General Public License is intended to
permit developers of non-free programs to use free libraries, while
preserving your freedom as a user of such programs to change the free
libraries that are incorporated in them. (We have not seen how to achieve
this as regards changes in header files, but we have achieved it as regards
changes in the actual functions of the Library.) The hope is that this
will lead to faster development of free libraries.
X
X The precise terms and conditions for copying, distribution and
modification follow. Pay close attention to the difference between a
"work based on the library" and a "work that uses the library". The
former contains code derived from the library, while the latter only
works together with the library.
X
X Note that it is possible for a library to be covered by the ordinary
General Public License rather than by this special one.
X
X GNU LIBRARY GENERAL PUBLIC LICENSE
X TERMS AND CONDITIONS FOR COPYING, DISTRIBUTION AND MODIFICATION
X
X 0. This License Agreement applies to any software library which
contains a notice placed by the copyright holder or other authorized
party saying it may be distributed under the terms of this Library
General Public License (also called "this License"). Each licensee is
addressed as "you".
X
X A "library" means a collection of software functions and/or data
prepared so as to be conveniently linked with application programs
(which use some of those functions and data) to form executables.
X
X The "Library", below, refers to any such software library or work
which has been distributed under these terms. A "work based on the
Library" means either the Library or any derivative work under
copyright law: that is to say, a work containing the Library or a
portion of it, either verbatim or with modifications and/or translated
straightforwardly into another language. (Hereinafter, translation is
included without limitation in the term "modification".)
X
X "Source code" for a work means the preferred form of the work for
making modifications to it. For a library, complete source code means
all the source code for all modules it contains, plus any associated
interface definition files, plus the scripts used to control compilation
and installation of the library.
X
X Activities other than copying, distribution and modification are not
covered by this License; they are outside its scope. The act of
running a program using the Library is not restricted, and output from
such a program is covered only if its contents constitute a work based
on the Library (independent of the use of the Library in a tool for
writing it). Whether that is true depends on what the Library does
and what the program that uses the Library does.
X
X 1. You may copy and distribute verbatim copies of the Library's
complete source code as you receive it, in any medium, provided that
you conspicuously and appropriately publish on each copy an
appropriate copyright notice and disclaimer of warranty; keep intact
all the notices that refer to this License and to the absence of any
warranty; and distribute a copy of this License along with the
Library.
X
X You may charge a fee for the physical act of transferring a copy,
and you may at your option offer warranty protection in exchange for a
fee.
X
X 2. You may modify your copy or copies of the Library or any portion
of it, thus forming a work based on the Library, and copy and
distribute such modifications or work under the terms of Section 1
above, provided that you also meet all of these conditions:
X
X a) The modified work must itself be a software library.
X
X b) You must cause the files modified to carry prominent notices
X stating that you changed the files and the date of any change.
X
X c) You must cause the whole of the work to be licensed at no
X charge to all third parties under the terms of this License.
X
X d) If a facility in the modified Library refers to a function or a
X table of data to be supplied by an application program that uses
X the facility, other than as an argument passed when the facility
X is invoked, then you must make a good faith effort to ensure that,
X in the event an application does not supply such function or
X table, the facility still operates, and performs whatever part of
X its purpose remains meaningful.
X
X (For example, a function in a library to compute square roots has
X a purpose that is entirely well-defined independent of the
X application. Therefore, Subsection 2d requires that any
X application-supplied function or table used by this function must
X be optional: if the application does not supply it, the square
X root function must still compute square roots.)
X
These requirements apply to the modified work as a whole. If
identifiable sections of that work are not derived from the Library,
and can be reasonably considered independent and separate works in
themselves, then this License, and its terms, do not apply to those
sections when you distribute them as separate works. But when you
distribute the same sections as part of a whole which is a work based
on the Library, the distribution of the whole must be on the terms of
this License, whose permissions for other licensees extend to the
entire whole, and thus to each and every part regardless of who wrote
it.
X
Thus, it is not the intent of this section to claim rights or contest
your rights to work written entirely by you; rather, the intent is to
exercise the right to control the distribution of derivative or
collective works based on the Library.
X
In addition, mere aggregation of another work not based on the Library
with the Library (or with a work based on the Library) on a volume of
a storage or distribution medium does not bring the other work under
the scope of this License.
X
X 3. You may opt to apply the terms of the ordinary GNU General Public
License instead of this License to a given copy of the Library. To do
this, you must alter all the notices that refer to this License, so
that they refer to the ordinary GNU General Public License, version 2,
instead of to this License. (If a newer version than version 2 of the
ordinary GNU General Public License has appeared, then you can specify
that version instead if you wish.) Do not make any other change in
these notices.
X
X Once this change is made in a given copy, it is irreversible for
that copy, so the ordinary GNU General Public License applies to all
subsequent copies and derivative works made from that copy.
X
X This option is useful when you wish to copy part of the code of
the Library into a program that is not a library.
X
X 4. You may copy and distribute the Library (or a portion or
derivative of it, under Section 2) in object code or executable form
under the terms of Sections 1 and 2 above provided that you accompany
it with the complete corresponding machine-readable source code, which
must be distributed under the terms of Sections 1 and 2 above on a
medium customarily used for software interchange.
X
X If distribution of object code is made by offering access to copy
from a designated place, then offering equivalent access to copy the
source code from the same place satisfies the requirement to
distribute the source code, even though third parties are not
compelled to copy the source along with the object code.
X
X 5. A program that contains no derivative of any portion of the
Library, but is designed to work with the Library by being compiled or
linked with it, is called a "work that uses the Library". Such a
work, in isolation, is not a derivative work of the Library, and
therefore falls outside the scope of this License.
X
X However, linking a "work that uses the Library" with the Library
creates an executable that is a derivative of the Library (because it
contains portions of the Library), rather than a "work that uses the
library". The executable is therefore covered by this License.
Section 6 states terms for distribution of such executables.
X
X When a "work that uses the Library" uses material from a header file
that is part of the Library, the object code for the work may be a
derivative work of the Library even though the source code is not.
Whether this is true is especially significant if the work can be
linked without the Library, or if the work is itself a library. The
threshold for this to be true is not precisely defined by law.
X
X If such an object file uses only numerical parameters, data
structure layouts and accessors, and small macros and small inline
functions (ten lines or less in length), then the use of the object
file is unrestricted, regardless of whether it is legally a derivative
work. (Executables containing this object code plus portions of the
Library will still fall under Section 6.)
X
X Otherwise, if the work is a derivative of the Library, you may
distribute the object code for the work under the terms of Section 6.
Any executables containing that work also fall under Section 6,
whether or not they are linked directly with the Library itself.
X
X 6. As an exception to the Sections above, you may also compile or
link a "work that uses the Library" with the Library to produce a
work containing portions of the Library, and distribute that work
under terms of your choice, provided that the terms permit
modification of the work for the customer's own use and reverse
engineering for debugging such modifications.
X
X You must give prominent notice with each copy of the work that the
Library is used in it and that the Library and its use are covered by
this License. You must supply a copy of this License. If the work
during execution displays copyright notices, you must include the
copyright notice for the Library among them, as well as a reference
directing the user to the copy of this License. Also, you must do one
of these things:
X
X a) Accompany the work with the complete corresponding
X machine-readable source code for the Library including whatever
X changes were used in the work (which must be distributed under
X Sections 1 and 2 above); and, if the work is an executable linked
X with the Library, with the complete machine-readable "work that
X uses the Library", as object code and/or source code, so that the
X user can modify the Library and then relink to produce a modified
X executable containing the modified Library. (It is understood
X that the user who changes the contents of definitions files in the
X Library will not necessarily be able to recompile the application
X to use the modified definitions.)
X
X b) Accompany the work with a written offer, valid for at
X least three years, to give the same user the materials
X specified in Subsection 6a, above, for a charge no more
X than the cost of performing this distribution.
X
X c) If distribution of the work is made by offering access to copy
X from a designated place, offer equivalent access to copy the above
X specified materials from the same place.
X
X d) Verify that the user has already received a copy of these
X materials or that you have already sent this user a copy.
X
X For an executable, the required form of the "work that uses the
Library" must include any data and utility programs needed for
reproducing the executable from it. However, as a special exception,
the source code distributed need not include anything that is normally
distributed (in either source or binary form) with the major
components (compiler, kernel, and so on) of the operating system on
which the executable runs, unless that component itself accompanies
the executable.
X
X It may happen that this requirement contradicts the license
restrictions of other proprietary libraries that do not normally
accompany the operating system. Such a contradiction means you cannot
use both them and the Library together in an executable that you
distribute.
X
X 7. You may place library facilities that are a work based on the
Library side-by-side in a single library together with other library
facilities not covered by this License, and distribute such a combined
library, provided that the separate distribution of the work based on
the Library and of the other library facilities is otherwise
permitted, and provided that you do these two things:
X
X a) Accompany the combined library with a copy of the same work
X based on the Library, uncombined with any other library
X facilities. This must be distributed under the terms of the
X Sections above.
X
X b) Give prominent notice with the combined library of the fact
X that part of it is a work based on the Library, and explaining
X where to find the accompanying uncombined form of the same work.
X
X 8. You may not copy, modify, sublicense, link with, or distribute
the Library except as expressly provided under this License. Any
attempt otherwise to copy, modify, sublicense, link with, or
distribute the Library is void, and will automatically terminate your
rights under this License. However, parties who have received copies,
or rights, from you under this License will not have their licenses
terminated so long as such parties remain in full compliance.
X
X 9. You are not required to accept this License, since you have not
signed it. However, nothing else grants you permission to modify or
distribute the Library or its derivative works. These actions are
prohibited by law if you do not accept this License. Therefore, by
modifying or distributing the Library (or any work based on the
Library), you indicate your acceptance of this License to do so, and
all its terms and conditions for copying, distributing or modifying
the Library or works based on it.
X
X 10. Each time you redistribute the Library (or any work based on the
Library), the recipient automatically receives a license from the
original licensor to copy, distribute, link with or modify the Library
subject to these terms and conditions. You may not impose any further
restrictions on the recipients' exercise of the rights granted herein.
You are not responsible for enforcing compliance by third parties to
this License.
X
X 11. If, as a consequence of a court judgment or allegation of patent
infringement or for any other reason (not limited to patent issues),
conditions are imposed on you (whether by court order, agreement or
otherwise) that contradict the conditions of this License, they do not
excuse you from the conditions of this License. If you cannot
distribute so as to satisfy simultaneously your obligations under this
License and any other pertinent obligations, then as a consequence you
may not distribute the Library at all. For example, if a patent
license would not permit royalty-free redistribution of the Library by
all those who receive copies directly or indirectly through you, then
the only way you could satisfy both it and this License would be to
refrain entirely from distribution of the Library.
X
If any portion of this section is held invalid or unenforceable under any
particular circumstance, the balance of the section is intended to apply,
and the section as a whole is intended to apply in other circumstances.
X
It is not the purpose of this section to induce you to infringe any
patents or other property right claims or to contest validity of any
such claims; this section has the sole purpose of protecting the
integrity of the free software distribution system which is
implemented by public license practices. Many people have made
generous contributions to the wide range of software distributed
through that system in reliance on consistent application of that
system; it is up to the author/donor to decide if he or she is willing
to distribute software through any other system and a licensee cannot
impose that choice.
X
This section is intended to make thoroughly clear what is believed to
be a consequence of the rest of this License.
X
X 12. If the distribution and/or use of the Library is restricted in
certain countries either by patents or by copyrighted interfaces, the
original copyright holder who places the Library under this License may add
an explicit geographical distribution limitation excluding those countries,
so that distribution is permitted only in or among countries not thus
excluded. In such case, this License incorporates the limitation as if
written in the body of this License.
X
X 13. The Free Software Foundation may publish revised and/or new
versions of the Library General Public License from time to time.
Such new versions will be similar in spirit to the present version,
but may differ in detail to address new problems or concerns.
X
Each version is given a distinguishing version number. If the Library
specifies a version number of this License which applies to it and
"any later version", you have the option of following the terms and
conditions either of that version or of any later version published by
the Free Software Foundation. If the Library does not specify a
license version number, you may choose any version ever published by
the Free Software Foundation.
X
X 14. If you wish to incorporate parts of the Library into other free
programs whose distribution conditions are incompatible with these,
write to the author to ask for permission. For software which is
copyrighted by the Free Software Foundation, write to the Free
Software Foundation; we sometimes make exceptions for this. Our
decision will be guided by the two goals of preserving the free status
of all derivatives of our free software and of promoting the sharing
and reuse of software generally.
X
X NO WARRANTY
X
X 15. BECAUSE THE LIBRARY IS LICENSED FREE OF CHARGE, THERE IS NO
WARRANTY FOR THE LIBRARY, TO THE EXTENT PERMITTED BY APPLICABLE LAW.
EXCEPT WHEN OTHERWISE STATED IN WRITING THE COPYRIGHT HOLDERS AND/OR
OTHER PARTIES PROVIDE THE LIBRARY "AS IS" WITHOUT WARRANTY OF ANY
KIND, EITHER EXPRESSED OR IMPLIED, INCLUDING, BUT NOT LIMITED TO, THE
IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
PURPOSE. THE ENTIRE RISK AS TO THE QUALITY AND PERFORMANCE OF THE
LIBRARY IS WITH YOU. SHOULD THE LIBRARY PROVE DEFECTIVE, YOU ASSUME
THE COST OF ALL NECESSARY SERVICING, REPAIR OR CORRECTION.
X
X 16. IN NO EVENT UNLESS REQUIRED BY APPLICABLE LAW OR AGREED TO IN
WRITING WILL ANY COPYRIGHT HOLDER, OR ANY OTHER PARTY WHO MAY MODIFY
AND/OR REDISTRIBUTE THE LIBRARY AS PERMITTED ABOVE, BE LIABLE TO YOU
FOR DAMAGES, INCLUDING ANY GENERAL, SPECIAL, INCIDENTAL OR
CONSEQUENTIAL DAMAGES ARISING OUT OF THE USE OR INABILITY TO USE THE
LIBRARY (INCLUDING BUT NOT LIMITED TO LOSS OF DATA OR DATA BEING
RENDERED INACCURATE OR LOSSES SUSTAINED BY YOU OR THIRD PARTIES OR A
FAILURE OF THE LIBRARY TO OPERATE WITH ANY OTHER SOFTWARE), EVEN IF
SUCH HOLDER OR OTHER PARTY HAS BEEN ADVISED OF THE POSSIBILITY OF SUCH
DAMAGES.
X
X END OF TERMS AND CONDITIONS
X
X How to Apply These Terms to Your New Libraries
X
X If you develop a new library, and you want it to be of the greatest
possible use to the public, we recommend making it free software that
everyone can redistribute and change. You can do so by permitting
redistribution under these terms (or, alternatively, under the terms of the
ordinary General Public License).
X
X To apply these terms, attach the following notices to the library. It is
safest to attach them to the start of each source file to most effectively
convey the exclusion of warranty; and each file should have at least the
"copyright" line and a pointer to where the full notice is found.
X
X
X Copyright (C)
X
X This library is free software; you can redistribute it and/or
X modify it under the terms of the GNU Library General Public
X License as published by the Free Software Foundation; either
X version 2 of the License, or (at your option) any later version.
X
X This library is distributed in the hope that it will be useful,
X but WITHOUT ANY WARRANTY; without even the implied warranty of
X MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
X Library General Public License for more details.
X
X You should have received a copy of the GNU Library General Public
X License along with this library; if not, write to the Free
X Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
X
Also add information on how to contact you by electronic and paper mail.
X
You should also get your employer (if you work as a programmer) or your
school, if any, to sign a "copyright disclaimer" for the library, if
necessary. Here is a sample; alter the names:
X
X Yoyodyne, Inc., hereby disclaims all copyright interest in the
X library `Frob' (a library for tweaking knobs) written by James Random Hacker.
X
X , 1 April 1990
X Ty Coon, President of Vice
X
That's all there is to it!
SHAR_EOF
$shar_touch -am 0723151198 'QuadProg/COPYING.LIB' &&
chmod 0644 'QuadProg/COPYING.LIB' ||
$echo 'restore of' 'QuadProg/COPYING.LIB' 'failed'
if ( md5sum --help 2>&1 | grep 'sage: md5sum \[' ) >/dev/null 2>&1 \
&& ( md5sum --version 2>&1 | grep -v 'textutils 1.12' ) >/dev/null; then
md5sum -c << SHAR_EOF >/dev/null 2>&1 \
|| $echo 'QuadProg/COPYING.LIB:' 'MD5 check failed'
bbbfea63f32a8f2a1ecbea22ce897766 QuadProg/COPYING.LIB
SHAR_EOF
else
shar_count="`LC_ALL= LC_CTYPE= LANG= wc -c < 'QuadProg/COPYING.LIB'`"
test 25266 -eq "$shar_count" ||
$echo 'QuadProg/COPYING.LIB:' 'original size' '25266,' 'current size' "$shar_count!"
fi
fi
# ============= QuadProg/README ==============
if test -f 'QuadProg/README' && test "$first_param" != -c; then
$echo 'x -' SKIPPING 'QuadProg/README' '(file already exists)'
else
$echo 'x -' extracting 'QuadProg/README' '(text)'
sed 's/^X//' << 'SHAR_EOF' > 'QuadProg/README' &&
Quadratic Programming functions for S-plus
X
solve.QP(...) Solve a quadratic problem (no compact storage)
solve.QP.compact(...) Solve a quadratic problem with constraints
X stored in compact form
X
SHAR_EOF
$shar_touch -am 0723151198 'QuadProg/README' &&
chmod 0644 'QuadProg/README' ||
$echo 'restore of' 'QuadProg/README' 'failed'
if ( md5sum --help 2>&1 | grep 'sage: md5sum \[' ) >/dev/null 2>&1 \
&& ( md5sum --version 2>&1 | grep -v 'textutils 1.12' ) >/dev/null; then
md5sum -c << SHAR_EOF >/dev/null 2>&1 \
|| $echo 'QuadProg/README:' 'MD5 check failed'
311ae94c12606cc463c1b85a08e9991b QuadProg/README
SHAR_EOF
else
shar_count="`LC_ALL= LC_CTYPE= LANG= wc -c < 'QuadProg/README'`"
test 231 -eq "$shar_count" ||
$echo 'QuadProg/README:' 'original size' '231,' 'current size' "$shar_count!"
fi
fi
# ============= QuadProg/README.1st ==============
if test -f 'QuadProg/README.1st' && test "$first_param" != -c; then
$echo 'x -' SKIPPING 'QuadProg/README.1st' '(file already exists)'
else
$echo 'x -' extracting 'QuadProg/README.1st' '(text)'
sed 's/^X//' << 'SHAR_EOF' > 'QuadProg/README.1st' &&
This library contains routines and documentation for solving quadratic
programming problems.
X
I implemented these routines (using the algorithm of Goldfarb and
Idnani, 1982, 1983) and decided to distribute them in the hope that
they may also be useful to others. Hence, this software is provided
"as is" without any expressed or implied warranty. Specifcally:
X
X This library is free software; you can redistribute it and/or
X modify it under the terms of the GNU Library General Public
X License as published by the Free Software Foundation; either
X version 2 of the License, or (at your option) any later version.
X
X This library is distributed in the hope that it will be useful,
X but WITHOUT ANY WARRANTY; without even the implied warranty of
X MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
X Library General Public License for more details.
X
X You should have received a copy of the GNU Library General Public
X License along with this library; if not, write to the Free Software
X Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA
X 02111-1307 USA
X
The complete distribution consists of the following files:
X
solve.QP.compact.f fortran routine to solve a quadratic programming
X problem where the matrix defining the constraints
X is stored in a compact format.
solve.QP.f fortran routine to solve a quadratic programming
X problem.
aind.f achck.f utility functions for the S-interface of
X solve.QP.compact.f.
util.f some fortran utility functions related to
X functions from the LINPACK package.
solve.QP.compact.d S-PLUS help file for solve.QP.compact.
solve.QP.d S-PLUS help file for solve.QP.compact.
QuadProg.q the S-Plus functions of this library.
Makefile a makefile which should facilitate the
X installation of this library.
COPYING.LIB GNU Library General Public License (mentioned above).
README short description of the routines of the library
X (necessary to make `library(help=QuadProg)' work
X in S-PLUS once the library is installed).
README.1st you are just reading it ;-).
README.INSTALL some information on how to install the library.
README.pd README file which should go into the parent
X directory of the QuadProg directory in which the
X library resides in the end (necessary so that
X `library()' will give the desired information in
X S-PLUS).
ChangeLog history of the development of this library.
X
If you have any troubles or if you find any bugs in this library,
please let me know and I will try to fix it. I would also be
interested to hear from you if you find this library useful and for
what kind of problems you are using it.
X
Best wishes,
X
X Berwin
X
============================ Full address =============================
X
Berwin A Turlach Tel.: +61 (8) 8303 5418 (Secr.)
Department of Statistics +61 (8) 8303 3218 (direct)
The University of Adelaide FAX : +61 (8) 8303 3696
Adelaide SA 5005 e-mail: bturlach@stats.adelaide.edu.au
Australia http://www.stats.adelaide.edu.au/people/bturlach
X
=========================== Former address ============================
X
Berwin A. Turlach
Centre for Mathematics and its Applications Tel.: +61 (6) 249 0706 (Secr.)
Statistical Science Program +61 (6) 249 0731 (direct)
The Australian National University FAX : +61 (6) 249 5549
Canberra ACT 0200 e-mail: berwin.turlach@anu.edu.au
Australia
X
X
X
SHAR_EOF
$shar_touch -am 0723151198 'QuadProg/README.1st' &&
chmod 0644 'QuadProg/README.1st' ||
$echo 'restore of' 'QuadProg/README.1st' 'failed'
if ( md5sum --help 2>&1 | grep 'sage: md5sum \[' ) >/dev/null 2>&1 \
&& ( md5sum --version 2>&1 | grep -v 'textutils 1.12' ) >/dev/null; then
md5sum -c << SHAR_EOF >/dev/null 2>&1 \
|| $echo 'QuadProg/README.1st:' 'MD5 check failed'
8ec248fd0e539eaa552a18d60d437c97 QuadProg/README.1st
SHAR_EOF
else
shar_count="`LC_ALL= LC_CTYPE= LANG= wc -c < 'QuadProg/README.1st'`"
test 3833 -eq "$shar_count" ||
$echo 'QuadProg/README.1st:' 'original size' '3833,' 'current size' "$shar_count!"
fi
fi
# ============= QuadProg/README.INSTALL ==============
if test -f 'QuadProg/README.INSTALL' && test "$first_param" != -c; then
$echo 'x -' SKIPPING 'QuadProg/README.INSTALL' '(file already exists)'
else
$echo 'x -' extracting 'QuadProg/README.INSTALL' '(text)'
sed 's/^X//' << 'SHAR_EOF' > 'QuadProg/README.INSTALL' &&
This library was put together using the CHAPTER utility of S-PLUS.
The Makefile created by this utility was modified so that dynamic
loading (see WHICH_LOAD below) is the default. The hope is that this
should make the installation process fairly simple. The development
was done in a UNIX environment (also notable by the filenames) and the
instructions given here assume that a UNIX environment is being used.
If the installation does not work then please let me know.
Alternatively, you can use your own method of installation.
X
INSTALLING THE LIBRARY:
X
Place the file QuadProg.shar in the directory of your choice.
X
For explanatory purposes, we will suppose that the directory has the
path name:
X ~/splus/libraries
X
Starting from this directory type the commands:
X
1) sh QuadProg.shar (This creates the directory QuadProg underneath
X the current directory and unpacks all the required
X files)
2) cd QuadProg
X
3) Splus SHOME (This tells you the location of the S-PLUS home
X directory. Next, edit the Makefile and make the SHOME
X variable point to the location of this directory.
X Also, change the WHICH_LOAD variable to your preferred
X version. The description on using the library given
X below assumes that you choose one of the dynamic
X loading options. If you choose `dyn.load.shared', you
X may have to change in the function .First.lib [located
X in QuadProg.q] the line `dyn.load(X)' to
X `dyn.load.shared(X)' [can anybody confirm this?])
X
4) make install (This installs the S-PLUS functions and help files)
X
5) make load (This makes the required object code)
X
6) make clean (This deletes unnecessary files)
X
7) Does the parent directory (i.e. ~/splus/libraries) have a README
X file?
X yes: Append the last line of README.pd to that file (and delete
X README.pd)
X no: type the following command: mv README.pd ../README
X
X
USING THE LIBRARY:
X
Enter S-PLUS from any directory. Assume that the absolute path for
X~/splus/mylibrary is /home/userxyz/splus/library. To make the
QuadProg functions available to the current session type:
X
X assign("lib.loc","/home/userxyz/splus/mylibrary",where=0)
X library(QuadProg)
X
To save having to type the first command each time you want to use
QuadProg it is suggested that you create a .First function containing
these commands. It is then activated each time you enter S-PLUS.
X
SHAR_EOF
$shar_touch -am 0723151198 'QuadProg/README.INSTALL' &&
chmod 0644 'QuadProg/README.INSTALL' ||
$echo 'restore of' 'QuadProg/README.INSTALL' 'failed'
if ( md5sum --help 2>&1 | grep 'sage: md5sum \[' ) >/dev/null 2>&1 \
&& ( md5sum --version 2>&1 | grep -v 'textutils 1.12' ) >/dev/null; then
md5sum -c << SHAR_EOF >/dev/null 2>&1 \
|| $echo 'QuadProg/README.INSTALL:' 'MD5 check failed'
f71401f3c746308ae5745e1a4d5448c9 QuadProg/README.INSTALL
SHAR_EOF
else
shar_count="`LC_ALL= LC_CTYPE= LANG= wc -c < 'QuadProg/README.INSTALL'`"
test 2553 -eq "$shar_count" ||
$echo 'QuadProg/README.INSTALL:' 'original size' '2553,' 'current size' "$shar_count!"
fi
fi
# ============= QuadProg/README.pd ==============
if test -f 'QuadProg/README.pd' && test "$first_param" != -c; then
$echo 'x -' SKIPPING 'QuadProg/README.pd' '(file already exists)'
else
$echo 'x -' extracting 'QuadProg/README.pd' '(text)'
sed 's/^X//' << 'SHAR_EOF' > 'QuadProg/README.pd' &&
SECTION BRIEF DESCRIPTION
X
QuadProg Functions to solve Quadratic Programming Problems.
X
X
SHAR_EOF
$shar_touch -am 0723151198 'QuadProg/README.pd' &&
chmod 0644 'QuadProg/README.pd' ||
$echo 'restore of' 'QuadProg/README.pd' 'failed'
if ( md5sum --help 2>&1 | grep 'sage: md5sum \[' ) >/dev/null 2>&1 \
&& ( md5sum --version 2>&1 | grep -v 'textutils 1.12' ) >/dev/null; then
md5sum -c << SHAR_EOF >/dev/null 2>&1 \
|| $echo 'QuadProg/README.pd:' 'MD5 check failed'
478a6b4894eece8df0221b970ca8a2d8 QuadProg/README.pd
SHAR_EOF
else
shar_count="`LC_ALL= LC_CTYPE= LANG= wc -c < 'QuadProg/README.pd'`"
test 104 -eq "$shar_count" ||
$echo 'QuadProg/README.pd:' 'original size' '104,' 'current size' "$shar_count!"
fi
fi
# ============= QuadProg/ChangeLog ==============
if test -f 'QuadProg/ChangeLog' && test "$first_param" != -c; then
$echo 'x -' SKIPPING 'QuadProg/ChangeLog' '(file already exists)'
else
$echo 'x -' extracting 'QuadProg/ChangeLog' '(text)'
sed 's/^X//' << 'SHAR_EOF' > 'QuadProg/ChangeLog' &&
Thu Jul 23 05:13:14 1998 Berwin Turlach
X
X * Release 1.4 ready.
X
X * First.lib.template:
X Corrected copyright notice. Changed to Release 1.4.
X
X * solve.QP.d, solve.QP.compact.d:
X Corrected copyright notice. Changed occurrence of .AT to .AG (pointed
X out by Brian Ripley).
X
X * solve.QP.compact.sf, solve.QP.sf: Corrected copyright notice.
X
X * solve.QP.f, solve.QP.compact.f, aind.f, achck.f:
X Corrected copyright notice, rearranged declaration of variable (it
X seems f2c and the windows compiler do not like the old sequence;
X pointed out by Andreas Weingessel and Brian Ripley, respectively)
X
Mon Jul 6 06:25:56 1998 Berwin Turlach
X
X * Release 1.3 ready.
X
X * First.lib.template:
X Changed to GNU Library Public License, changed start up message and
X changed to Release 1.3.
X
X * solve.QP.compact.sf, solve.QP.sf:
X Changed to GNU Library Public License and added code to check at the
X end whether constraints are really satisfied.
X
X * achck.f: Initial revision
X
X * solve.QP.f, solve.QP.compact.f:
X Changed to GNU Library Public License and corrected a typo in the
X documentation pointed out by Peter Chen
X
X * solve.QP.compact.d, solve.QP.d:
X Added GNU Library Public License notice
X
X * aind.f: Changed to GNU Library Public License
X
Tue Dec 30 07:24:26 1997 Berwin Turlach
X
X * Release 1.2 ready.
X
X * First.lib.template:
X Changed to Release 1.2 after recent bug fix in Fortran code and
X changes to S-plus interface and help files.
X
X * solve.QP.compact.sf:
X Added a Fortran routine that checks whether all elements of Aind which
X are accessed are valid. Prior (S-plus coded) test was too crude and
X failed if constraints had different lengths and shorter constraints
X were filled up with zeros.
X
X * aind.f: Initial revision
X
Sun Dec 28 07:49:18 1997 Berwin Turlach
X
X * solve.QP.sf, solve.QP.compact.sf:
X Deleted `miter' from parameter list as underlying Fortran code doesn't
X have it anymore. Deleted all storage.mode commands and defined mode
X in call to Fortran routine (my preferred style these days). Added
X checks whether dimensions of arguments are compatible and sensible.
X
X * solve.QP.d, solve.QP.compact.d:
X Reworded the help file and brought it closer to the style described in
X the complements of V&R 2nd ed.
X
X * solve.QP.f, solve.QP.compact.f:
X Deleted `miter' from parameter list and program as it is not
X necessary, was there for debugging reasons but annoys now.
X
Wed Nov 12 08:52:22 1997 Berwin Turlach
X
X * solve.QP.f, solve.QP.compact.f:
X Fixed a bug, if an equality constraint is dropped and the sign of the
X normal vector changes, bvec() wasn't changed.
X
Tue Jul 1 07:41:19 1997 Berwin Turlach
X
X * solve.QP.compact.sf, solve.QP.sf:
X Changed error message from "is singular" to "is not positive definite"
X
Mon Jun 30 01:22:38 1997 Berwin Turlach
X
X * solve.QP.compact.f, solve.QP.f, solve.QP.sf, First.lib.template, solve.QP.compact.sf:
X Corrected address of Free Software Foundation
X
Mon Jun 23 09:09:28 1997 Berwin Turlach
X
X * Release 1.1 ready.
SHAR_EOF
$shar_touch -am 0723151198 'QuadProg/ChangeLog' &&
chmod 0644 'QuadProg/ChangeLog' ||
$echo 'restore of' 'QuadProg/ChangeLog' 'failed'
if ( md5sum --help 2>&1 | grep 'sage: md5sum \[' ) >/dev/null 2>&1 \
&& ( md5sum --version 2>&1 | grep -v 'textutils 1.12' ) >/dev/null; then
md5sum -c << SHAR_EOF >/dev/null 2>&1 \
|| $echo 'QuadProg/ChangeLog:' 'MD5 check failed'
3ea799e511a1db6855a535c97f89bb57 QuadProg/ChangeLog
SHAR_EOF
else
shar_count="`LC_ALL= LC_CTYPE= LANG= wc -c < 'QuadProg/ChangeLog'`"
test 3234 -eq "$shar_count" ||
$echo 'QuadProg/ChangeLog:' 'original size' '3234,' 'current size' "$shar_count!"
fi
fi
rm -fr _sh14571
exit 0