/[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 5148 - (show annotations)
Mon Sep 15 01:25:23 2014 UTC (4 years, 11 months ago) by caltinay
File size: 22558 byte(s)
Merging ripley diagonal storage + CUDA support into trunk.
Options file version has been incremented due to new options
'cuda' and 'nvccflags'.

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

Properties

Name Value
svn:mergeinfo /branches/amg_from_3530/paso/src/SparseMatrix_MatrixVector.cpp:3531-3826 /branches/diaplayground/paso/src/SparseMatrix_MatrixVector.cpp:4940-5147 /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