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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 4829 - (show annotations)
Thu Apr 3 04:02:53 2014 UTC (5 years, 5 months ago) by caltinay
File size: 22446 byte(s)
checkpointing some SparseMatrix cleanup.

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: raw scaled vector update operation:
21 * out = alpha * A * in + beta * out
22 ****************************************************************************/
23
24 /* Author: Lutz Gross, l.gross@uq.edu.au */
25
26 /****************************************************************************/
27
28 #include "SparseMatrix.h"
29
30 namespace paso {
31
32 // forward declaration
33 void SparseMatrix_MatrixVector_CSR_OFFSET0_stripe(double alpha, dim_t nRows,
34 dim_t row_block_size,
35 dim_t col_block_size,
36 const index_t* ptr,
37 const index_t* index,
38 const double* val,
39 const double* in,
40 double beta, double* out);
41
42 /* CSC format with offset 0 */
43 void SparseMatrix_MatrixVector_CSC_OFFSET0(double alpha,
44 const_SparseMatrix_ptr A,
45 const double* in,
46 double beta, double* out)
47 {
48 if (ABS(beta) > 0) {
49 if (beta != 1.) {
50 #pragma omp parallel for schedule(static)
51 for (index_t irow=0; irow < A->numRows*A->row_block_size; irow++) {
52 out[irow] *= beta;
53 }
54 }
55 } else {
56 #pragma omp parallel for schedule(static)
57 for (index_t irow=0; irow < A->numRows*A->row_block_size; irow++) {
58 out[irow] = 0;
59 }
60 }
61
62 if (A->pattern->isEmpty())
63 return;
64
65 if (ABS(alpha) > 0) {
66 if (A->col_block_size==1 && A->row_block_size==1) {
67 /* TODO: parallelize (good luck!) */
68 for (index_t icol=0; icol < A->pattern->numOutput; ++icol) {
69 #pragma ivdep
70 for (index_t iptr=A->pattern->ptr[icol];
71 iptr < A->pattern->ptr[icol+1]; ++iptr) {
72 out[A->pattern->index[iptr]] += alpha * A->val[iptr] * in[icol];
73 }
74 }
75 } else if (A->col_block_size==2 && A->row_block_size==2) {
76 /* TODO: parallelize */
77 for (index_t ic=0; ic < A->pattern->numOutput; ic++) {
78 #pragma ivdep
79 for (index_t iptr=A->pattern->ptr[ic];
80 iptr < A->pattern->ptr[ic+1]; iptr++) {
81 const index_t ir = 2*A->pattern->index[iptr];
82 out[ ir] += alpha * (A->val[iptr*4 ]*in[ic] + A->val[iptr*4+2]*in[1+ic]);
83 out[1+ir] += alpha * (A->val[iptr*4+1]*in[ic] + A->val[iptr*4+3]*in[1+ic]);
84 }
85 }
86 } else if (A->col_block_size==3 && A->row_block_size==3) {
87 /* TODO: parallelize */
88 for (index_t ic=0; ic < A->pattern->numOutput; ic++) {
89 #pragma ivdep
90 for (index_t iptr=A->pattern->ptr[ic];
91 iptr < A->pattern->ptr[ic+1]; iptr++) {
92 const index_t ir = 3*A->pattern->index[iptr];
93 out[ ir] += alpha * (A->val[iptr*9 ]*in[ic] + A->val[iptr*9+3]*in[1+ic] + A->val[iptr*9+6]*in[2+ic]);
94 out[1+ir] += alpha * (A->val[iptr*9+1]*in[ic] + A->val[iptr*9+4]*in[1+ic] + A->val[iptr*9+7]*in[2+ic]);
95 out[2+ir] += alpha * (A->val[iptr*9+2]*in[ic] + A->val[iptr*9+5]*in[1+ic] + A->val[iptr*9+8]*in[2+ic]);
96 }
97 }
98 } else {
99 /* TODO: parallelize */
100 for (index_t ic=0; ic < A->pattern->numOutput; ic++) {
101 for (index_t iptr=A->pattern->ptr[ic];
102 iptr < A->pattern->ptr[ic+1]; iptr++) {
103 for (index_t irb=0; irb < A->row_block_size; irb++) {
104 const index_t irow = irb+A->row_block_size*A->pattern->index[iptr];
105 #pragma ivdep
106 for (index_t icb=0; icb < A->col_block_size; icb++) {
107 const index_t icol=icb+A->col_block_size*ic;
108 out[irow] += alpha * A->val[iptr*A->block_size+irb+A->row_block_size*icb] * in[icol];
109 }
110 }
111 }
112 }
113 } // blocksizes
114 } // alpha > 0
115 }
116
117 /* CSC format with offset 1 */
118 void SparseMatrix_MatrixVector_CSC_OFFSET1(double alpha,
119 const_SparseMatrix_ptr A,
120 const double* in,
121 double beta, double* out)
122 {
123 if (ABS(beta) > 0) {
124 if (beta != 1.) {
125 #pragma omp parallel for schedule(static)
126 for (index_t irow=0; irow < A->numRows*A->row_block_size; irow++) {
127 out[irow] *= beta;
128 }
129 }
130 } else {
131 #pragma omp parallel for schedule(static)
132 for (index_t irow=0; irow < A->numRows*A->row_block_size; irow++) {
133 out[irow] = 0;
134 }
135 }
136
137 if (ABS(alpha) > 0) {
138 if (A->col_block_size==1 && A->row_block_size==1) {
139 /* TODO: parallelize (good luck!) */
140 for (index_t icol=0; icol < A->pattern->numOutput; icol++) {
141 #pragma ivdep
142 for (index_t iptr=A->pattern->ptr[icol]-1;
143 iptr<A->pattern->ptr[icol+1]-1; iptr++) {
144 out[A->pattern->index[iptr]-1]+= alpha * A->val[iptr] * in[icol];
145 }
146 }
147 } else if (A->col_block_size==2 && A->row_block_size==2) {
148 /* TODO: parallelize */
149 for (index_t ic=0; ic < A->pattern->numOutput; ic++) {
150 for (index_t iptr=A->pattern->ptr[ic]-1;
151 iptr < A->pattern->ptr[ic+1]-1; iptr++) {
152 const index_t ir=2*(A->pattern->index[iptr]-1);
153 out[ ir] += alpha * (A->val[iptr*4 ]*in[ic] + A->val[iptr*4+2]*in[1+ic]);
154 out[1+ir] += alpha * (A->val[iptr*4+1]*in[ic] + A->val[iptr*4+3]*in[1+ic]);
155 }
156 }
157 } else if (A->col_block_size==3 && A->row_block_size==3) {
158 /* TODO: parallelize */
159 for (index_t ic=0; ic < A->pattern->numOutput; ic++) {
160 #pragma ivdep
161 for (index_t iptr=A->pattern->ptr[ic]-1;
162 iptr < A->pattern->ptr[ic+1]-1; iptr++) {
163 const index_t ir=3*(A->pattern->index[iptr]-1);
164 out[ ir] += alpha * (A->val[iptr*9 ]*in[ic] + A->val[iptr*9+3]*in[1+ic] + A->val[iptr*9+6]*in[2+ic]);
165 out[1+ir] += alpha * (A->val[iptr*9+1]*in[ic] + A->val[iptr*9+4]*in[1+ic] + A->val[iptr*9+7]*in[2+ic]);
166 out[2+ir] += alpha * (A->val[iptr*9+2]*in[ic] + A->val[iptr*9+5]*in[1+ic] + A->val[iptr*9+8]*in[2+ic]);
167 }
168 }
169 } else {
170 /* TODO: parallelize */
171 for (index_t ic=0; ic < A->pattern->numOutput; ic++) {
172 for (index_t iptr=A->pattern->ptr[ic]-1;
173 iptr < A->pattern->ptr[ic+1]-1; iptr++) {
174 for (index_t irb=0; irb < A->row_block_size; irb++) {
175 const index_t irow=irb+A->row_block_size*(A->pattern->index[iptr]-1);
176 #pragma ivdep
177 for (index_t icb=0; icb < A->col_block_size; icb++) {
178 const index_t icol=icb+A->col_block_size*ic;
179 out[irow] += alpha * A->val[iptr*A->block_size+irb+A->row_block_size*icb] * in[icol];
180 }
181 }
182 }
183 }
184 } // blocksizes
185 } // alpha > 0
186 }
187
188 /* CSR format with offset 1 */
189 void SparseMatrix_MatrixVector_CSR_OFFSET1(double alpha,
190 const_SparseMatrix_ptr A,
191 const double* in,
192 double beta, double* out)
193 {
194 if (ABS(beta) > 0) {
195 if (beta != 1.) {
196 #pragma omp parallel for schedule(static)
197 for (index_t irow=0; irow < A->numRows*A->row_block_size; irow++) {
198 out[irow] *= beta;
199 }
200 }
201 } else {
202 #pragma omp parallel for schedule(static)
203 for (index_t irow=0; irow < A->numRows*A->row_block_size; irow++) {
204 out[irow] = 0;
205 }
206 }
207
208 if (ABS(alpha) > 0) {
209 if (A->col_block_size==1 && A->row_block_size==1) {
210 #pragma omp parallel for schedule(static)
211 for (index_t irow=0; irow < A->pattern->numOutput; irow++) {
212 double reg=0.;
213 #pragma ivdep
214 for (index_t iptr=A->pattern->ptr[irow]-1;
215 iptr < A->pattern->ptr[irow+1]-1; ++iptr) {
216 reg += A->val[iptr] * in[A->pattern->index[iptr]-1];
217 }
218 out[irow] += alpha * reg;
219 }
220 } else if (A->col_block_size==2 && A->row_block_size==2) {
221 #pragma omp parallel for schedule(static)
222 for (index_t ir=0; ir < A->pattern->numOutput; ir++) {
223 double reg1=0.;
224 double reg2=0.;
225 #pragma ivdep
226 for (index_t iptr=A->pattern->ptr[ir]-1;
227 iptr < A->pattern->ptr[ir+1]-1; iptr++) {
228 const index_t ic=2*(A->pattern->index[iptr]-1);
229 reg1 += A->val[iptr*4 ]*in[ic] + A->val[iptr*4+2]*in[1+ic];
230 reg2 += A->val[iptr*4+1]*in[ic] + A->val[iptr*4+3]*in[1+ic];
231 }
232 out[ 2*ir] += alpha * reg1;
233 out[1+2*ir] += alpha * reg2;
234 }
235 } else if (A->col_block_size==3 && A->row_block_size==3) {
236 #pragma omp parallel for schedule(static)
237 for (index_t ir=0; ir < A->pattern->numOutput; ir++) {
238 double reg1=0.;
239 double reg2=0.;
240 double reg3=0.;
241 #pragma ivdep
242 for (index_t iptr=A->pattern->ptr[ir]-1;
243 iptr < A->pattern->ptr[ir+1]-1; iptr++) {
244 const index_t ic=3*(A->pattern->index[iptr]-1);
245 reg1 += A->val[iptr*9 ]*in[ic] + A->val[iptr*9+3]*in[1+ic] + A->val[iptr*9+6]*in[2+ic];
246 reg2 += A->val[iptr*9+1]*in[ic] + A->val[iptr*9+4]*in[1+ic] + A->val[iptr*9+7]*in[2+ic];
247 reg3 += A->val[iptr*9+2]*in[ic] + A->val[iptr*9+5]*in[1+ic] + A->val[iptr*9+8]*in[2+ic];
248 }
249 out[ 3*ir] += alpha * reg1;
250 out[1+3*ir] += alpha * reg2;
251 out[2+3*ir] += alpha * reg3;
252 }
253 } else {
254 #pragma omp parallel for schedule(static)
255 for (index_t ir=0; ir < A->pattern->numOutput; ir++) {
256 for (index_t iptr=A->pattern->ptr[ir]-1;
257 iptr < A->pattern->ptr[ir+1]-1; iptr++) {
258 for (index_t irb=0; irb < A->row_block_size; irb++) {
259 double reg=0.;
260 #pragma ivdep
261 for (index_t icb=0; icb < A->col_block_size; icb++) {
262 const index_t icol=icb+A->col_block_size*(A->pattern->index[iptr]-1);
263 reg += A->val[iptr*A->block_size+irb+A->row_block_size*icb] * in[icol];
264 }
265 const index_t irow=irb+A->row_block_size*ir;
266 out[irow] += alpha * reg;
267 }
268 }
269 }
270 } // blocksizes
271 } // alpha > 0
272 }
273
274 /* CSR format with offset 0 */
275 void SparseMatrix_MatrixVector_CSR_OFFSET0(double alpha,
276 const_SparseMatrix_ptr A,
277 const double* in,
278 double beta, double* out)
279 {
280 /*#define PASO_DYNAMIC_SCHEDULING_MVM */
281
282 #if defined PASO_DYNAMIC_SCHEDULING_MVM && defined _OPENMP
283 #define USE_DYNAMIC_SCHEDULING
284 #endif
285
286 const dim_t nrow=A->numRows;
287 const dim_t np=omp_get_max_threads();
288 const dim_t len=nrow/np;
289
290 #ifdef USE_DYNAMIC_SCHEDULING
291 dim_t chunk_size=1;
292 char* chksz_chr=getenv("PASO_CHUNK_SIZE_MVM");
293 if (chksz_chr!=NULL) sscanf(chksz_chr, "%d",&chunk_size);
294 chunk_size=MIN(MAX(1,chunk_size),len);
295 dim_t n_chunks=nrow/chunk_size;
296 if (n_chunks*chunk_size<nrow) n_chunks+=1;
297
298 #pragma omp parallel for schedule(dynamic,1)
299 for (dim_t p=0; p < n_chunks; p++) {
300 const dim_t irow=chunk_size*p;
301 const dim_t local_n=MIN(chunk_size,nrow-chunk_size*p);
302 SparseMatrix_MatrixVector_CSR_OFFSET0_stripe(alpha, local_n,
303 A->row_block_size, A->col_block_size, &(A->pattern->ptr[irow]),
304 A->pattern->index, A->val, in, beta,
305 &out[irow*A->row_block_size]);
306 }
307
308 #else // static scheduling
309 const dim_t rest=nrow-len*np;
310
311 #pragma omp parallel for
312 for (dim_t p=0; p < np; p++) {
313 const dim_t irow=len*p+MIN(p,rest);
314 const dim_t local_n=len+(p<rest ? 1 :0 );
315 SparseMatrix_MatrixVector_CSR_OFFSET0_stripe(alpha, local_n,
316 A->row_block_size, A->col_block_size, &(A->pattern->ptr[irow]),
317 A->pattern->index, A->val, in, beta,
318 &out[irow*A->row_block_size]);
319 }
320 #endif // scheduling
321 }
322
323 /* CSR format with offset 0 */
324 void SparseMatrix_MatrixVector_CSR_OFFSET0_stripe(const double alpha,
325 const dim_t nRows, const dim_t row_block_size,
326 const dim_t col_block_size, const index_t* ptr, const index_t* index,
327 const double* val, const double* in, const double beta, double* out)
328 {
329 if (ABS(beta) > 0) {
330 if (beta != 1.) {
331 for (index_t irow=0; irow < nRows*row_block_size; irow++)
332 out[irow] *= beta;
333 }
334 } else {
335 for (index_t irow=0; irow < nRows*row_block_size; irow++)
336 out[irow] = 0;
337 }
338
339 if (ABS(alpha) > 0) {
340 if (col_block_size==1 && row_block_size ==1) {
341 for (index_t irow=0; irow < nRows; ++irow) {
342 double reg=0.;
343 #pragma ivdep
344 for (index_t iptr=ptr[irow]; iptr<ptr[irow+1]; ++iptr) {
345 reg += val[iptr] * in[index[iptr]];
346 }
347 out[irow] += alpha * reg;
348 }
349 } else if (col_block_size==2 && row_block_size==2) {
350 for (index_t ir=0; ir < nRows; ir++) {
351 double reg1=0.;
352 double reg2=0.;
353 #pragma ivdep
354 for (index_t iptr=ptr[ir]; iptr<ptr[ir+1]; iptr++) {
355 const index_t ic=2*index[iptr];
356 const index_t Aiptr=iptr*4;
357 const double in1=in[ic];
358 const double in2=in[1+ic];
359 const double A00=val[Aiptr ];
360 const double A10=val[Aiptr+1];
361 const double A01=val[Aiptr+2];
362 const double A11=val[Aiptr+3];
363 reg1 += A00*in1 + A01*in2;
364 reg2 += A10*in1 + A11*in2;
365 }
366 out[ 2*ir] += alpha * reg1;
367 out[1+2*ir] += alpha * reg2;
368 }
369 } else if (col_block_size==3 && row_block_size==3) {
370 for (index_t ir=0; ir < nRows; ir++) {
371 double reg1=0.;
372 double reg2=0.;
373 double reg3=0.;
374 #pragma ivdep
375 for (index_t iptr=ptr[ir]; iptr < ptr[ir+1]; iptr++) {
376 const index_t ic=3*index[iptr];
377 const index_t Aiptr=iptr*9;
378 const double in1=in[ic];
379 const double in2=in[1+ic];
380 const double in3=in[2+ic];
381 const double A00=val[Aiptr ];
382 const double A10=val[Aiptr+1];
383 const double A20=val[Aiptr+2];
384 const double A01=val[Aiptr+3];
385 const double A11=val[Aiptr+4];
386 const double A21=val[Aiptr+5];
387 const double A02=val[Aiptr+6];
388 const double A12=val[Aiptr+7];
389 const double A22=val[Aiptr+8];
390 reg1 += A00*in1 + A01*in2 + A02*in3;
391 reg2 += A10*in1 + A11*in2 + A12*in3;
392 reg3 += A20*in1 + A21*in2 + A22*in3;
393 }
394 out[ 3*ir] += alpha * reg1;
395 out[1+3*ir] += alpha * reg2;
396 out[2+3*ir] += alpha * reg3;
397 }
398 } else { // blocksizes > 3
399 const dim_t block_size=col_block_size*row_block_size;
400 for (index_t ir=0; ir < nRows; ir++) {
401 for (index_t iptr=ptr[ir]; iptr < ptr[ir+1]; iptr++) {
402 for (index_t irb=0; irb < row_block_size; irb++) {
403 double reg=0.;
404 #pragma ivdep
405 for (index_t icb=0; icb < col_block_size; icb++) {
406 const index_t icol=icb+col_block_size*index[iptr];
407 reg += val[iptr*block_size+irb+row_block_size*icb] * in[icol];
408 }
409 const index_t irow=irb+row_block_size*ir;
410 out[irow] += alpha * reg;
411 }
412 }
413 }
414 }
415 }
416 }
417
418 /* CSR format with offset 0 (diagonal only) */
419 void SparseMatrix_MatrixVector_CSR_OFFSET0_DIAG(double alpha,
420 const_SparseMatrix_ptr A,
421 const double* in,
422 double beta, double* out)
423 {
424 if (ABS(beta) > 0) {
425 if (beta != 1.) {
426 #pragma omp parallel for schedule(static)
427 for (index_t irow=0; irow < A->numRows*A->row_block_size; irow++)
428 out[irow] *= beta;
429 }
430 } else {
431 #pragma omp parallel for schedule(static)
432 for (index_t irow=0; irow < A->numRows*A->row_block_size; irow++)
433 out[irow] = 0;
434 }
435
436 if (ABS(alpha) > 0) {
437 if (A->block_size == 1) {
438 #pragma omp parallel for schedule(static)
439 for (index_t irow=0; irow < A->pattern->numOutput; irow++) {
440 double reg=0.;
441 #pragma ivdep
442 for (index_t iptr=A->pattern->ptr[irow];
443 iptr < A->pattern->ptr[irow+1]; ++iptr) {
444 reg += A->val[iptr] * in[A->pattern->index[iptr]];
445 }
446 out[irow] += alpha * reg;
447 }
448 } else if (A->block_size == 2) {
449 #pragma omp parallel for schedule(static)
450 for (index_t ir=0; ir < A->pattern->numOutput; ir++) {
451 double reg1=0.;
452 double reg2=0.;
453 #pragma ivdep
454 for (index_t iptr=A->pattern->ptr[ir];
455 iptr < A->pattern->ptr[ir+1]; iptr++) {
456 const index_t ic = 2*A->pattern->index[iptr];
457 reg1 += A->val[iptr*2 ]*in[ ic];
458 reg2 += A->val[iptr*2+1]*in[1+ic];
459 }
460 out[ 2*ir] += alpha * reg1;
461 out[1+2*ir] += alpha * reg2;
462 }
463 } else if (A->block_size == 3) {
464 #pragma omp parallel for schedule(static)
465 for (index_t ir=0; ir < A->pattern->numOutput; ir++) {
466 double reg1=0.;
467 double reg2=0.;
468 double reg3=0.;
469 #pragma ivdep
470 for (index_t iptr=A->pattern->ptr[ir];
471 iptr < A->pattern->ptr[ir+1]; iptr++) {
472 const index_t ic = 3*A->pattern->index[iptr];
473 reg1 += A->val[iptr*3 ]*in[ ic];
474 reg2 += A->val[iptr*3+1]*in[1+ic];
475 reg3 += A->val[iptr*3+2]*in[2+ic];
476 }
477 out[ 3*ir] += alpha * reg1;
478 out[1+3*ir] += alpha * reg2;
479 out[2+3*ir] += alpha * reg3;
480 }
481 } else if (A->block_size == 4) {
482 #pragma omp parallel for schedule(static)
483 for (index_t ir=0; ir < A->pattern->numOutput; ir++) {
484 double reg1=0.;
485 double reg2=0.;
486 double reg3=0.;
487 double reg4=0.;
488 #pragma ivdep
489 for (index_t iptr=A->pattern->ptr[ir];
490 iptr < A->pattern->ptr[ir+1]; iptr++) {
491 const index_t ic = 4 * A->pattern->index[iptr];
492 reg1 += A->val[iptr*4 ]*in[ ic];
493 reg2 += A->val[iptr*4+1]*in[1+ic];
494 reg3 += A->val[iptr*4+2]*in[2+ic];
495 reg4 += A->val[iptr*4+3]*in[3+ic];
496 }
497 out[ 4*ir] += alpha * reg1;
498 out[1+4*ir] += alpha * reg2;
499 out[2+4*ir] += alpha * reg3;
500 out[3+4*ir] += alpha * reg4;
501 }
502 } else { // blocksize > 4
503 #pragma omp parallel for schedule(static)
504 for (index_t ir=0; ir < A->pattern->numOutput; ir++) {
505 for (index_t iptr=A->pattern->ptr[ir];
506 iptr < A->pattern->ptr[ir+1]; iptr++) {
507 for (index_t ib=0; ib < A->block_size; ib++) {
508 const index_t irow = ib+A->row_block_size*ir;
509 const index_t icol = ib+A->col_block_size*A->pattern->index[iptr];
510 out[irow] += alpha * A->val[iptr*A->block_size+ib] * in[icol];
511 }
512 }
513 }
514 } // blocksize
515 } // alpha > 0
516 }
517
518 } // namespace paso
519

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26