Optimization: principles and algorithms, by Michel Bierlaire
Go to the documentation of this file.
1 %> \file
2 %> Algorithm 9.2: Conjugate gradient algorithm for quadratic problems. Implementation of algorithm 9.2 of \cite Bier15-book
3 %>
4 %> @author <a href="http://people.epfl.ch/michel.bierlaire">Michel Bierlaire</a>
5 %> @date Sat Apr 5 23:32:26 2014
6 %> @ingroup Algorithms
7 %> @ingroup chap09
8
9 %> @note Tested with \ref run0908cg.m
10 %> @note Called by \ref newtonLocalQuadratic
11
12 %> Applies the conjugate gradient method to solve \f[\min_x \frac{1}{2} x^T Q x + b^T x\f] where \f$Q \in \mathbb{R}^n\times\mathbb{R}^n \f$ and \f$b \in \mathbb{R}^n\f$.
13 %> @param Q matrix of size \f$n \times n \f$.
14 %> @param b vector of size \f$n\f$.
15 %> @param x0 starting point
16 %> @param printlevel if different from 0, informations are printed at each iteration (Default: 0)
17 %> @return [D, solution]
18 %> @return D: matrix gathering all directions generated by the algorithms
19 %> @return solution: optimal solution
20 function [D, solution] = conjugateGradient(Q,b,x0,printlevel=0)
21  n = size(x0,1) ;
22  xk = x0 ;
23  gk = Q * xk + b ;
24  dk = -gk ;
25  D = dk ;
26  betak = 0 ;
27  for k=1:n
28  xprint = xk ;
29  gprint = gk ;
30  dprint = dk ;
31  if (printlevel != 0)
32  printf("%3d\t%+10.5e\t%+10.5e\t%+10.5e",k,xprint(1),gprint(1),dprint(1)) ;
33  endif
34  denom = dk' * Q * dk ;
35  if (denom <= 0)
36  error("The matrix must de positive definite") ;
37  endif
38  alphak = - (dk' * gk) / denom;
39  xk = xk + alphak * dk ;
40  gkp1 = Q * xk + b ;
41  if (printlevel != 0)
42  if (betak == 0)
43  printf("\t%+10.5e\n",alphak) ;
44  else
45  printf("\t%+10.5e\t%+10.5e\n",alphak,betak) ;
46  endif
47  for i=2:n
48  printf("\t%+10.5e\t%+10.5e\t%+10.5e\t\t\n",xprint(i),gprint(i),dprint(i)) ;
49  endfor
50  endif
51  betak = (gkp1' * gkp1) / (gk' * gk) ;
52  dk = -gkp1 + betak * dk ;
53  D = [D dk] ;
54  gk = gkp1 ;
55 % printf("%d %f\n",k,norm(gk)) ;
56  endfor
57  if (printlevel != 0)
58  printf("%3d\t%+10.5e\t%+10.5e\t\t\t\n",n+1,xk(i),gk(i)) ;
59  for i=2:n
60  printf("\t%+10.5e\t%+10.5e\t\t\t\n",xk(i),gk(i)) ;
61  endfor
62  endif
63  solution = xk ;
64 endfunction
function newtonLocalQuadratic(in obj, in x0, in eps, in cg, in maxiter)
Applies local Newton algorithm to solve where is the gradient of the objective function.
function conjugateGradient(in Q, in b, in x0, in printlevel)
Applies the conjugate gradient method to solve where and .