# [SciPy-user] [Fwd: Re: Computing eigenvalues by trace minimization]

dmitrey openopt@ukr....
Fri Jun 29 04:40:12 CDT 2007

1) you don't need to use .T in dot, it's executing automatically:
f = lambda x: dot(x.T,dot(A,x)) # => dot(x,dot(A,x))
h = lambda x: dot(x.T,dot(B,x))-1.0 # => dot(x,dot(B,x))-1.0
2) I have no symeig module, I will try to install it now
3) I didn't decide yet what to do if stop creteria funtol or xtol report
true, but constraints are not satisfied yet (despite alg can
successfully continue to decriase them, but I can't know will he do the
trick or no and how many iterations it will require).
So now you should decide, what contol is ok for your problem, and then
set appropriate contol, funtol, xtol
note that funtol and xtol are stop criteria, not desired tolerance of
solution.
So, if you need contol=1e-6, as it is set by default, you need to reduce
funtol and xtol (also, maybe, gradtol, if it will stop your problem
before desired contol is achieved) from default values (1e-6) to
something like 1e-8.
4) please wait some time before I will implement some changes about df,
dc and dh. I will inform you about a hour later (as well as about my
symeig module installation results) to update svn.
Use only func values for now.
5) Some weeks later (I hope) there will be excellent native OO QPQC
solver (for positive-defined matrices only).

HTH, D.

Nils Wagner wrote:
>  Dmitrey,
>
> You asked me to move my inquiry to scipy-user.
> How can I improve my script wrt. the results ?
> How can I provide the gradients ?
> Is it possible to extend the code to rectangular matrices  instead of
> vectors ?
>
> python -i test_qpqc.py
> starting solver lincher (BSD  license)  with problem  unnamed
> itn 0: Fk= 0.285714285714 maxResidual= 0
> itn 10 : Fk= 0.0820810434319 maxResidual= 0.00260656014331
> N= 1.0 alpha= 0.978713763748
> itn 20 : Fk= 0.0810194693773 maxResidual= 1.71023808193e-05
> N= 1.0 alpha= 0.132923656348
> itn 22 : Fk= 0.0810154618808 maxResidual= 3.51005361554e-06
> N= 1.0 alpha= 0.354395906532
> solver lincher finished solving the problem unnamed
> istop:  4 (|| F[k] - F[k-1] || < funtol)
> Solver:   Time elapsed = 0.24   CPU Time Elapsed = 0.24
> NO FEASIBLE SOLUTION is obtained (max residual = 3.51005361554e-06)
>
> x_opt: [ 0.11973864  0.22980377  0.32143896  0.38732479  0.42178464
> 0.4223829
>   0.38829845  0.32320301  0.2311393   0.12059048]
> f_opt: [ 0.08101546]
>
>
> Smallest eigenvalue by symeig 0.081014052771
>
> Nils
>
>
> ------------------------------------------------------------------------
>
> from scikits.openopt import NLP
> from scipy import diag, ones, identity, rand, trace, dot, shape, zeros
> from scipy.linalg import solve, norm
> from symeig import symeig
>
> #f = lambda x: trace(dot(x.T,dot(A,x))) # Version for m > 1
> #h = lambda x: dot(x.T,dot(B,x))-identity(m)
> f = lambda x: dot(x.T,dot(A,x))
> h = lambda x: dot(x.T,dot(B,x))-1.0
>
> n = 10 # order of A and B
> m = 1  # subspace dimension # m > 1 doesn't work
>
> A = diag(2*ones(n))-diag(ones(n-1),-1)-diag(ones(n-1),1)
> B = identity(n)
>
> #x0 = rand(n,m)          # Initial subspace
> b = zeros(n)
> b[-1] = 1.
> x0 = solve(A,b)
> x0 = x0/norm(x0)
> p = NLP(f, x0, h=h)
>
> #p.df = lambda x: 2*dot(A,x)
> #p.dh = lambda x: 2*dot(B,x)
>
> r = p.solve('lincher')
> print
> print 'x_opt:',r.xf
> print 'f_opt:',r.ff
> print
>
> w,v = symeig(A,B)
> print
> print 'Smallest eigenvalue by symeig',w[0]
>
>
> ------------------------------------------------------------------------
>
> Subject:
> Re: Computing eigenvalues by trace minimization
> From:
> dmitrey <openopt@ukr.net>
> Date:
> Thu, 28 Jun 2007 17:48:21 +0300
> To:
> Nils Wagner <nwagner@iam.uni-stuttgart.de>
>
> To:
> Nils Wagner <nwagner@iam.uni-stuttgart.de>
>
>
> Hi Nils,
> constraints.
> Currently you can yuse only lincher to solve the problem.
> f = lambda x: trace (X^T A X)
>
> h = lambda x X^T B X - I_p
>
> p = NLP(f, x0, h=h)
> r = p.solve('lincher')
>
> p.df = ...
> p.dh = ...
>
> However, lincher can fail to solve the problem if something like
> ill-conditioned matrix will be encountered.
> Our department's chiefman Petro I. Stetsyuk said me he can provide an
> excellent QPQC algorithm for my Python QPQC solver (based on ralg),
> but currently he is busy with other problems. Maybe, it will be ready
> some weeks or months later.
> HTH, D.
>
>
>
> Nils Wagner wrote:
>> Hi Dmitrey,
>>
>> This inquiry is not related to my previous problem.
>> Can I use openopt to solve the quadratic minimization problem with
>> openopt
>>
>> minimize trace (X^T A X)
>>
>> subjected to the constraints
>>
>> X^T B X = I_p,
>>
>> where I_p denotes the identity matrix of order p. A is symmetric
>> positive semidefinite
>> and B is symmetric positive definite. A, B \in \mathds{R}^{n \times n}.
>> A small example would be appreciated.
>>
>>
>>
>> Nils
>>
>>
>> Reference: A. Sameh and J. Wisniewski
>> A trace minimization algorithm for the generalized eigenvalue problem
>> SIAM J. Numer. Anal. Vol. 19 (1982) pp. 1243-1259
>>
>>
>>
>>
>
> ------------------------------------------------------------------------
>
> _______________________________________________
> SciPy-user mailing list
> SciPy-user@scipy.org
> http://projects.scipy.org/mailman/listinfo/scipy-user
>