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

Diff of /trunk/paso/src/Pattern_coupling.c

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 1883 by artak, Tue Oct 14 04:46:27 2008 UTC revision 1890 by artak, Fri Oct 17 00:14:22 2008 UTC
# Line 42  Line 42 
42  #define IS_IN_MIS -3  #define IS_IN_MIS -3
43  #define IS_CONNECTED_TO_MIS -4  #define IS_CONNECTED_TO_MIS -4
44    
45    
46    
47  void Paso_Pattern_coup(Paso_SparseMatrix* A, index_t* mis_marker) {  void Paso_Pattern_coup(Paso_SparseMatrix* A, index_t* mis_marker) {
48    
49    index_t index_offset=(A->pattern->type & PATTERN_FORMAT_OFFSET1 ? 1:0);    index_t index_offset=(A->pattern->type & PATTERN_FORMAT_OFFSET1 ? 1:0);
# Line 91  void Paso_Pattern_coup(Paso_SparseMatrix Line 93  void Paso_Pattern_coup(Paso_SparseMatrix
93       /* swap to TRUE/FALSE in mis_marker */       /* swap to TRUE/FALSE in mis_marker */
94       #pragma omp parallel for private(i) schedule(static)       #pragma omp parallel for private(i) schedule(static)
95       for (i=0;i<n;i++) mis_marker[i]=(mis_marker[i]==IS_IN_MIS);       for (i=0;i<n;i++) mis_marker[i]=(mis_marker[i]==IS_IN_MIS);
96    }
97    
98    /*
99     *
100     * Return a strength of connection mask using the classical
101     * strength of connection measure by Ruge and Stuben.
102     *
103     * Specifically, an off-diagonal entry A[i.j] is a strong
104     * connection if:
105     *  
106     *      -A[i,j] >= theta * max( -A[i,k] )   where k != i
107     *
108     * Otherwise, the connection is weak.
109     *  
110     */  
111    void Paso_Pattern_RS(Paso_SparseMatrix* A, index_t* mis_marker, double theta)
112    {
113      index_t index_offset=(A->pattern->type & PATTERN_FORMAT_OFFSET1 ? 1:0);
114      dim_t i,j;
115      index_t naib,iptr;
116      double threshold,min_offdiagonal;
117      bool_t flag;
118      dim_t n=A->pattern->numOutput;
119      if (A->pattern->type & PATTERN_FORMAT_SYM) {
120        Paso_setError(TYPE_ERROR,"Paso_Pattern_mis: symmetric matrix pattern is not supported yet");
121        return;
122      }
123    
124    /* is there any vertex available ?*/
125         while (Paso_Util_isAny(n,mis_marker,IS_AVAILABLE)) {
126    
127          #pragma omp parallel for private(i,iptr,min_offdiagonal,threshold) schedule(static)
128          for (i=0;i<n;++i) {
129            min_offdiagonal = A->val[A->pattern->ptr[i]];
130            for (iptr=A->pattern->ptr[i]-index_offset;iptr<A->pattern->ptr[i+1]-index_offset; ++iptr) {
131                if(A->pattern->index[iptr] != i){
132                    min_offdiagonal = MIN(min_offdiagonal,A->val[iptr]);
133                }
134            }
135    
136            threshold = theta*min_offdiagonal;
137            #pragma omp parallel for private(iptr) schedule(static)
138            for (iptr=A->pattern->ptr[i]-index_offset;iptr<A->pattern->ptr[i+1]-index_offset; ++iptr) {
139                mis_marker[i]=IS_CONNECTED_TO_MIS;
140                if(A->val[iptr] <= threshold){
141                    if(A->pattern->index[iptr] != i){
142                        mis_marker[i]=IS_IN_MIS;
143                    }
144                }
145            }
146          }
147        }
148         /* swap to TRUE/FALSE in mis_marker */
149         #pragma omp parallel for private(i) schedule(static)
150         for (i=0;i<n;i++) mis_marker[i]=(mis_marker[i]==IS_IN_MIS);
151  }  }
152  #undef IS_AVAILABLE  #undef IS_AVAILABLE
153  #undef IS_IN_MIS_NOW  #undef IS_IN_MIS_NOW
Legend:
Removed from v.1883  
changed lines
  Added in v.1890
  ViewVC Help
Powered by ViewVC 1.1.26