/[escript]/trunk/paso/src/Smoother.cpp
ViewVC logotype

Contents of /trunk/paso/src/Smoother.cpp

Parent Directory Parent Directory | Revision Log Revision Log


Revision 4803 - (show annotations)
Wed Mar 26 06:52:28 2014 UTC (5 years, 10 months ago) by caltinay
File size: 18041 byte(s)
Removed obsolete wrappers for malloc and friends.
Paso_Pattern -> paso::Pattern

1
2 /*****************************************************************************
3 *
4 * Copyright (c) 2003-2014 by University of Queensland
5 * http://www.uq.edu.au
6 *
7 * Primary Business: Queensland, Australia
8 * Licensed under the Open Software License version 3.0
9 * http://www.opensource.org/licenses/osl-3.0.php
10 *
11 * Development until 2012 by Earth Systems Science Computational Center (ESSCC)
12 * Development 2012-2013 by School of Earth Sciences
13 * Development from 2014 by Centre for Geoscience Computing (GeoComp)
14 *
15 *****************************************************************************/
16
17
18 /************************************************************************************/
19
20 /* Paso: Gauss-Seidel */
21
22 /************************************************************************************/
23
24 /* Author: artak@uq.edu.au */
25
26 /************************************************************************************/
27
28 #include "Paso.h"
29 #include "Preconditioner.h"
30 #include "PasoUtil.h"
31 #include "BlockOps.h"
32
33 #include <stdio.h>
34
35
36 /************************************************************************************/
37
38 /* free all memory used by Smoother */
39
40 void Paso_Preconditioner_Smoother_free(Paso_Preconditioner_Smoother * in) {
41 if (in!=NULL) {
42 Paso_Preconditioner_LocalSmoother_free(in->localSmoother);
43 delete in;
44 }
45 }
46 void Paso_Preconditioner_LocalSmoother_free(Paso_Preconditioner_LocalSmoother * in) {
47 if (in!=NULL) {
48 delete[] in->diag;
49 delete[] in->pivot;
50 delete[] in->buffer;
51 delete in;
52 }
53 }
54 /************************************************************************************/
55
56 /* constructs the symmetric Gauss-Seidel preconditioner */
57
58 Paso_Preconditioner_Smoother* Paso_Preconditioner_Smoother_alloc(Paso_SystemMatrix * A_p, const bool jacobi, const bool is_local, const bool verbose)
59 {
60
61 /* allocations: */
62 Paso_Preconditioner_Smoother* out=new Paso_Preconditioner_Smoother;
63 if (! Esys_checkPtr(out)) {
64 out->localSmoother=Paso_Preconditioner_LocalSmoother_alloc(A_p->mainBlock,jacobi,verbose);
65 out->is_local=is_local;
66 }
67 if (Esys_MPIInfo_noError(A_p->mpi_info)) {
68 return out;
69 } else {
70 Paso_Preconditioner_Smoother_free(out);
71 return NULL;
72 }
73 }
74 Paso_Preconditioner_LocalSmoother* Paso_Preconditioner_LocalSmoother_alloc(paso::SparseMatrix * A_p, const bool jacobi, bool verbose)
75 {
76
77 const dim_t n=A_p->numRows;
78 const dim_t n_block=A_p->row_block_size;
79 const dim_t block_size=A_p->block_size;
80
81 double time0=Esys_timer();
82 /* allocations: */
83 Paso_Preconditioner_LocalSmoother* out=new Paso_Preconditioner_LocalSmoother;
84 if (! Esys_checkPtr(out)) {
85
86 out->diag=new double[((size_t) n) * ((size_t) block_size)];
87 out->pivot=new index_t[ ((size_t) n) * ((size_t) n_block)];
88 out->buffer=new double[((size_t) n) * ((size_t) n_block)];
89 out->Jacobi=jacobi;
90
91 if ( ! ( Esys_checkPtr(out->diag) || Esys_checkPtr(out->pivot) ) ) {
92 paso::SparseMatrix_invMain(A_p, out->diag, out->pivot );
93 }
94
95 }
96 time0=Esys_timer()-time0;
97
98 if (Esys_noError()) {
99 return out;
100 } else {
101 Paso_Preconditioner_LocalSmoother_free(out);
102 return NULL;
103 }
104 }
105
106 /*
107
108 performs a few sweeps of the form
109
110 S (x_{k} - x_{k-1}) = b - A x_{k-1}
111
112 where x_{0}=0 and S provides some approximation of A.
113
114 Under MPI S is built using A_p->mainBlock only.
115 If Smoother is local the defect b - A x_{k-1} is calculated using A_p->mainBlock only.
116
117 The MPI-versioned smoother works in the following way:
118 (1) initialize x
119 (2) calculate residual r_n = b - A * x_n
120 (3) store residual r_n in b_new and pass to Paso_Preconditioner_LocalSmoother_Sweep() as input
121 (4) /delta x_{n+1} is returned (stored in b_new) as the output of Paso_Preconditioner_LocalSmoother_Sweep()
122 (5) recover x_{n+1} by /delta x_{n+1} + x_n
123 (6) repeat steps 2-5 for until sweeps number reduce to 0
124 */
125
126 void Paso_Preconditioner_Smoother_solve(Paso_SystemMatrix* A_p, Paso_Preconditioner_Smoother * smoother, double * x, const double * b,
127 const dim_t sweeps, const bool x_is_initial)
128 {
129 const dim_t n= (A_p->mainBlock->numRows) * (A_p->mainBlock->row_block_size);
130
131 double *b_new = smoother->localSmoother->buffer;
132 dim_t nsweeps=sweeps;
133 if (smoother->is_local) {
134 Paso_Preconditioner_LocalSmoother_solve(A_p->mainBlock,smoother->localSmoother,x,b,sweeps,x_is_initial);
135 } else {
136 if (! x_is_initial) {
137 Paso_Copy(n, x, b);
138
139 Paso_Preconditioner_LocalSmoother_Sweep(A_p->mainBlock,smoother->localSmoother,x);
140 nsweeps--;
141 }
142 while (nsweeps > 0 ) {
143 Paso_Copy(n, b_new, b);
144 Paso_SystemMatrix_MatrixVector_CSR_OFFSET0((-1.), A_p, x, 1., b_new); /* b_new = b - A*x */
145 Paso_Preconditioner_LocalSmoother_Sweep(A_p->mainBlock,smoother->localSmoother,b_new);
146 Paso_AXPY(n, x, 1., b_new);
147 nsweeps--;
148 }
149
150 }
151 }
152
153 err_t Paso_Preconditioner_Smoother_solve_byTolerance(Paso_SystemMatrix* A_p, Paso_Preconditioner_Smoother * smoother,
154 double * x, const double * b,
155 const double atol, dim_t *sweeps, const bool x_is_initial)
156 {
157 const dim_t n= (A_p->mainBlock->numRows) * (A_p->mainBlock->row_block_size);
158 double *b_new = smoother->localSmoother->buffer;
159 const dim_t max_sweeps=*sweeps;
160 dim_t s=0;
161 double norm_dx = atol * 2.;
162 err_t errorCode = PRECONDITIONER_NO_ERROR;
163
164 if (! x_is_initial) {
165 Paso_Copy(n, x, b);
166 Paso_Preconditioner_LocalSmoother_Sweep(A_p->mainBlock,smoother->localSmoother,x);
167 norm_dx=Paso_lsup(n,x,A_p->mpi_info);
168 s++;
169 }
170 while (norm_dx > atol) {
171 Paso_Copy(n, b_new, b);
172 Paso_SystemMatrix_MatrixVector((-1.), A_p, x, 1., b_new); /* b_new = b - A*x */
173 Paso_Preconditioner_LocalSmoother_Sweep(A_p->mainBlock,smoother->localSmoother,b_new);
174 norm_dx=Paso_lsup(n,b_new,A_p->mpi_info);
175 Paso_AXPY(n, x, 1., b_new);
176 if (s >= max_sweeps) {
177 errorCode = PRECONDITIONER_MAXITER_REACHED;
178 break;
179 }
180 s++;
181 }
182 *sweeps=s;
183 return errorCode;
184 }
185
186 void Paso_Preconditioner_LocalSmoother_solve(paso::SparseMatrix* A_p, Paso_Preconditioner_LocalSmoother * smoother, double * x, const double * b,
187 const dim_t sweeps, const bool x_is_initial)
188 {
189 const dim_t n= (A_p->numRows) * (A_p->row_block_size);
190 double *b_new = smoother->buffer;
191 dim_t nsweeps=sweeps;
192
193 if (! x_is_initial) {
194 Paso_Copy(n, x, b);
195 Paso_Preconditioner_LocalSmoother_Sweep(A_p,smoother,x);
196 nsweeps--;
197 }
198
199 while (nsweeps > 0 ) {
200 Paso_Copy(n, b_new, b);
201
202 paso::SparseMatrix_MatrixVector_CSR_OFFSET0((-1.), A_p, x, 1., b_new); /* b_new = b - A*x */
203 Paso_Preconditioner_LocalSmoother_Sweep(A_p,smoother,b_new);
204 Paso_AXPY(n, x, 1., b_new);
205 nsweeps--;
206 }
207 }
208
209 /*
210 Gauss-Seidel sweep to calculate /delta x_{n+1} from matrix A and
211 residual r_n. It has two steps: forward substitution and backward
212 substitution.
213 Forward substitution: (D+L) * /delta x_{n+1/2} = r_n
214 Backward substitution: (D+U) * /delta x_{n+1} = r_n - L * /delta x_{n+1/2}
215 = D * /delta x_{n+1/2}
216 where A = D + L + U
217 /delta x_{n+1/2} = x_{n+1/2} - x_n
218 /delta x_{n+1} = x_{n+1} - x_n
219
220 Input: Matrix A
221 residual r_n (in x)
222 Output: /delta x_{n+1} (in x)
223 */
224 void Paso_Preconditioner_LocalSmoother_Sweep(paso::SparseMatrix* A, Paso_Preconditioner_LocalSmoother * smoother, double * x)
225 {
226 const dim_t nt=omp_get_max_threads();
227 if (smoother->Jacobi) {
228 Paso_BlockOps_solveAll(A->row_block_size,A->numRows,smoother->diag,smoother->pivot,x);
229 } else {
230 if (nt < 2) {
231 Paso_Preconditioner_LocalSmoother_Sweep_sequential(A,smoother,x);
232 } else {
233 Paso_Preconditioner_LocalSmoother_Sweep_colored(A,smoother,x);
234 }
235 }
236 }
237
238 /* inplace Gauss-Seidel sweep in sequential mode: */
239
240 void Paso_Preconditioner_LocalSmoother_Sweep_sequential(paso::SparseMatrix* A_p, Paso_Preconditioner_LocalSmoother * smoother, double * x)
241 {
242 const dim_t n=A_p->numRows;
243 const dim_t n_block=A_p->row_block_size;
244 double *diag = smoother->diag;
245 index_t* pivot = smoother->pivot;
246 const dim_t block_len=A_p->block_size;
247
248
249 register dim_t i,k;
250 register index_t iptr_ik, mm;
251 register double rtmp;
252 int failed = 0;
253
254 const index_t* ptr_main = paso::SparseMatrix_borrowMainDiagonalPointer(A_p);
255
256 (void)pivot; /* silence warning from var being unused by macros */
257 (void)block_len;
258 /* forward substitution */
259
260 if (n_block==1) {
261 x[0]*=diag[0];
262 for (i = 1; i < n; ++i) {
263 mm=ptr_main[i];
264 /* x_i=x_i-a_ik*x_k (with k<i) */
265 rtmp=x[i];
266 for (iptr_ik=A_p->pattern->ptr[i];iptr_ik<mm; ++iptr_ik) {
267 k=A_p->pattern->index[iptr_ik];
268 rtmp-=A_p->val[iptr_ik]*x[k];
269 }
270 x[i]=rtmp*diag[i];
271 }
272 } else if (n_block==2) {
273 Paso_BlockOps_MViP_2(&diag[0], &x[0]);
274 for (i = 1; i < n; ++i) {
275 mm=ptr_main[i];
276 for (iptr_ik=A_p->pattern->ptr[i];iptr_ik<mm; ++iptr_ik) {
277 k=A_p->pattern->index[iptr_ik];
278 Paso_BlockOps_SMV_2(&x[2*i], &A_p->val[4*iptr_ik], &x[2*k]);
279 }
280 Paso_BlockOps_MViP_2(&diag[4*i], &x[2*i]);
281 }
282 } else if (n_block==3) {
283 Paso_BlockOps_MViP_3(&diag[0], &x[0]);
284 for (i = 1; i < n; ++i) {
285 mm=ptr_main[i];
286 for (iptr_ik=A_p->pattern->ptr[i];iptr_ik<mm; ++iptr_ik) {
287 k=A_p->pattern->index[iptr_ik];
288 Paso_BlockOps_SMV_3(&x[3*i], &A_p->val[9*iptr_ik], &x[3*k]);
289 }
290 Paso_BlockOps_MViP_3(&diag[9*i], &x[3*i]);
291 }
292 } else {
293 Paso_BlockOps_solve_N(n_block, &x[0], &diag[0], &pivot[0], &failed);
294 for (i = 1; i < n; ++i) {
295 mm=ptr_main[i];
296 for (iptr_ik=A_p->pattern->ptr[i];iptr_ik<mm; ++iptr_ik) {
297 k=A_p->pattern->index[iptr_ik];
298 Paso_BlockOps_SMV_N(n_block, &x[n_block*i], &A_p->val[block_len*iptr_ik], &x[n_block*k]);
299 }
300 Paso_BlockOps_solve_N(n_block, &x[n_block*i], &diag[block_len*i], &pivot[n_block*i], &failed);
301 }
302 }
303
304 /* backward sweeps */
305 if (n_block==1) {
306 for (i = n-2; i > -1; --i) {
307 mm=ptr_main[i];
308 rtmp=x[i]*A_p->val[mm];
309 for (iptr_ik=mm+1; iptr_ik < A_p->pattern->ptr[i+1]; ++iptr_ik) {
310 k=A_p->pattern->index[iptr_ik];
311 rtmp-=A_p->val[iptr_ik]*x[k];
312 }
313 x[i]=diag[i]*rtmp;
314 }
315
316 } else if (n_block==2) {
317 for (i = n-2; i > -1; --i) {
318 mm=ptr_main[i];
319 Paso_BlockOps_MViP_2(&A_p->val[4*mm], &x[2*i]);
320 for (iptr_ik=mm+1; iptr_ik < A_p->pattern->ptr[i+1]; ++iptr_ik) {
321 k=A_p->pattern->index[iptr_ik];
322 Paso_BlockOps_SMV_2(&x[2*i], &A_p->val[4*iptr_ik], &x[2*k]);
323 }
324 Paso_BlockOps_MViP_2(&diag[i*4], &x[2*i]);
325
326 }
327 } else if (n_block==3) {
328 for (i = n-2; i > -1; --i) {
329 mm=ptr_main[i];
330 Paso_BlockOps_MViP_3(&A_p->val[9*mm], &x[3*i]);
331 for (iptr_ik=mm+1; iptr_ik < A_p->pattern->ptr[i+1]; ++iptr_ik) {
332 k=A_p->pattern->index[iptr_ik];
333 Paso_BlockOps_SMV_3(&x[3*i], &A_p->val[9*iptr_ik], &x[3*k]);
334 }
335 Paso_BlockOps_MViP_3(&diag[i*9], &x[3*i]);
336 }
337 } else {
338 double *y=new double[n_block];
339 for (i = n-2; i > -1; --i) {
340 mm=ptr_main[i];
341 Paso_BlockOps_MV_N(n_block, &y[0], &A_p->val[block_len*mm], &x[n_block*i]);
342 for (iptr_ik=mm+1; iptr_ik < A_p->pattern->ptr[i+1]; ++iptr_ik) {
343 k=A_p->pattern->index[iptr_ik];
344 Paso_BlockOps_SMV_N(n_block, &y[0], &A_p->val[block_len*iptr_ik], &x[n_block*k]);
345 }
346 Paso_BlockOps_Cpy_N(n_block ,&x[n_block*i], &y[0]);
347 Paso_BlockOps_solve_N(n_block, &x[n_block*i], &diag[i*block_len], &pivot[i*n_block], &failed);
348 }
349 delete[] y;
350 }
351
352 if (failed > 0) {
353 Esys_setError(ZERO_DIVISION_ERROR, "Paso_Preconditioner_LocalSmoother_Sweep_sequential: non-regular main diagonal block.");
354 }
355
356 return;
357 }
358
359 void Paso_Preconditioner_LocalSmoother_Sweep_colored(paso::SparseMatrix* A_p, Paso_Preconditioner_LocalSmoother * smoother, double * x)
360 {
361 const dim_t n=A_p->numRows;
362 const dim_t n_block=A_p->row_block_size;
363 double *diag = smoother->diag;
364 index_t* pivot = smoother->pivot;
365 const dim_t block_len=A_p->block_size;
366 double *y;
367
368 register dim_t i,k;
369 register index_t color,iptr_ik, mm;
370 register double rtmp;
371 int failed = 0;
372
373 const index_t* coloring = paso::Pattern_borrowColoringPointer(A_p->pattern);
374 const dim_t num_colors = paso::Pattern_getNumColors(A_p->pattern);
375 const index_t* ptr_main = paso::SparseMatrix_borrowMainDiagonalPointer(A_p);
376
377 (void)pivot; /* These vars are dropped by some macros*/
378 (void)block_len;
379
380 #pragma omp parallel private(mm, i,iptr_ik,k,rtmp, color, y)
381 {
382 if (n_block>3) {
383 y=new double[n_block];
384 } else {
385 y=NULL;
386 }
387 /* forward substitution */
388
389 /* color = 0 */
390 if (n_block==1) {
391 #pragma omp for schedule(static)
392 for (i = 0; i < n; ++i) {
393 if (coloring[i]== 0 ) x[i]*=diag[i];
394 }
395 } else if (n_block==2) {
396 #pragma omp for schedule(static)
397 for (i = 0; i < n; ++i) {
398 if (coloring[i]== 0 ) Paso_BlockOps_MViP_2(&diag[i*4], &x[2*i]);
399 }
400 } else if (n_block==3) {
401 #pragma omp for schedule(static)
402 for (i = 0; i < n; ++i) {
403 if (coloring[i]==0) Paso_BlockOps_MViP_3(&diag[i*9], &x[3*i]);
404 }
405 } else {
406 #pragma omp for schedule(static)
407 for (i = 0; i < n; ++i) {
408 if (coloring[i]==0) Paso_BlockOps_solve_N(n_block, &x[n_block*i], &diag[block_len*i], &pivot[n_block*i], &failed);
409 }
410 }
411
412
413 for (color=1;color<num_colors;++color) {
414
415 if (n_block==1) {
416 #pragma omp for schedule(static)
417 for (i = 0; i < n; ++i) {
418 if (coloring[i]==color) {
419 /* x_i=x_i-a_ik*x_k */
420 rtmp=x[i];
421 for (iptr_ik=A_p->pattern->ptr[i];iptr_ik<A_p->pattern->ptr[i+1]; ++iptr_ik) {
422 k=A_p->pattern->index[iptr_ik];
423 if (coloring[k]<color) rtmp-=A_p->val[iptr_ik]*x[k];
424 }
425 x[i]=diag[i]*rtmp;
426 }
427 }
428 } else if (n_block==2) {
429 #pragma omp for schedule(static)
430 for (i = 0; i < n; ++i) {
431 if (coloring[i]==color) {
432 for (iptr_ik=A_p->pattern->ptr[i];iptr_ik<A_p->pattern->ptr[i+1]; ++iptr_ik) {
433 k=A_p->pattern->index[iptr_ik];
434 if (coloring[k]<color) Paso_BlockOps_SMV_2(&x[2*i], &A_p->val[4*iptr_ik], &x[2*k]);
435 }
436 Paso_BlockOps_MViP_2(&diag[4*i], &x[2*i]);
437 }
438 }
439 } else if (n_block==3) {
440 #pragma omp for schedule(static)
441 for (i = 0; i < n; ++i) {
442 if (coloring[i]==color) {
443 for (iptr_ik=A_p->pattern->ptr[i];iptr_ik<A_p->pattern->ptr[i+1]; ++iptr_ik) {
444 k=A_p->pattern->index[iptr_ik];
445 if (coloring[k]<color) Paso_BlockOps_SMV_3(&x[3*i], &A_p->val[9*iptr_ik], &x[3*k]);
446 }
447 Paso_BlockOps_MViP_3(&diag[9*i], &x[3*i]);
448 }
449 }
450 } else {
451 #pragma omp for schedule(static)
452 for (i = 0; i < n; ++i) {
453 if (coloring[i] == color) {
454 for (iptr_ik=A_p->pattern->ptr[i];iptr_ik<A_p->pattern->ptr[i+1]; ++iptr_ik) {
455 k=A_p->pattern->index[iptr_ik];
456 if (coloring[k]<color) Paso_BlockOps_SMV_N(n_block, &x[n_block*i], &A_p->val[block_len*iptr_ik], &x[n_block*k]);
457 }
458 Paso_BlockOps_solve_N(n_block, &x[n_block*i], &diag[block_len*i], &pivot[n_block*i], &failed);
459 }
460 }
461 }
462 } /* end of coloring loop */
463
464 /* backward substitution */
465 for (color=(num_colors)-2 ;color>-1;--color) { /* Note: color=(num_colors)-1 is not required */
466 if (n_block==1) {
467 #pragma omp for schedule(static)
468 for (i = 0; i < n; ++i) {
469 if (coloring[i]==color) {
470 mm=ptr_main[i];
471 rtmp=A_p->val[mm]*x[i];
472 for (iptr_ik=A_p->pattern->ptr[i];iptr_ik<A_p->pattern->ptr[i+1]; ++iptr_ik) {
473 k=A_p->pattern->index[iptr_ik];
474 if (coloring[k]>color) rtmp-=A_p->val[iptr_ik]*x[k];
475 }
476 x[i]= rtmp*diag[i];
477 }
478 }
479 } else if (n_block==2) {
480 #pragma omp for schedule(static)
481 for (i = 0; i < n; ++i) {
482 if (coloring[i]==color) {
483 mm=ptr_main[i];
484 Paso_BlockOps_MViP_2(&A_p->val[4*mm], &x[2*i]);
485 for (iptr_ik=A_p->pattern->ptr[i];iptr_ik<A_p->pattern->ptr[i+1]; ++iptr_ik) {
486 k=A_p->pattern->index[iptr_ik];
487 if (coloring[k]>color) Paso_BlockOps_SMV_2(&x[2*i], &A_p->val[4*iptr_ik], &x[2*k]);
488 }
489 Paso_BlockOps_MViP_2(&diag[4*i], &x[2*i]);
490 }
491 }
492 } else if (n_block==3) {
493 #pragma omp for schedule(static)
494 for (i = 0; i < n; ++i) {
495 if (coloring[i]==color) {
496 mm=ptr_main[i];
497 Paso_BlockOps_MViP_3(&A_p->val[9*mm], &x[3*i]);
498 for (iptr_ik=A_p->pattern->ptr[i];iptr_ik<A_p->pattern->ptr[i+1]; ++iptr_ik) {
499 k=A_p->pattern->index[iptr_ik];
500 if (coloring[k]>color) Paso_BlockOps_SMV_3(&x[3*i], &A_p->val[9*iptr_ik], &x[3*k]);
501 }
502 Paso_BlockOps_MViP_3(&diag[9*i], &x[3*i]);
503 }
504 }
505 } else {
506 #pragma omp for schedule(static)
507 for (i = 0; i < n; ++i) {
508 if (coloring[i]==color) {
509 mm=ptr_main[i];
510 Paso_BlockOps_MV_N(n_block, &y[0], &A_p->val[block_len*mm], &x[n_block*i]);
511 for (iptr_ik=A_p->pattern->ptr[i];iptr_ik<A_p->pattern->ptr[i+1]; ++iptr_ik) {
512 k=A_p->pattern->index[iptr_ik];
513 if (coloring[k]>color) Paso_BlockOps_SMV_N(n_block, &y[0], &A_p->val[block_len*iptr_ik], &x[n_block*k]);
514 }
515 Paso_BlockOps_Cpy_N(n_block ,&x[n_block*i], &y[0]);
516 Paso_BlockOps_solve_N(n_block, &x[n_block*i], &diag[i*block_len], &pivot[i*n_block], &failed);
517 }
518 }
519 }
520 }
521 delete[] y;
522 }
523 if (failed > 0) {
524 Esys_setError(ZERO_DIVISION_ERROR, "Paso_Preconditioner_LocalSmoother_Sweep_colored: non-regular main diagonal block.");
525 }
526 return;
527 }

Properties

Name Value
svn:mergeinfo /branches/amg_from_3530/paso/src/Smoother.cpp:3531-3826 /branches/lapack2681/paso/src/Smoother.cpp:2682-2741 /branches/pasowrap/paso/src/Smoother.cpp:3661-3674 /branches/py3_attempt2/paso/src/Smoother.cpp:3871-3891 /branches/restext/paso/src/Smoother.cpp:2610-2624 /branches/ripleygmg_from_3668/paso/src/Smoother.cpp:3669-3791 /branches/stage3.0/paso/src/Smoother.cpp:2569-2590 /branches/symbolic_from_3470/paso/src/Smoother.cpp:3471-3974 /branches/symbolic_from_3470/ripley/test/python/paso/src/Smoother.cpp:3517-3974 /release/3.0/paso/src/Smoother.cpp:2591-2601 /trunk/paso/src/Smoother.cpp:4257-4344 /trunk/ripley/test/python/paso/src/Smoother.cpp:3480-3515

  ViewVC Help
Powered by ViewVC 1.1.26