Mercurial > octave-nkf

comparison scripts/sparse/qmr.m @ 19867:9fc020886ae9

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
maint: Clean up m-files to follow Octave coding conventions. Try to trim long lines to < 80 chars. Use '##' for single line comments. Use '(...)' around tests for if/elseif/switch/while. Abut cell indexing operator '{' next to variable. Abut array indexing operator '(' next to variable. Use space between negation operator '!' and following expression. Use two newlines between endfunction and start of %!test or %!demo code. Remove unnecessary parens grouping between short-circuit operators. Remove stray extra spaces (typos) between variables and assignment operators. Remove stray extra spaces from ends of lines.
author Rik <rik@octave.org>
date Mon, 23 Feb 2015 14:54:39 -0800
parents 5f2c0ca0ef51
children 0a3ca546d7fc
comparison
equal deleted inserted replaced
19866:a1acca0c2216 19867:9fc020886ae9
184 v = vt / rho0; 184 v = vt / rho0;
185 y = y / rho0; 185 y = y / rho0;
186 w = wt / xi1; 186 w = wt / xi1;
187 z = z / xi1; 187 z = z / xi1;
188 188
189 delta1 = z' * y; ## If delta1 == 0, method fails. 189 delta1 = z' * y; # If delta1 == 0, method fails.
190 190
191 yt = M2m1x (y); 191 yt = M2m1x (y);
192 zt = M1tm1x (z); 192 zt = M1tm1x (z);
193 193
194 if (iter == 1) 194 if (iter == 1)
198 p = yt - (xi1*delta1/eps0) * p; 198 p = yt - (xi1*delta1/eps0) * p;
199 q = zt - (rho0*delta1/eps0) * q; 199 q = zt - (rho0*delta1/eps0) * q;
200 endif 200 endif
201 pt = Ax (p); 201 pt = Ax (p);
202 202
203 eps0 = q' * pt; ## If eps0 == 0, method fails. 203 eps0 = q' * pt; # If eps0 == 0, method fails.
204 beta1 = eps0 / delta1; ## If beta1 == 0, method fails. 204 beta1 = eps0 / delta1; # If beta1 == 0, method fails.
205 vt = pt - beta1 * v; 205 vt = pt - beta1 * v;
206 206
207 y = M1m1x (vt); 207 y = M1m1x (vt);
208 rho1 = norm(y); 208 rho1 = norm(y);
209 wt = Atx (q) - beta1 * w; 209 wt = Atx (q) - beta1 * w;
210 z = M2tm1x (wt); 210 z = M2tm1x (wt);
211 211
212 xi1 = norm(z); 212 xi1 = norm(z);
213 theta1 = rho1 / (gamma0 * abs(beta1)); 213 theta1 = rho1 / (gamma0 * abs(beta1));
214 gamma1 = 1 / sqrt(1 + theta1^2); ## If gamma1 == 0, method fails. 214 gamma1 = 1 / sqrt(1 + theta1^2); # If gamma1 == 0, method fails.
215 eta1 = -eta0 * rho0 * gamma1^2 / (beta1 * gamma0^2); 215 eta1 = -eta0 * rho0 * gamma1^2 / (beta1 * gamma0^2);
216 216
217 if (iter == 1) 217 if (iter == 1)
218 d = eta1 * p; 218 d = eta1 * p;
219 s = eta1 * pt; 219 s = eta1 * pt;