Revista de Matema´tica: Teor´ıa y Aplicaciones 2008 15(1) : 1–12
cimpa – ucr issn: 1409-2433
me´todos de punto interior para optimizacio´n
cuadra´tica convexa con matrices
no definidas positivas
Gonzalo Palencia F.∗ Rosina Hing C.† Mariledy Rojas C.‡
Denysde Medina S.§
Recibido/Received: 2 Feb 2007; Aceptado/Accepted: 31 Oct 2007
Resumen
En este art´ıculo se obtiene una modificacio´n del algoritmo recursivo de Cholesky
que permite la factorizacio´n de matrices semidefinidas positivas, au´n cuando e´stas
no sean definidas positivas, sin incrementar el costo computacional. Gracias a esta
factorizacio´n se transforman los Problemas de Programacio´n Cuadra´tica Convexa en
Problemas Co´nicos de Segundo Orden, los cuales se resuelven con la ayuda de la
generalizacio´n del algoritmo predictor-corrector de Menhrotra para dichos problemas.
Se realizan experimentos nume´ricos para validar los resultados.
Palabras clave: programacio´n cuadra´tica convexa, conos de segundo orden, me´todos de
punto interior.
Abstract
In this article a modification of the recursive algorithm of Cholesky is obtained
that allows the factorization of Semi Definite Positive Matrices, even though these
are not positive defined, without increasing the computational cost. Thanks to this
factorization Convex Quadratic Programming Problems are transformed into Second
Order Conical Problems, which are solved with the aid of the generalization of the
Predictor-Corrector algorithm of Mehrotra for these problems. There are carried out
numeric experiments for validating the results.
∗Departamento de Matema´tica, Facultad de Matema´tica, F´ısica y Computacio´n, Universidad Central
“Marta Abreu” de Las Villas, Carretera a Camajuan´ı 5.5 Kms, Santa Clara, Cuba. E-Mail:
†Misma direccio´n que G. Palencia. E-Mail: rhing@uclv.edu.cu
‡Departamento de Computacio´n, misma direccio´n que G. Palencia.
§Misma direccio´n que M. Rojas. E-Mail: denysde@uclv.edu.cu
1
2 G. Palencia – R. Hing – M. Rojas – D. Medina Rev.Mate.Teor.Aplic. (2008) 15(1)
Keywords: convex quadratic programming, second-order cones, interior point methods.
Mathematics Subject Classification: 90C52, 90C25.
1 Introduccio´n
En este trabajo consideramos el siguiente Problema de Programacio´n Cuadra´tica Convexa
(QCP):
1
2x
tQx+ ctx→ min
sujeto a: Ax = b
x ≥ 0
c, x ∈ Rn, Q ∈M(n, n), A ∈M(m,n), b ∈ Rn
rankA = m,Qt = Q,Q ≥ 0
(1)
Este es uno de los modelos ma´s utilizados para describir problemas relacionados con la
ingenier´ıa, el control, las finanzas, la optimizacio´n robusta y la optimizacio´n combinatoria.
Referencias a algunas de estas aplicaciones aparecen en [5] y [3].
Aunque se trata de un problema no lineal, el hecho de que tanto las restricciones como
las condiciones de Karush-Kuhn-Tucker se describen mediante sistemas lineales, permite
elaborar me´todos de solucio´n ma´s simples que los generales de la Programacio´n No-Lineal.
De hecho, en el caso en que Q sea definida positiva, estos problemas no son ma´s dif´ıciles
de resolver que los de Programacio´n Lineal.
En 1959 aparecio´ el primer algoritmo Simplex Cuadra´tico formulado por P. Wolfe, que
fue el ma´s ventajoso en esa e´poca, pero que resultaba realmente engorroso por el incre-
mento del nu´mero de variables que produc´ıa y por la necesidad de aplicar un procedimiento
de Fase I para determinar una solucio´n inicial factible.
En la actualidad para resolver el problema (1) se utilizan me´todos basados en la
condicio´n necesaria de primer orden que parten de un estimado inicial de la solucio´n
[13] y que requieren la resolucio´n del sistema de Karush-Kuhn-Tucker:[ −p
λ
]
=
[
Q AT
A 0
] [ −p
λ
]
=
[
g
0
]
, (2)
El paso de mayor complejidad computacional de este algoritmo es la solucio´n del sis-
tema (2). Si la matriz A es de rango completo y Q es definida positiva, es decir Z∗QZ > 0,
donde las columnas de Z constituyen una base del subespacio nulo de A, existen me´todos
para resolverlo que son eficientes para problemas pequen˜os o de mediana escala. Usual-
mente esto se hace factorizando directamente la matriz del sistema utilizando me´todos de
factorizacio´n de matrices indefinidas como el de Bunch-Parlett, invirtiendo la matriz Q y
luego factorizando ATkQ
−1Ak, o utilizando el me´todo del subespacio nulo que requiere el
ca´lculo de Z y la factorizacio´n de la matriz Z∗QZ [13]. Si Z∗QZ no es definida positiva, el
problema de determinar cuando un punto factible es minimizador global es un problema
NP-hard [9]
Aunque el algoritmo anterior ha sido descrito a partir de un punto inicial factible, el
mismo permite comenzar con un estimado inicial de la solucio´n del problema, que no sea
necesariamente factible [13].
optimizacio´n cuadra´tica convexa con matrices no definidas positivas 3
La importancia de buscar me´todos cada vez ma´s eficientes para resolver los QCP se
ha incrementado en las u´ltimas de´cadas debido a que uno de las te´cnicas ma´s usadas para
resolver problemas generales de programacio´n no lineal es la Programacio´n Cuadra´tica Se-
cuencial, que requiere resolver en cada iteracio´n un problema cuadra´tico. Un gran nu´mero
de paquetes computacionales como el NPSOL, NOPQL, OPSYC, OPTIMA, MATLAB y
SQP esta´n basados en este enfoque. En su forma ma´s pura este algoritmo reemplaza la
funcio´n objetivo en el punto actual, por su aproximacio´n cuadra´tica dada por:
q (xk) = ∇f (xk)T d+ 12d
T∇2xxL (xk, λk) d,
y reemplaza las restricciones por su aproximacio´n lineal. Por tanto para determinar el
paso desde el punto actual debe resolverse el siguiente problema cuadr´atico:
min q (xk)
s.a ci (xk) +∇ci (xk)T d = 0, i ∈Wk
Existen diferentes variantes de este algoritmo, por ejemplo, se sustituye la matriz
Hessiana del Lagrangiano por una aproximacio´n BFGS de e´ste, se utiliza bu´squeda lineal
para mejorar la convergencia o se utiliza el Lagrangiano aumentado con una funcio´n de
me´rito o una funcio´n de penalidad para resolver el compromiso entre el decrecimiento de
la funcio´n objetivo y la verificacio´n de las restricciones.
Debido al fabuloso desarrollo alcanzado por la Programacio´n Lineal (LP) a partir
de la aplicacio´n de algoritmos basados en los Me´todos del Punto Interior (IPM), que han
permitido resolver problemas con millones de variables y millones de restricciones en pocos
minutos [15] se comienzan a aplicar estas te´cnicas a los QP convexos. Efectivamente, si
escribimos las condiciones de KKT al problema (1), incluyendo las restricciones de igualdad
obtenemos que el minimizador debe satisfacer el sistema:
Qx+ c−ATλ− s = 0
Ax− b = 0
sixi = 0, i = 1, . . . , n
s.x ≥ 0
muy parecido al que se obtiene en la LP, por tanto podemos aplicar las te´cnicas de los
IPM y desarrollar algoritmos similares por ejemplo al Predictor-Corrector de Mehrotra.
La u´nica diferencia es que en este caso los incrementos de las variables primales y duales
no son independientes, ya que esta´n relacionados por la primera ecuacio´n. En este caso, al
igual que en los algoritmos que usan el conjunto activo, puede apreciarse que pueden surgir
dificultades computacionales si la matriz Q no es definida positiva. La implementacio´n es
muy eficiente para los problemas estrictamente convexos y aunque las iteraciones son ma´s
costosas que en los Me´todos del Conjunto Activo, se logra la convergencia con un nu´mero
menor de iteraciones.
El desarrollo de la teor´ıa de los IPM ha dominado el campo de la Optimizacio´n en los
u´ltimos 15 an˜os y debido fundamentalmente a la influencia de los trabajos de Nesterov y
Nemirovski [8] ha dado lugar al surgimiento de la Programacio´n Co´nica (CP) y en partic-
ular de la Programacio´n Co´nica de Segundo Orden (QCQP), donde se trata el problema
4 G. Palencia – R. Hing – M. Rojas – D. Medina Rev.Mate.Teor.Aplic. (2008) 15(1)
de la minimizacio´n de una funcio´n lineal sobre la interseccio´n de una variedad lineal af´ın,
con el producto cartesiano de conos de segundo orden. Estos problemas se caracterizan
por ser muy tratables computacionalmente ya que en [8] se prueba que los algoritmos de
IPM para estos problemas alcanzan una complejidad de
√
r donde r es el nu´mero de conos
de segundo orden que aparecen en el problema y Nemirovski y Scheinberg [10] prueban
que el Me´todo Primal Dual de Punto interior de la Programacio´n Lineal puede ser apli-
cado palabra por palabra a los QCQPs; el desarrollo de estos algoritmos comienza con los
trabajos de Nesterov y Todd [11] y [12].
Un problema QCQP se describe de la forma siguiente:
min c∗x
s.a. Ax = b
x ∈ K
(3)
donde K es el producto de conos convexos Ki cerrados de segundo orden,
Ki =
{(
x0i , xˆi
) | ‖xˆi‖ ≤ x0i , x0i ≥ 0, xˆi ∈ Rni} , i = 1, . . . , N . En particular el octante no
negativo puede ser considerado como un producto de conos de segundo orden unidimen-
sionales [14].
En los trabajos de Alizadeh y Goldfarb [2] demuestran que los QCP pueden ser expre-
sados como QCQP, de la siguiente forma:
minu0
s.a Q1/2x− uˆ = 12Q−1/2c
Ax = b
x ≥ 0, (u0, uˆ) ∈ K
K = {(u0, uˆ) | u0 ≥ 0, ‖uˆ‖ ≤ u0} .
(4)
La obtencio´n de esta transformacio´n desde el punto de vista pra´ctico implica el ca´lculo
de la matriz Q1/2 = λ1/2T T , donde λ1/2 = diag
(√
λ1, . . . ,
√
λn
)
y T es una matriz ortog-
onal cuyas columnas son los vectores propios de la matriz Q; lo que a su vez requiere el
ca´lculo de los valores propios, ca´lculo de vectores propios y ortonormalizacio´n de e´stos.
Tambie´n se requiere invertir la matriz Q1/2, lo cual no es posible si la matriz Q no es
definida positiva.
Dentro de este contexto el objetivo de nuestro trabajo es buscar un algoritmo eficiente
para transformar el QCP a un QCQP, que sea eficiente au´n cuando la matriz Q no sea
definida positiva y aplicar para su solucio´n un algoritmo Primal-Dual de Punto Interior.
2 Factorizacio´n de Cholesky de matrices semidefinidas po-
sitivas
El algoritmo que proponemos esta´ basado en la aplicacio´n de la teor´ıa de las matrices
semidefinidas, por lo que para facilitar la exposicio´n aclareremos algunas notaciones que
sera´n utilizadas.
Notacio´n 1: Menores de orden p de la matriz Q
optimizacio´n cuadra´tica convexa con matrices no definidas positivas 5
Q
(
i1, . . . , ip
k1, . . . , kp
)
,
1 ≤ i1 < i2 < . . . < ip ≤ n
1 ≤ j1 < j2 < . . . < jp ≤ n. , p ≤ n.
donde los ı´ndices i representan las filas y los ı´ndices j las columnas, cuya interseccio´n
determina el menor.
Notacio´n 2: Menores principales de orden p de la matriz Q.
Q
(
i1, . . . , ip
i1, . . . , ip
)
, 1 ≤ i1 < i2 < . . . < ip ≤ n, p ≤ n.
Notacio´n 3: Menores principales diagonales de orden p de la matriz Q.
Q
(
1, . . . , p
1, . . . , p
)
, p ≤ n.
Un argumento importante que sera´ utilizado en nuestro trabajo es la caracterizacio´n
de las matrices semidefinidas positivas, expresadas en el teorema siguiente [6].
Teorema 1 Una condicio´n necesaria y suficiente para que la matriz Q sea semidefinida
positiva es que todos sus menores principales sean no negativos, es decir:
Q
(
i1, . . . , ip
i1, . . . , ip
)
≥ 0, 1 ≤ i1 < i2 < . . . < ip ≤ n, p = 1, . . . , n. (5)
Nuestro primer paso sera´ llevar la forma cuadra´tica contenida en la funcio´n objetivo a
su forma cano´nica sin utilizar las matrices que aparecen en (4). En su lugar utilizaremos
la factorizac´ıon de Cholesky de la matriz Q.
Para las matrices definidas positivas se usa el me´todo de factorizacio´n de Cholesky,
que descompone la matriz Q en el producto siguiente:
Q = LL∗,
donde L es una matriz triangular inferior. El algoritmo para calcular los elementos de L
no requiere pivotacio´n, requiere solamente n
3
3 operaciones y en su forma recursiva, para
matrices Q definidas positivas [16], utiliza las fo´rmulas siguientes:
Q =
[
q11 Q
∗
21
Q21 Q22
]
=
[
l11 0
L21 L22
] [
l11 L
∗
21
0 L∗22
]
Q =
[
l211 l11L
∗
21
l11L21 L21L
∗
21 + L22L
∗
22.
]
.
Para determinar la primera columna de L se usan las igualdades:
l11 =
√
q11, l11L21 = Q21. (6)
De las fo´rmulas anteriores es fa´cil ver que L22L∗22, que es tambie´n definida positiva [16],
es la factorizacio´n de la matriz ya conocida Q22 − L21L∗21. Por tanto, usando fo´rmulas
6 G. Palencia – R. Hing – M. Rojas – D. Medina Rev.Mate.Teor.Aplic. (2008) 15(1)
similares a la anteriores se halla la primera columna de L22 y consecuentemente la segunda
columna de L. Repitiendo recursivamente el proceso, al cabo de n iteraciones se obtiene
la matriz L buscada.
Es evidente que si tratamos de aplicar el algoritmo anterior a una matriz que sea semi
definida positiva pero no definida positiva, falla en algu´n paso en el que lˆ11 = 0 por lo que
es necesario hacer algunas modificaciones a e´ste.
Se han propuesto algunos algoritmos que primeramente determinan el rango de la
matriz y seleccionan una submatriz DP a la cual se le aplica el algoritmo anterior, pero
estos ca´lculos son considerados inestables [7]. Se proponen tambie´n algoritmos basados en
el conocimiento de la estructura del subespacio nulo de la matriz Q para una computacio´n
muy precisa de la factorizacio´n de Cholesky, muy eficientes para el caso en que tanto la
matriz Q, como la base del subespacio nulo conocida, sean suficientemente poco densas [4].
En nuestro trabajo proponemos un algoritmo que no requiere ninguna informacio´n previa
sobre la matriz semidefinida positiva, para lograrlo demostraremos la siguiente propiedad.
Proposicio´n 1 Si Q es una matriz semidefinida positiva y qii = 0 entonces qij = 0,
j = 1, . . . , n.
Demostracio´n: Sea qii = 0, para todo j = 1, . . . , n se tiene que los menores de segundo
orden que se obtienen de la interseccio´n de las filas i y j con las columnas i y j son iguales
a: qiiqjj − q2ij = −q2ij, y aplicando (5) se obtiene la igualdad a 0 de todos los elementos de
la fila y la columna j-e´simas.
Con este resultado podemos analizar las igualdades (6) del algoritmo recursivo de
Cholesky y obtener que si q11 = 0 entonces l11 = 0 pero como tambie´n la columna
Q21 = Q12 = 0, la segunda igualdad se cumple para cualquier valor de L21, por lo que la
descomposicio´n no es u´nica. Para nuestro algoritmo tomamos la que simplifica ma´s los
ca´lculos, esto es L21 = 0, lo que implica que L22L∗22 = Q22. Esto significa que en nuestro
algoritmo modificado no se incrementan los ca´lculos, por el contrario pueden disminuir
cuando disminuye el rango Q y el nu´mero de operaciones en el peorde los casos sigue
siendo n
3
3 .
Por ejemplo apliquemos la modificacio´n propuesta a la matriz 1 0 10 0 0
1 0 3
 , con valores propios : √2 + 2, 2−√2, 0.
La matriz es semidefinida positiva ya que sus valores propios son no negativos, pero
no es definida positiva ya que uno de sus valores propios es nulo. En el primer paso del
algoritmo de Cholesky se obtiene:.
l11 =
√
1 = 1, l11L∗21 = 1
[
l21 l31
]
=
[
0 1
]
,
L22L
∗
22 =
[
l222 l22l32
l22l32 l
2
32 + l
2
33
]
=
[
0 0
0 2
]
= Qˆ
En el segundo paso debemos aplicar la modificacio´n propuesta:
l22 =
√
0 = 0, l22L32 = 0l32 = 0
optimizacio´n cuadra´tica convexa con matrices no definidas positivas 7
La solucio´n es indeterminada, si hacemos l32 = 0, entonces:[
l22 l32
]
=
[
0 0
]
, l233 = qˆ33 = 2, l32 =
√
2
Por tanto la factorizacio´n obtenida sera´: 1 0 00 0 0
1 0
√
2
 1 0 10 0 0
0 0
√
2
 =
 1 0 10 0 0
1 0 3

Para comprobar que la descomposicio´n no es u´nica podemos hacer l32 = 1, entonces se
obtiene que: [
l22 l32
]
=
[
0 1
]
, l233 = qˆ33 − l232 = 2− 1 = 1, l33 = 1
Por tanto la factorizacio´n obtenida sera´: 1 0 00 0 0
1 1 1
 1 0 10 0 1
0 0 1
 =
 1 0 10 0 0
1 0 3

Para un segundo ejemplo trabajamos con la matriz factorizada en [4], utilizando la
informacio´n sobre su subespacio nulo. Con la modificacio´n propuesta al algoritmo recursivo
de Cholesky, sin ninguna informacio´n adicional, obtuvimos el mismo resultado:
1 0 1 1 3
0 9 3 9 9
1 3 3 6 8
1 9 6 14 16
3 9 8 16 22
 =

1 0 0 0 0
0 3 0 0 0
1 1 1 0 0
1 3 2 0 0
3 3 2 0 0
 ·

1 0 1 1 3
0 3 1 3 3
0 0 1 2 2
0 0 0 0 0
0 0 0 0 0

Hacemos notar que la factorizacio´n de Cholesky de esta matriz no fue hallada con los
paquetes de ca´lculo disponibles, por no ser definida positiva.
3 Transformacio´n del QCP a un QCQP
Como resultado de aplicar nuestro algoritmo modificado de Cholesky a la matriz Q de
rango r, se obtiene la matriz triangular inferior L =
[
l1 l2 . . . . ln
]
, donde exacta-
mente n − r vectores columnas son nulos. Si r < n esta matriz es singular, por lo que
construimos una nueva matriz de la forma siguiente:
Sea I = {i | li 6= 0} y J = {j | lj = 0} , entonces hacemos
L˜ =
[
l˜1 l˜2 . . . l˜n
]
, l˜k =
{
lk, k ∈ I
ek k ∈ J
donde ek es el k-e´simo vector cano´nico.
8 G. Palencia – R. Hing – M. Rojas – D. Medina Rev.Mate.Teor.Aplic. (2008) 15(1)
Utilizamos ahora la transformacio´n lineal ξ = L˜∗x, es decir ξk = l˜∗kx, k ∈ I y
ξk = xk, k ∈ J , con la que la forma cuadra´tica definida por Q se transforma en:
x∗Qx = x∗LL∗x = (L∗x)∗ L∗x =
∑
k∈I
ξ2k.
Usando la matriz de permutacio´n:
P ∗ =
[
ei1 . . . eirr ej1 . . . ejn−r
]
donde los ı´ndices i pertenecen a I y los j a J y la transformacio´n z = Pξ = PL˜∗x, logramos
que las primeras r componentes de z sean iguales a las componentes de ξ cuyos ı´ndices se
encuentran en I.
Denotando por
(
PL˜∗
)−1
= H,y G = [hik]i=1,...r;k=1...n el problema (1) se puede ree-
scribir como: ∑r
i=1 z
2
i + c
∗Hz → min
s.a AHz = b
Gz ≥ 0
zj ≥ 0, j = r + 1, . . . , n.
Introduciendo variables superplus y las notaciones: ζ, c˜ ∈ Rn+r el problema se transforma
en: ∑r
i=1 ζ
2
i + c˜
∗ζ → min
s.a A˜ζ = b˜
ζj ≥ 0, j = r + 1, n+ r.
(7)
Se introduce ahora una nueva variable unidimensional y, y se considera a ζ1 como un
vector r-dimensional que contiene las primeras r componentes de ζ, el problema anterior
puede reescribirse como:
y1 − y2 → min
s.a. A˜ζ = b˜∑r
i=1 ζ
2
i + c˜
∗ζ ≤ y1 − y2
ζj ≥ 0, j = r + 1, n+ r
yi ≥ 0, i = 1, 2.
(8)
Se ha tenido en cuenta que la parte lineal de la restriccio´n cuadra´tica puede ser nega-
tiva, por lo que se acota con la diferencia de dos variables no negativas.
Para transformar la restriccio´n cuadra´tica de (8) en un cono de segundo orden usare-
mos una te´cnica similar a la utilizada en [3] para transformar una restriccio´n cuadra´tica
convexa en la interseccio´n de dos restricciones lineales afines y un cono de segundo orden.
Primeramente hacemos y1−y2− c˜∗ζ = θ y an˜adimos la restriccio´n lineal t = 1. Con lo que
la restriccio´n cuadra´tica se transforma en
∑r
i=1 ζ
2
i ≤ θt, y1 − y2 − c˜∗ζ = θ, y t = 1.
optimizacio´n cuadra´tica convexa con matrices no definidas positivas 9
Haciendo ahora θ = u+ v y t = u− v, se obtiene finalmente el QCQP siguiente:
y1 − y2 → min
s.a. A˜ζ = b˜
u+ v − y1 + y2 + c˜∗ζ = 0
u− v = 1∑r
i=1 ζ
2
i + v
2 ≤ u2, u ≥ 0
ζj ≥ 0, j = r + 1, n+ r, y1, y2 ≥ 0.
(9)
Es evidente que (9) es un problema de QCQP, ya que la funcio´n objetivo es lineal y las
filas 2a., 3a. y 4a. contienen las restricciones lineales. La penu´ltima restriccio´n representa
un cono de segundo orden de dimensio´n r+ 2 y la u´ltima el producto cartesiano de n+ 2
conos de segundo orden unidimensionales.
El ca´lculo de mayor complejidad que debemos hacer es la inversio´n de la matriz L˜∗,
que por ser triangular so´lo requiere un nu´mero de operaciones que es equivalente tambie´n
a n
3
3 , por lo que incluyendo la factorizacio´n de Q la transformacio´n requiere un nu´mero
de operaciones que es aproximadamente 2n
3
3 .
4 Resultados experimentales
Los QCQP son analizados detalladamente en [1], donde se formulan condiciones de de-
generacio´n y de complementariedad estricta y adema´s se elabora para ellos un algoritmo
primal-dual de punto interior, eficiente y nume´ricamente estable bajo condiciones de no
degeneracio´n de los problemas primal y dual y condiciones de complementariedad estricta.
En este algoritmo se implementa una generalizacio´n del me´todo precdictor-corrector de
Mehrotra, adaptado para ser aplicado a CP definidos con conos cuadra´ticos de segundo
(QCQP) orden o conos de matrices semidefinidas positivas (SDP).
En nuestro trabajo se elabora e implementa un algoritmo para QCP con dos fases,
en la primera se transforma el problema en un QCQP utilizando el me´todo propuesto
en el ep´ıgrafe anterior y en la segunda se resuelve el problema transformado usando el
algoritmo propuesto en [1]. Este algoritmo comienza con un punto inicial no factible por
lo que utilizamos la misma heur´ıstica propuesta en [1] para generar este punto inicial dentro
del cono, pero no factible. Se usan los vectores rp y rd para representar las violaciones
de la factibilidad en el primal y el dual respectivamente. El algoritmo termina cuando el
error = gap+ ‖rp‖+ ‖rd‖ ≤ 10−7, o cuando ‖x‖ o ‖z‖ exceden a un valor definido por un
para´metro, en este u´ltimo caso decimos que el primal o el dual no tienen solucio´n acotada.
Otras condiciones de parada esta´n dadas por valores muy pequen˜os de los taman˜os del paso
en el primal o el dual y por un nu´mero ma´ximo de iteraciones. Finalmente se multiplica
por la matriz H el vector de las primeras n componentes de la solucio´n del QCQP y se
tiene la solucio´n del QCP. Este co´digo fue denominado Larx.
Para validar las dos fases de nuestro algoritmo, comparamos los resultados que se
obtienen al aplicarlo a 14 QCP de diferentes dimensiones, generados aleatoriamente, con
los resultados que se obtienen al resolver estos mismos problemas con las funciones de
optimizacio´n del Mathematica 5.5. Los resultados aparecen en la tabla 1, donde m es
10 G. Palencia – R. Hing – M. Rojas – D. Medina Rev.Mate.Teor.Aplic. (2008) 15(1)
el nu´mero de restricciones y n es el nu´mero de variables. Todos los problemas tienen
matrices semidefinidas positivas, en la tabla DP significa que dicha matriz es adema´s
definida positiva y NDP que no lo es. Para cada uno de los co´digos aparece el valor de la
funcio´n objetivo obtenido y el tiempo de ejecucio´n expresado en segundos.
Problema Larx Mathematica
m n Tipo Objetivo Seg. Objetivo Seg.
3 7 NDP 31.394 0.2 31.394 1.375
NDP 36.963 0.2 36.953 0.547
NDP 31.969 0.0 31.969 1.187
8 13 NDP 56.034 0.4 56.034 4.782
NDP 225.680 0.8 230.598 4.422
4 8 DP 37.718 0.4 37.716 1.266
DP 51.166 0.2 51.166 1.625
12 17 DP 538.030 2.7 541.314 4.907
537.420 1.4 538.189 4.953
424.130 1.9 424.133 4.187
10 20 DP 345.110 3.2 345.037 13.641
DP 729.280 4.0 729.270 18.938
15 25 NDP 730.840 11 730.838 28.703
NDP 927.200 12 927.193 34.86
NDP 792.880 11 792.880 23.688
Tabla 1: Resultados comparativos entre el nuevo me´todo Larx y Mathematica 5.5 sobre
14 QCP generados aleatoriamente (DP: matriz definida positiva, NDP: matriz no definida
positiva).
Como se puede apreciar en la tabla 1, los valores de las soluciones obtenidas son bas-
tante similares, hay 5 pequen˜as diferencias que favorecen al Mathematica y 4 diferencias,
no todas tan pequen˜as, que favorecen al Larx. Casi todas debidas a errores de aproxi-
macio´n al calcular el valor final de la funcio´n objetivo, pues los vectores soluciones son
iguales.
Los resultados obtenidos al resolver con el Larx otros 45 QCP generados aleatoriamente
se resumen en la tabla 2.
Las tres primeras columnas tienen el mismo significado que en la tabla 1 y caracterizan
a los problemas cuyos resultados aparecen en cada fila. N es el nu´mero de instancias del
problema.analizado en la fila. En R.E. aparecen los lugares decimales de la primera cifra
no nula en el mayor y el menor error obtenido en los problemas de la fila. M.E. es el
error promedio, M.I. es el nu´mero promedio de iteraciones y MT el tiempo promedio de
ejecucio´n de los problemas de la fila. Las columnas ME y MI reflejan la calidad de las
soluciones y la eficiencia del Larx..
optimizacio´n cuadra´tica convexa con matrices no definidas positivas 11
m n Tipo N R.E. M.E. M.I. M.T.
3 7 NDP 3 E-6,E-9 7.5E-7 27 0.03
4 8 DP 3 E-7,E-12 2.7E-7 73.6 0.2
8 9 DP 4 E-7,E-6 1.7E-7 28.5 0.225
8 12 DP 3 E-8,E-13 3.3E-8 132.3 0.67
8 13 NDP 3 E-10,E-12 2.3E-10 142 1.43
9 13 DP 5 E-7,E-10 8.4E-8 186.6 1.08
9 14 DP 5 E-6,E-11 6,5E-7 156,6 1.14
9 14 NDP 3 E-9,E-11 2.9E-9 151,7 1.06
10 15 NDP 5 E-7,E-12 4.3E-8 318.2 9.9
10 20 DP 5 E-8,E-11 5.1E-9 246.6 5.46
12 19 NDP 5 E-7,E-12 1.2E-7 253 4.38
15 25 NDP 1 E-5, 4.9E-5 312 33
Tabla 2: Resultados obtenidos con Larx sobre 45 QCP generados aleatoriamente.
5 Conclusiones
En este trabajo se obtiene una modificacio´n del algoritmo recursivo de Cholesky, que
permite factorizar matrices semidefinidas positivas sin incrementar el costo computacional
en las no definidas positivas, ya que no requiere informacio´n previa sobre el subespacio
nulo de e´stas. Con ayuda de esta factorizacio´n se logra transformar cualquier problema de
Programacio´n Cuadr´atica Convexa en un Problema de Programacio´n Co´nica de Segundo
Orden, au´n cuando la matriz de la forma cuadra´tica no sea definida positiva. Lo anterior
permite aplicar una generalizacio´n del algoritmo predictor-corrector de Mehrotra de Punto
Interior, para resolver el Problema Cuadra´tico Convexo. Estos resultados fueron validados
por las pruebas nume´ricas realizadas.
Aunque los experimentos fueron realizados con problemas a pequen˜a escala, los resul-
tados obtenidos nos permiten conjeturar que este procedimiento de transformacio´n nos
ofrece una alternativa para poder aplicar los eficientes Me´todos de Punto Interior a la
solucio´n de Problemas Cuadra´ticos Convexos a gran escala, incluyendo los casos en que
estos problemas no sean estrictamente convexos.
Referencias
[1] Alizadeh, F.; Schmieta, S.H. (1997) “Optimization with semidefinite, quadratic and
linear constraints”, Rutcor Research Report, RRR 23–97, November.
[2] Alizadeh, F.; Goldfarb, D. (2001) “Second-Order Cone Programming, Rutcor Re-
search Report.
[3] Andersen, E.D.; Ross, C.; Terlaky, T. (2002) “On implementing a primal-dual
interior-point method for conic quadratic optimization”, Springer Verlag.
12 G. Palencia – R. Hing – M. Rojas – D. Medina Rev.Mate.Teor.Aplic. (2008) 15(1)
[4] Arbenz, P.; Dramac, Z. (2000) “On semipositive definite matrices whith null space
known”, ETH Zu¨rich, Computer Science Departament, Technical Report Research
#352, November.
[5] Ben-Tal, A.; Nemirovski, A. (2001) Lectures on Modern Convex Optimization: Anal-
ysis, Algorithms and Engineering Applications. MPS/Siam, Series on Optimization,
SIAM.
[6] Gantmacher, F.R. (1959) The Theory of Matrices. Chelsea, New York.
[7] Higham, N. (1996) “The symmetric indefinite pactorization: stability an applications
in optimization”.
[8] Lustig, J.; Marsten, R.E,; Shanno, D.F.(1994) “Interior point methods for linear
programming: computational state of the art”, ORSA J. on Comptutation.
[9] Murty, Kabadi (1987) “Some NP-complete problems in quadratic and nonlinear pro-
gramming”, Mathematical Programming 39.
[10] Neemirovski, A.; Cherinberg, K. (1996) “Extension of Karmarkar algoritm onto con-
vex quadratically constrained quadratic programming”, Mathematical Programming
72.
[11] Nesterov, Y.E.; Todd, M.J. (1997) “Self-scaled barriers and interior-point methods
for convex programming”, Math. of Oper. Res. 22(1): 1–42.
[12] Nesterov, Y.E.; Todd, M.J. (1998) “Primal-dual interior-point methods for self-scaled
cones”, SIAM J. Optim. 8: 324–364.
[13] Nocedal, J.; Wright, S.J. (1999) Numerical Optimization. Springer Series in Opera-
tions Research, Springer Verlag, New York.
[14] Renegar, J. (2001) A Mathematical View of Interior-Point. Methods in Convex Op-
timization, MPS-SIAM Series en Optimization.
[15] Ross, C.; Terlaky, T.; Vial, J. (1997) Theory and Algoritms for Linear Optimization
an Interior Point Approach. John Wiley and Sons, New York.
[16] Vandenberghe, L. (2002) “The Cholesky factorization”, EE103 Winter, Germany.