/[escript]/trunk/paso/src/SparseMatrix_nullifyRowsAndCols.c
ViewVC logotype

Annotation of /trunk/paso/src/SparseMatrix_nullifyRowsAndCols.c

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1564 - (hide annotations)
Thu May 22 09:31:33 2008 UTC (11 years, 9 months ago) by gross
File MIME type: text/plain
File size: 6208 byte(s)
some openmp dynamic scheduling for MVM.
1 ksteube 1315
2     /* $Id: SparseMatrix_nullifyRowsAndCols.c 1306 2007-09-18 05:51:09Z ksteube $ */
3    
4     /*******************************************************
5     *
6     * Copyright 2003-2007 by ACceSS MNRF
7     * Copyright 2007 by University of Queensland
8     *
9     * http://esscc.uq.edu.au
10     * Primary Business: Queensland, Australia
11     * Licensed under the Open Software License version 3.0
12     * http://www.opensource.org/licenses/osl-3.0.php
13     *
14     *******************************************************/
15    
16     /**************************************************************/
17    
18     /* Paso: SparseMatrix */
19    
20     /* nullify rows and columns in the matrix */
21    
22     /* the rows and columns are marked by positive values in */
23     /* mask_row and mask_col. Values on the main diagonal */
24     /* which are marked to set to zero by both mask_row and */
25     /* mask_col are set to main_diagonal_value */
26    
27    
28     /**************************************************************/
29    
30     /* Copyrights by ACcESS Australia 2003, 2007 */
31     /* Author: gross@access.edu.au */
32    
33     /**************************************************************/
34    
35     #include "Paso.h"
36     #include "SparseMatrix.h"
37    
38     void Paso_SparseMatrix_nullifyRowsAndCols_CSC_BLK1(Paso_SparseMatrix* A, double* mask_row, double* mask_col, double main_diagonal_value) {
39     index_t index_offset=(A->type & MATRIX_FORMAT_OFFSET1 ? 1:0);
40     index_t irow, iptr, icol;
41     #pragma omp parallel for private(irow, iptr,icol) schedule(static)
42     for (icol=0;icol< A->pattern->numOutput;icol++) {
43 gross 1552 #pragma ivdep
44 ksteube 1315 for (iptr=A->pattern->ptr[icol]-index_offset;iptr<A->pattern->ptr[icol+1]-index_offset; iptr++) {
45     irow=A->pattern->index[iptr]-index_offset;
46     if (mask_col[icol]>0. || mask_row[irow]>0. ) {
47     if (irow==icol) {
48     A->val[iptr]=main_diagonal_value;
49     } else {
50     A->val[iptr]=0;
51     }
52     }
53     }
54     }
55     }
56     void Paso_SparseMatrix_nullifyRowsAndCols_CSR_BLK1(Paso_SparseMatrix* A, double* mask_row, double* mask_col, double main_diagonal_value) {
57     index_t index_offset=(A->type & MATRIX_FORMAT_OFFSET1 ? 1:0);
58     index_t irow, iptr, icol;
59     #pragma omp parallel for private(irow, iptr,icol) schedule(static)
60     for (irow=0;irow< A->pattern->numOutput;irow++) {
61     /* TODO: test mask_row here amd not inside every loop */
62 gross 1552 #pragma ivdep
63 ksteube 1315 for (iptr=A->pattern->ptr[irow]-index_offset;iptr<A->pattern->ptr[irow+1]-index_offset; iptr++) {
64     icol=A->pattern->index[iptr]-index_offset;
65     if (mask_col[icol]>0. || mask_row[irow]>0. ) {
66     if (irow==icol) {
67     A->val[iptr]=main_diagonal_value;
68     } else {
69     A->val[iptr]=0;
70     }
71     }
72     }
73     }
74     }
75     void Paso_SparseMatrix_nullifyRowsAndCols_CSC(Paso_SparseMatrix* A, double* mask_row, double* mask_col, double main_diagonal_value) {
76     index_t index_offset=(A->type & MATRIX_FORMAT_OFFSET1 ? 1:0);
77     index_t ir,icol,iptr,icb,irb,irow,ic,l;
78     #pragma omp parallel for private(l,irow, iptr,icol,ic,irb,icb) schedule(static)
79     for (ic=0;ic< A->pattern->numOutput;ic++) {
80     for (iptr=A->pattern->ptr[ic]-index_offset;iptr<A->pattern->ptr[ic+1]-index_offset; iptr++) {
81     for (irb=0;irb< A->row_block_size;irb++) {
82     irow=irb+A->row_block_size*(A->pattern->index[iptr]-index_offset);
83 gross 1552 #pragma ivdep
84 ksteube 1315 for (icb=0;icb< A->col_block_size;icb++) {
85     icol=icb+A->col_block_size*ic;
86     if (mask_col[icol]>0. || mask_row[irow]>0. ) {
87     l=iptr*A->block_size+irb+A->row_block_size*icb;
88     if (irow==icol) {
89     A->val[l]=main_diagonal_value;
90     } else {
91     A->val[l]=0;
92     }
93     }
94     }
95     }
96     }
97     }
98     }
99     void Paso_SparseMatrix_nullifyRowsAndCols_CSR(Paso_SparseMatrix* A, double* mask_row, double* mask_col, double main_diagonal_value) {
100     index_t index_offset=(A->type & MATRIX_FORMAT_OFFSET1 ? 1:0);
101     index_t ir,icol,iptr,icb,irb,irow,ic,l;
102     #pragma omp parallel for private(l,irow, iptr,icol,ir,irb,icb) schedule(static)
103     for (ir=0;ir< A->pattern->numOutput;ir++) {
104     for (iptr=A->pattern->ptr[ir]-index_offset;iptr<A->pattern->ptr[ir+1]-index_offset; iptr++) {
105     for (irb=0;irb< A->row_block_size;irb++) {
106     irow=irb+A->row_block_size*ir;
107 gross 1552 #pragma ivdep
108 ksteube 1315 for (icb=0;icb< A->col_block_size;icb++) {
109     icol=icb+A->col_block_size*(A->pattern->index[iptr]-index_offset);
110     if (mask_col[icol]>0. || mask_row[irow]>0. ) {
111     l=iptr*A->block_size+irb+A->row_block_size*icb;
112     if (irow==icol) {
113     A->val[l]=main_diagonal_value;
114     } else {
115     A->val[l]=0;
116     }
117     }
118     }
119     }
120     }
121     }
122     }
123     void Paso_SparseMatrix_nullifyRows_CSR_BLK1(Paso_SparseMatrix* A, double* mask_row, double main_diagonal_value) {
124     index_t index_offset=(A->type & MATRIX_FORMAT_OFFSET1 ? 1:0);
125     index_t irow, iptr, icol;
126     #pragma omp parallel for private(irow, iptr,icol) schedule(static)
127     for (irow=0;irow< A->pattern->numOutput;irow++) {
128     if (mask_row[irow]>0.) {
129 gross 1552 #pragma ivdep
130 ksteube 1315 for (iptr=A->pattern->ptr[irow]-index_offset;iptr<A->pattern->ptr[irow+1]-index_offset; iptr++) {
131     icol=A->pattern->index[iptr]-index_offset;
132     if (irow==icol) {
133     A->val[iptr]=main_diagonal_value;
134     } else {
135     A->val[iptr]=0;
136     }
137     }
138     }
139     }
140     }
141     void Paso_SparseMatrix_nullifyRows_CSR(Paso_SparseMatrix* A, double* mask_row, double main_diagonal_value) {
142     index_t index_offset=(A->type & MATRIX_FORMAT_OFFSET1 ? 1:0);
143     index_t ir,icol,iptr,icb,irb,irow,ic,l;
144     #pragma omp parallel for private(l,irow, iptr,icol,ir,irb,icb) schedule(static)
145     for (ir=0;ir< A->pattern->numOutput;ir++) {
146     for (iptr=A->pattern->ptr[ir]-index_offset;iptr<A->pattern->ptr[ir+1]-index_offset; iptr++) {
147     for (irb=0;irb< A->row_block_size;irb++) {
148     irow=irb+A->row_block_size*ir;
149     if (mask_row[irow]>0. ) {
150 gross 1552 #pragma ivdep
151 ksteube 1315 for (icb=0;icb< A->col_block_size;icb++) {
152     icol=icb+A->col_block_size*(A->pattern->index[iptr]-index_offset);
153     l=iptr*A->block_size+irb+A->row_block_size*icb;
154     if (irow==icol) {
155     A->val[l]=main_diagonal_value;
156     } else {
157     A->val[l]=0;
158     }
159     }
160     }
161     }
162     }
163     }
164     }

  ViewVC Help
Powered by ViewVC 1.1.26