# Path: nlpy/tools/norms.py
# def norm1(x):
# def norm2(x):
# def normp(x,p):
# def norm_infty(x):
# def normest(A, tol=1.0e-6, maxits=100):
#         A = np.random.randn(n,m)
#         A = np.random.rand(n,n)
#         A = .5*(A.T + A)


from nlpy.model.amplpy import AmplModel
from nlpy.krylov.linop import SimpleLinearOperator
from nlpy.krylov.ppcg import ProjectedCG as Ppcg
from nlpy.tools import norms
import numpy
import sys

if len(sys.argv) < 2:
        print 'This demo reads the nl file of an AMPL-decoded model and'
        print 'demonstrates the projected conjugate gradient by solving the'
        print 'following quadratic program with trust-region constraint:'
        print '     min   <g, x> + 1/2 <x, H x>'
        print '     s.t.  J d = 0,  |d| <= radius,'
        print
        print ' where g = grad(x0), H = Hess(x0), J = jac(x0), radius = |g|/10'
        print 'Please supply a constrained problem as input'
        sys.exit(0)

problem = sys.argv[1]
nlp = AmplModel(problem)

if nlp.m == 0:
    sys.stderr.write( 'Problem should have general constraints\n' )
    nlp.close()
    sys.exit(1)

x0 = numpy.ones(nlp.n) #nlp.x0            # Initial point
g  = numpy.random.random(nlp.n) #nlp.grad(x0)      # Objective gradient
J  = nlp.jac(x0)       # Obtain constraints Jacobian

# Path: nlpy/tools/norms.py
# def norm1(x):
# def norm2(x):
# def normp(x,p):
# def norm_infty(x):
# def normest(A, tol=1.0e-6, maxits=100):
#         A = np.random.randn(n,m)
#         A = np.random.rand(n,n)
#         A = .5*(A.T + A)


# Path: examples/demo_ppcg.py
from nlpy.model.amplpy import AmplModel
from nlpy.krylov.linop import SimpleLinearOperator
from nlpy.krylov.ppcg import ProjectedCG as Ppcg
from nlpy.tools import norms
import numpy
import sys


if len(sys.argv) < 2:
        print 'This demo reads the nl file of an AMPL-decoded model and'
        print 'demonstrates the projected conjugate gradient by solving the'
        print 'following quadratic program with trust-region constraint:'
        print '     min   <g, x> + 1/2 <x, H x>'
        print '     s.t.  J d = 0,  |d| <= radius,'
        print
        print ' where g = grad(x0), H = Hess(x0), J = jac(x0), radius = |g|/10'
        print 'Please supply a constrained problem as input'
        sys.exit(0)

problem = sys.argv[1]
nlp = AmplModel(problem)

if nlp.m == 0:
    sys.stderr.write( 'Problem should have general constraints\n' )
    nlp.close()
    sys.exit(1)

x0 = numpy.ones(nlp.n) #nlp.x0            # Initial point
g  = numpy.random.random(nlp.n) #nlp.grad(x0)      # Objective gradient
J  = nlp.jac(x0)       # Obtain constraints Jacobian

delta = 1000 #max(1.0, 100 * norms.norm_infty(g))