# [SciPy-user] Usage of fmin_ncg

Nils Wagner nwagner@iam.uni-stuttgart...
Thu Mar 8 10:06:38 CST 2007

```Hi all,

I have a question wrt optimize.fmin_ncg(F,v_0,Fp, fhess_p = Fhess_p).
especially the last argument is of interest. I am not sure if I applied
it correctly.

The help function yields

fhess_p -- a function to compute the Hessian of f times an
arbitrary vector: fhess_p (x, p, *args)

Is p the arbitrary vector ?
Any hint will be appreciated.

Nils

P.S.: The Hessian is not directly available.

from scipy import *
from pylab import plot, show, semilogy, xlabel, ylabel, legend, savefig,
title, figure, ylim
#
#
def R(v):
""" Rayleigh quotient """
return dot(v.T,A*v)/dot(v.T,B*v)

def F(v):
""" Objective function """
Bv = bsolve(v)
res = linalg.norm(A*Bv-R(Bv)*B*Bv)/linalg.norm(B*Bv)
h =    0.25*dot(bsolve(v).T,v)**2+0.5*dot(sigma_solve(v).T,v)
data.append(res)
return h

def Fp(v):
""" Gradient of the objective function """
return dot(bsolve(v).T,v)*bsolve(v)+sigma_solve(v)

def Fhess_p(v,p):
""" Hessian applied to an arbitrary vector p """
return
dot(bsolve(v).T,v)*bsolve(p)+sigma_solve(p)+2*dot(outer(bsolve(v),bsolve(v)),p)

#
# Test matrices from MatrixMarket
#
n = A.shape[0]

sigma = 47.5                                       # Preassigned shift
bsolve = linsolve.splu(B).solve                    # inv(B)*( )
sigma_solve = linsolve.splu(A - sigma*B).solve     # inv(A-sigma*B)*( )

random.seed(10)
v_0 = random.rand(n)                               # Initial vector

data=[]
v = optimize.fmin_ncg(F,v_0,Fp, fhess_p = Fhess_p)
semilogy(arange(0,len(data)),data)
print  'Newton CG',-1./(2*sqrt(-F(v)))+sigma
print R(bsolve(v))
print

```