impls/pounders/gqt.c

90 c     positive number delta, this subroutine determines a vector
97 c           norm(x) <= delta.
102 c            norm(x) <= (1+rtol)*delta,
106 c            abs(norm(x) - delta) <= rtol*delta.
115 c       subroutine gqt(n,a,lda,b,delta,rtol,atol,itmax,
137 c       delta is a double precision variable.
138 c         On entry delta is a bound on the Euclidean norm of x.
139 c         On exit delta is unchanged.
153 c              norm(x) <= (1 + rtol)*delta
165 c            multiplier for the constraint norm(x) <= delta.
214 PetscErrorCode gqt(PetscInt n, PetscReal *a, PetscInt lda, PetscReal *b, PetscReal delta, PetscReal rtol, PetscReal atol, PetscInt itmax, PetscReal *retpar, PetscReal *retf, PetscReal *x, PetscInt *retinfo, PetscInt *retits, PetscReal *z, PetscReal *wa1, PetscReal *wa2)
216   PetscReal    f = 0.0, p001 = 0.001, p5 = 0.5, minusone = -1, delta2 = delta * delta;
262   parl = PetscMax(bnorm / delta - parl, pars);
264   paru = PetscMax(0.0, bnorm / delta + paru);
311       if (PetscAbs(xnorm - delta) <= rtol * delta || (par == 0 && xnorm <= (1.0 + rtol) * delta)) info = 1;
318       /* Compute a negative curvature solution of the form x + alpha*z,  where norm(x+alpha*z)==delta */
321       if (xnorm < delta) {
323         PetscCallBLAS("BLASdot", prod = BLASdot_(&blasn, z, &blas1, x, &blas1) / delta);
324         temp  = (delta - xnorm) * ((delta + xnorm) / delta);
325         alpha = temp / (PetscAbs(prod) + PetscSqrtScalar(prod * prod + temp / delta));
350         parc = (xnorm - delta) / (delta * temp * temp);
354       if (xnorm > delta) {
356       } else if (xnorm < delta) {
413         f = -p5 * (rxnorm * rxnorm + par * delta * delta - rznorm * rznorm);