Search All of the Math Forum:

Views expressed in these public forums are not endorsed by NCTM or The Math Forum.

Notice: We are no longer accepting new posts, but the forums will continue to be readable.

Replies: 4   Last Post: Aug 12, 2014 12:27 PM

 Messages: [ Previous | Next ]
 Todd McCollough Posts: 150 Registered: 7/17/09
Posted: Mar 16, 2011 8:18 PM

I have been trying to implement the nonlinear conjugate gradient algorithm and have been having some difficulty and am not sure what is going wrong.

On wikipedia there is a discussion of the conjugate gradient algorithm and a numerical example. http://en.wikipedia.org/wiki/Conjugate_gradient_method

This algorithm works in MATLAB as I modified the code provided by Xin-She Yang while using the the Polak&#8211;Ribière formula

A=[4,1;1,3];
b=[1;2];

% Solve the system by the conjugate gradient method
%u=rand(n,1); % initial guess u0
u=[2;1];
r=b-A*u; % intial r0
d=r; % intial d0=r0
delta=10^(-5); % accuracy

while max(abs(r))>delta
alpha=r'*r/(d'*A*d); % alpha_n
K=r'*r; % temporary value
u=u+alpha*d; % u_{n+1}
last_r=r;
r=r-alpha*A*d; % r_{n+1}
%beta=r'*r/K; % beta_n Fletcher Reeves
beta=(r'*(r-last_r))/K; % beta_n Polak Ribiere
d=r+beta*d; % d_{n+1}
u' % disp u on the screen
end

The nonlinear conjugate gradient method is what I am interested in since it allows a generic function to be minimized instead of Ax=b. Discussion on Wikipedia is here http://en.wikipedia.org/wiki/Nonlinear_conjugate_gradient

There is also discussion here on the modifications from going from linear to nonlinear conjugate gradient. http://www.ee.ucla.edu/~vandenbe/236C/lectures/cg.pdf

My attempt at the nonlinear conjugate gradient method is

A=[4,1;1,3];
b=[1;2];

x=[2;1];

F=A*x-b;

d=-g;

for k=1:10

a=linesrch(d,x,b,A);

x = x + a*d;

last_g=g;
K=g'*g;

F=A*x-b;

beta=norm((g'*(g-last_g)))^2/norm(K)^2; % beta_n Polak Ribiere

d=-g+beta*d;

x' % disp x on the screen

end

---

function a = linesrch(d,x,b,A)

a = 1e-4;
while myfun(a,d,x,b,A) < 0
a = 1.6*a;
end

a = fminunc(@(a) myfun(a,d,x,b,A),1);

----

function [f] = myfun(a,p,x,b,A)

Right now the results blow up. I am rying to replicate the results of the Ax=b example using nonlinear conjugate gradient so I can move on to other functions.

I see some recent related discussion here as well http://www.mathworks.com/matlabcentral/newsreader/view_thread/304151 which I also tried. Any suggestions would be welcome. Thanks.

Date Subject Author
3/16/11 Todd McCollough