IFPACK  Development
 All Classes Namespaces Files Functions Variables Enumerations Friends Pages
SortedList_dh.c
1 /*@HEADER
2 // ***********************************************************************
3 //
4 // Ifpack: Object-Oriented Algebraic Preconditioner Package
5 // Copyright (2002) Sandia Corporation
6 //
7 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8 // license for use of this work by or on behalf of the U.S. Government.
9 //
10 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
38 //
39 // ***********************************************************************
40 //@HEADER
41 */
42 
43 #include "SortedList_dh.h"
44 #include "Mem_dh.h"
45 #include "Parser_dh.h"
46 #include "Hash_i_dh.h"
47 #include "SubdomainGraph_dh.h"
48 
49 
51 {
52  int m; /* number of local rows */
53  int row; /* local number of row being factored */
54  int beg_row; /* global number of first locally owned row, wrt A */
55  int beg_rowP; /* global number of first locally owned row, wrt F */
56  int count; /* number of items entered in the list,
57  plus 1 (for header node)
58  */
59  int countMax; /* same as count, but includes number of items that my have
60  been deleted from calling SortedList_dhEnforceConstraint()
61  */
62  int *o2n_local; /* not owned! */
63  Hash_i_dh o2n_external; /* not owned! */
64 
65  SRecord *list; /* the sorted list */
66  int alloc; /* allocated length of list */
67  int getLower; /* index used for returning lower tri elts */
68  int get; /* index of returning all elts; */
69 
70  bool debug;
71 };
72 
73 static void lengthen_list_private (SortedList_dh sList);
74 
75 
76 #undef __FUNC__
77 #define __FUNC__ "SortedList_dhCreate"
78 void
79 SortedList_dhCreate (SortedList_dh * sList)
80 {
81  START_FUNC_DH
82  struct _sortedList_dh *tmp =
83  (struct _sortedList_dh *) MALLOC_DH (sizeof (struct _sortedList_dh));
84  CHECK_V_ERROR;
85  *sList = tmp;
86  tmp->m = 0;
87  tmp->row = -1;
88  tmp->beg_row = 0;
89  tmp->count = 1;
90  tmp->countMax = 1;
91  tmp->o2n_external = NULL;
92  tmp->o2n_local = NULL;
93 
94  tmp->get = 0;
95  tmp->getLower = 0;
96  tmp->alloc = 0;
97  tmp->list = NULL;
98  tmp->debug = Parser_dhHasSwitch (parser_dh, "-debug_SortedList");
99 END_FUNC_DH}
100 
101 
102 #undef __FUNC__
103 #define __FUNC__ "SortedList_dhDestroy"
104 void
105 SortedList_dhDestroy (SortedList_dh sList)
106 {
107  START_FUNC_DH if (sList->list != NULL)
108  {
109  FREE_DH (sList->list);
110  CHECK_V_ERROR;
111  }
112  FREE_DH (sList);
113  CHECK_V_ERROR;
114 END_FUNC_DH}
115 
116 
117 #undef __FUNC__
118 #define __FUNC__ "SortedList_dhInit"
119 void
120 SortedList_dhInit (SortedList_dh sList, SubdomainGraph_dh sg)
121 {
122  START_FUNC_DH sList->o2n_local = sg->o2n_col;
123  sList->m = sg->m;
124  sList->beg_row = sg->beg_row[myid_dh];
125  sList->beg_rowP = sg->beg_rowP[myid_dh];
126  sList->count = 1; /* "1" is for the header node */
127  sList->countMax = 1; /* "1" is for the header node */
128  sList->o2n_external = sg->o2n_ext;
129 
130  /* heuristic: "m" should be a good number of nodes */
131  sList->alloc = sList->m + 5;
132  sList->list = (SRecord *) MALLOC_DH (sList->alloc * sizeof (SRecord));
133  sList->list[0].col = INT_MAX;
134  sList->list[0].next = 0;
135 END_FUNC_DH}
136 
137 
138 #undef __FUNC__
139 #define __FUNC__ "SortedList_dhReset"
140 void
141 SortedList_dhReset (SortedList_dh sList, int row)
142 {
143  START_FUNC_DH sList->row = row;
144  sList->count = 1;
145  sList->countMax = 1;
146  sList->get = 0;
147  sList->getLower = 0;
148  sList->list[0].next = 0;
149 END_FUNC_DH}
150 
151 
152 #undef __FUNC__
153 #define __FUNC__ "SortedList_dhReadCount"
154 int
155 SortedList_dhReadCount (SortedList_dh sList)
156 {
157 START_FUNC_DH END_FUNC_VAL (sList->count - 1)}
158 
159 #undef __FUNC__
160 #define __FUNC__ "SortedList_dhResetGetSmallest"
161 void
162 SortedList_dhResetGetSmallest (SortedList_dh sList)
163 {
164  START_FUNC_DH sList->getLower = 0;
165  sList->get = 0;
166 END_FUNC_DH}
167 
168 #undef __FUNC__
169 #define __FUNC__ "SortedList_dhGetSmallest"
170 SRecord *
171 SortedList_dhGetSmallest (SortedList_dh sList)
172 {
173  START_FUNC_DH SRecord * node = NULL;
174  SRecord *list = sList->list;
175  int get = sList->get;
176 
177  get = list[get].next;
178 
179  if (list[get].col < INT_MAX)
180  {
181  node = &(list[get]);
182  sList->get = get;
183  }
184 END_FUNC_VAL (node)}
185 
186 #undef __FUNC__
187 #define __FUNC__ "SortedList_dhGetSmallestLowerTri"
188 SRecord *
189 SortedList_dhGetSmallestLowerTri (SortedList_dh sList)
190 {
191  START_FUNC_DH SRecord * node = NULL;
192  SRecord *list = sList->list;
193  int getLower = sList->getLower;
194  int globalRow = sList->row + sList->beg_rowP;
195 
196  getLower = list[getLower].next;
197 
198  if (list[getLower].col < globalRow)
199  {
200  node = &(list[getLower]);
201  sList->getLower = getLower;
202  }
203 END_FUNC_VAL (node)}
204 
205 
206 #undef __FUNC__
207 #define __FUNC__ "SortedList_dhPermuteAndInsert"
208 bool
209 SortedList_dhPermuteAndInsert (SortedList_dh sList, SRecord * sr,
210  double thresh)
211 {
212  START_FUNC_DH bool wasInserted = false;
213  int col = sr->col;
214  double testVal = fabs (sr->val);
215  int beg_row = sList->beg_row, end_row = beg_row + sList->m;
216  int beg_rowP = sList->beg_rowP;
217 
218  /* insertion of local indices */
219  if (col >= beg_row && col < end_row)
220  {
221  /* convert to local indexing and permute */
222  col -= beg_row;
223  col = sList->o2n_local[col];
224 
225  /* sparsification */
226  if (testVal > thresh || col == sList->row)
227  {
228  col += beg_rowP;
229  }
230  else
231  {
232  col = -1;
233 /*
234 fprintf(logFile, "local row: %i DROPPED: col= %i val= %g (thresh= %g)\n",
235  sList->row+1, sr->col+1, testVal, thresh);
236 */
237  }
238  }
239 
240 
241  /* insertion of external indices */
242  else
243  {
244  /* sparsification for external indices */
245  if (testVal < thresh)
246  goto END_OF_FUNCTION;
247 
248  /* permute column index */
249  if (sList->o2n_external == NULL)
250  {
251  col = -1;
252  }
253  else
254  {
255  int tmp = Hash_i_dhLookup (sList->o2n_external, col);
256  CHECK_ERROR (-1);
257  if (tmp == -1)
258  {
259  col = -1;
260  }
261  else
262  {
263  col = tmp;
264  }
265  }
266  }
267 
268  if (col != -1)
269  {
270  sr->col = col;
271  SortedList_dhInsert (sList, sr);
272  CHECK_ERROR (-1);
273  wasInserted = true;
274  }
275 
276 END_OF_FUNCTION:;
277 
278 END_FUNC_VAL (wasInserted)}
279 
280 
281 #undef __FUNC__
282 #define __FUNC__ "SortedList_dhInsertOrUpdate"
283 void
284 SortedList_dhInsertOrUpdate (SortedList_dh sList, SRecord * sr)
285 {
286  START_FUNC_DH SRecord * node = SortedList_dhFind (sList, sr);
287  CHECK_V_ERROR;
288 
289  if (node == NULL)
290  {
291  SortedList_dhInsert (sList, sr);
292  CHECK_V_ERROR;
293  }
294  else
295  {
296  node->level = MIN (sr->level, node->level);
297  }
298 END_FUNC_DH}
299 
300 
301 /* note: this does NOT check to see if item was already inserted! */
302 #undef __FUNC__
303 #define __FUNC__ "SortedList_dhInsert"
304 void
305 SortedList_dhInsert (SortedList_dh sList, SRecord * sr)
306 {
307  START_FUNC_DH int prev, next;
308  int ct, col = sr->col;
309  SRecord *list = sList->list;
310 
311  /* lengthen list if out of space */
312  if (sList->countMax == sList->alloc)
313  {
314  lengthen_list_private (sList);
315  CHECK_V_ERROR;
316  list = sList->list;
317  }
318 
319  /* add new node to end of list */
320  ct = sList->countMax;
321  sList->countMax += 1;
322  sList->count += 1;
323 
324  list[ct].col = col;
325  list[ct].level = sr->level;
326  list[ct].val = sr->val;
327 
328  /* splice new node into list */
329  prev = 0;
330  next = list[0].next;
331  while (col > list[next].col)
332  {
333  prev = next;
334  next = list[next].next;
335  }
336  list[prev].next = ct;
337  list[ct].next = next;
338 END_FUNC_DH}
339 
340 
341 #undef __FUNC__
342 #define __FUNC__ "SortedList_dhFind"
343 SRecord *
344 SortedList_dhFind (SortedList_dh sList, SRecord * sr)
345 {
346  START_FUNC_DH int i, count = sList->countMax;
347  int c = sr->col;
348  SRecord *s = sList->list;
349  SRecord *node = NULL;
350 
351  /* no need to traverse list in sorted order */
352  for (i = 1; i < count; ++i)
353  { /* start at i=1, since i=0 would be header node */
354 
355  if (s[i].col == c)
356  {
357  node = &(s[i]);
358  break;
359  }
360  }
361 
362 END_FUNC_VAL (node)}
363 
364 #undef __FUNC__
365 #define __FUNC__ "lengthen_list_private"
366 void
367 lengthen_list_private (SortedList_dh sList)
368 {
369  START_FUNC_DH SRecord * tmp = sList->list;
370  int size = sList->alloc = 2 * sList->alloc;
371 
372  SET_INFO ("lengthening list");
373 
374  sList->list = (SRecord *) MALLOC_DH (size * sizeof (SRecord));
375  memcpy (sList->list, tmp, sList->countMax * sizeof (SRecord));
376  SET_INFO ("doubling size of sList->list");
377  FREE_DH (tmp);
378  CHECK_V_ERROR;
379 END_FUNC_DH}
380 
381 
382 /*=====================================================================
383  * functions for enforcing subdomain constraint
384  *=====================================================================*/
385 
386 
387 static bool check_constraint_private (SubdomainGraph_dh sg,
388  int thisSubdomain, int col);
389 void delete_private (SortedList_dh sList, int col);
390 
391 #undef __FUNC__
392 #define __FUNC__ "SortedList_dhEnforceConstraint"
393 void
394 SortedList_dhEnforceConstraint (SortedList_dh sList, SubdomainGraph_dh sg)
395 {
396  START_FUNC_DH int thisSubdomain = myid_dh;
397  int col, count;
398  int beg_rowP = sList->beg_rowP;
399  int end_rowP = beg_rowP + sList->m;
400  bool debug = false;
401 
402  if (Parser_dhHasSwitch (parser_dh, "-debug_SortedList"))
403  debug = true;
404 
405  if (debug)
406  {
407  fprintf (logFile, "SLIST ======= enforcing constraint for row= %i\n",
408  1 + sList->row);
409 
410  fprintf (logFile, "\nSLIST ---- before checking: ");
411  count = SortedList_dhReadCount (sList);
412  CHECK_V_ERROR;
413  while (count--)
414  {
415  SRecord *sr = SortedList_dhGetSmallest (sList);
416  CHECK_V_ERROR;
417  fprintf (logFile, "%i ", sr->col + 1);
418  }
419  fprintf (logFile, "\n");
420  sList->get = 0;
421  }
422 
423  /* for each column index in the list */
424  count = SortedList_dhReadCount (sList);
425  CHECK_V_ERROR;
426 
427  while (count--)
428  {
429  SRecord *sr = SortedList_dhGetSmallest (sList);
430  CHECK_V_ERROR;
431  col = sr->col;
432 
433  if (debug)
434  {
435  fprintf (logFile, "SLIST next col= %i\n", col + 1);
436  }
437 
438 
439  /* if corresponding row is nonlocal */
440  if (col < beg_rowP || col >= end_rowP)
441  {
442 
443  if (debug)
444  {
445  fprintf (logFile, "SLIST external col: %i ; ", 1 + col);
446  }
447 
448  /* if entry would violate subdomain constraint, discard it
449  (snip it out of the list)
450  */
451  if (check_constraint_private (sg, thisSubdomain, col))
452  {
453  delete_private (sList, col);
454  CHECK_V_ERROR;
455  sList->count -= 1;
456 
457  if (debug)
458  {
459  fprintf (logFile, " deleted\n");
460  }
461  }
462  else
463  {
464  if (debug)
465  {
466  fprintf (logFile, " kept\n");
467  }
468  }
469  }
470  }
471  sList->get = 0;
472 
473  if (debug)
474  {
475  fprintf (logFile, "SLIST---- after checking: ");
476  count = SortedList_dhReadCount (sList);
477  CHECK_V_ERROR;
478  while (count--)
479  {
480  SRecord *sr = SortedList_dhGetSmallest (sList);
481  CHECK_V_ERROR;
482  fprintf (logFile, "%i ", sr->col + 1);
483  }
484  fprintf (logFile, "\n");
485  fflush (logFile);
486  sList->get = 0;
487  }
488 
489 END_FUNC_DH}
490 
491 
492 /* this is similar to a function in ilu_seq.c */
493 #undef __FUNC__
494 #define __FUNC__ "check_constraint_private"
495 bool
496 check_constraint_private (SubdomainGraph_dh sg, int p1, int j)
497 {
498  START_FUNC_DH bool retval = false;
499  int i, p2;
500  int *nabors, count;
501 
502  p2 = SubdomainGraph_dhFindOwner (sg, j, true);
503 
504  nabors = sg->adj + sg->ptrs[p1];
505  count = sg->ptrs[p1 + 1] - sg->ptrs[p1];
506 
507  for (i = 0; i < count; ++i)
508  {
509  if (nabors[i] == p2)
510  {
511  retval = true;
512  break;
513  }
514  }
515 
516 END_FUNC_VAL (!retval)}
517 
518 #undef __FUNC__
519 #define __FUNC__ "delete_private"
520 void
521 delete_private (SortedList_dh sList, int col)
522 {
523  START_FUNC_DH int curNode = 0;
524  SRecord *list = sList->list;
525  int next;
526 
527  /* find node preceeding the node to be snipped out */
528  /* 'list[curNode].next' is array index of the next node in the list */
529 
530  while (list[list[curNode].next].col != col)
531  {
532  curNode = list[curNode].next;
533  }
534 
535  /* mark node to be deleted as inactive (needed for Find()) */
536  next = list[curNode].next;
537  list[next].col = -1;
538 
539  /* snip */
540  next = list[next].next;
541  list[curNode].next = next;
542 END_FUNC_DH}