113 int i, j, len, count, col, idx = 0;
114 int *list, *marker, *fill, *tmpFill;
115 int temp, m, from = ctx->from, to = ctx->to;
116 int *n2o_row, *o2n_col, beg_row, beg_rowP;
140 n2o_row = ctx->sg->n2o_row;
141 o2n_col = ctx->sg->o2n_col;
142 beg_row = ctx->sg->beg_row[
myid_dh];
143 beg_rowP = ctx->sg->beg_rowP[
myid_dh];
146 list = (
int *)
MALLOC_DH ((m + 1) *
sizeof (int));
148 marker = (
int *)
MALLOC_DH (m *
sizeof (
int));
150 tmpFill = (
int *)
MALLOC_DH (m *
sizeof (
int));
152 for (i = 0; i < m; ++i)
156 for (i = 0; i < m; ++i)
165 for (i = from; i < to; ++i)
167 int row = n2o_row[i];
168 int globalRow = row + beg_row;
176 "ILU_seq ================================= starting local row: %i, (global= %i) level= %i\n",
194 len, CVAL, AVAL, o2n_col, ctx, debug);
198 if (idx + count > F->
alloc)
202 SET_INFO (
"REALLOCATED from ilu_seq");
213 fill[idx] = tmpFill[col];
225 while (cval[temp] != i)
242 fprintf (
logFile,
"ILU_seq: ");
243 for (j = rp[i]; j < rp[i + 1]; ++j)
248 fprintf (
logFile,
"%i,%i,%g ; ", 1 + cval[j], fill[j], aval[j]);
255 for (j = rp[i]; j < rp[i + 1]; ++j)
266 sprintf (
msgBuf_dh,
"zero diagonal in local row %i", i + 1);
281 int start = rp[from];
283 for (i = start; i < stop; ++i)
290 for (i = to + 1; i < m; ++i)
303 int h, i, j, len, count, col, idx = 0;
304 int *list, *marker, *fill, *tmpFill;
306 int *n2o_row, *o2n_col, *beg_rowP, *n2o_sub, blocks;
307 int *row_count, *dummy = NULL, dummy2[1];
312 bool bj =
false, constrained =
false;
323 if (!strcmp (ctx->algo_par,
"bj"))
348 dummy = (
int *)
MALLOC_DH (m *
sizeof (
int));
350 for (i = 0; i < m; ++i)
363 list = (
int *)
MALLOC_DH ((m + 1) *
sizeof (int));
365 marker = (
int *)
MALLOC_DH (m *
sizeof (
int));
367 tmpFill = (
int *)
MALLOC_DH (m *
sizeof (
int));
369 for (i = 0; i < m; ++i)
373 for (i = 0; i < m; ++i)
378 for (h = 0; h < blocks; ++h)
381 int curBlock = n2o_sub[h];
382 int first_row = beg_rowP[curBlock];
383 int end_row = first_row + row_count[curBlock];
387 fprintf (
logFile,
"\n\nILU_seq BLOCK: %i @@@@@@@@@@@@@@@ \n",
391 for (i = first_row; i < end_row; ++i)
393 int row = n2o_row[i];
399 "ILU_seq global: %i local: %i =================================\n",
400 1 + gr, 1 + i - first_row);
420 len, CVAL, AVAL, o2n_col, ctx, debug);
424 if (idx + count > F->
alloc)
428 SET_INFO (
"REALLOCATED from ilu_seq");
440 if (constrained && !bj)
442 if (col >= first_row && col < end_row)
445 fill[idx] = tmpFill[col];
453 fill[idx] = tmpFill[col];
467 if (col >= first_row && col < end_row)
470 fill[idx] = tmpFill[col];
484 fill[idx] = tmpFill[col];
495 while (cval[temp] != i)
509 fprintf (
logFile,
"ILU_seq: ");
510 for (j = rp[i]; j < rp[i + 1]; ++j)
515 fprintf (
logFile,
"%i,%i,%g ; ", 1 + cval[j], fill[j],
524 for (j = rp[i]; j < rp[i + 1]; ++j)
535 sprintf (
msgBuf_dh,
"zero diagonal in local row %i", i + 1);
568 int *list,
int *marker,
int *tmpFill,
569 int len,
int *CVAL,
double *AVAL,
573 int *cval = ctx->F->cval, *diag = ctx->F->diag, *rp = ctx->F->rp;
574 int *fill = ctx->F->fill;
576 int j, node, tmp, col, head;
577 int fill1, fill2, beg_row;
579 double thresh = ctx->sparseTolA;
582 scale = ctx->scale[localRow];
584 beg_row = ctx->sg->beg_row[
myid_dh];
591 for (j = 0; j < len; ++j)
600 if (fabs (val) > thresh || col == localRow)
603 while (col > list[tmp])
605 list[col] = list[tmp];
608 marker[col] = localRow;
613 if (marker[localRow] != localRow)
616 while (localRow > list[tmp])
618 list[localRow] = list[tmp];
619 list[tmp] = localRow;
620 tmpFill[localRow] = 0;
621 marker[localRow] = localRow;
630 while (list[head] < localRow)
633 fill1 = tmpFill[node];
637 fprintf (
logFile,
"ILU_seq sf updating from row: %i\n",
643 for (j = diag[node] + 1; j < rp[node + 1]; ++j)
646 fill2 = fill1 + fill[j] + 1;
653 if (marker[col] < localRow)
656 marker[col] = localRow;
657 tmpFill[col] = fill2;
658 while (col > list[tmp])
660 list[col] = list[tmp];
669 (fill2 < tmpFill[col]) ? fill2 : tmpFill[col];
684 int len,
int *CVAL,
double *AVAL,
689 int *rp = ctx->F->rp, *cval = ctx->F->cval;
690 int *diag = ctx->F->diag;
693 REAL_DH *aval = ctx->F->aval, scale;
695 scale = ctx->scale[localRow];
696 beg_row = ctx->sg->beg_row[
myid_dh];
700 for (j = rp[localRow]; j < rp[localRow + 1]; ++j)
708 for (j = 0; j < len; ++j)
714 work[col] = val * scale;
723 for (j = rp[localRow]; j < diag[localRow]; ++j)
729 pv = aval[diag[row]];
737 if (pc != 0.0 && pv != 0.0)
739 multiplier = pc / pv;
740 work[row] = multiplier;
745 "ILU_seq nf updating from row: %i; multiplier= %g\n",
746 1 + row, multiplier);
749 for (k = diag[row] + 1; k < rp[row + 1]; ++k)
752 work[col] -= (multiplier * aval[k]);
760 "ILU_seq nf NO UPDATE from row %i; pc = %g; pv = %g\n",
768 if (fabs (work[i]) <= pivotTol)
771 aval[diag[i]] = pivotFix;
791 int i, len, count, col, idx = 0;
793 int temp, m, from, to;
794 int *n2o_row, *o2n_col, beg_row, beg_rowP;
795 double *AVAL, droptol;
813 droptol = ctx->droptol;
820 n2o_row = ctx->sg->n2o_row;
821 o2n_col = ctx->sg->o2n_col;
822 beg_row = ctx->sg->beg_row[
myid_dh];
823 beg_rowP = ctx->sg->beg_rowP[
myid_dh];
827 list = (
int *)
MALLOC_DH ((m + 1) *
sizeof (int));
829 marker = (
int *)
MALLOC_DH (m *
sizeof (
int));
831 for (i = 0; i < m; ++i)
836 for (i = 0; i < m; ++i)
840 for (i = from; i < to; ++i)
842 int row = n2o_row[i];
843 int globalRow = row + beg_row;
853 len, CVAL, AVAL, work, ctx, debug);
860 if (idx + count > F->
alloc)
864 SET_INFO (
"REALLOCATED from ilu_seq");
877 if (col == i || fabs (val) > droptol)
891 while (cval[temp] != i)
898 sprintf (
msgBuf_dh,
"zero diagonal in local row %i", i + 1);
906 int start = rp[from];
908 for (i = start; i < stop; ++i)
921 int len,
int *CVAL,
double *AVAL,
925 int j, col, m = ctx->
m, *rp = F->rp, *cval = F->cval;
926 int tmp, *diag = F->diag;
928 int count = 0, beg_row;
930 double mult, *aval = F->aval;
931 double scale, pv, pc;
932 double droptol = ctx->droptol;
933 double thresh = ctx->sparseTolA;
935 scale = ctx->scale[localRow];
937 beg_row = ctx->sg->beg_row[
myid_dh];
945 for (j = 0; j < len; ++j)
954 if (fabs (val) > thresh || col == localRow)
957 while (col > list[tmp])
959 list[col] = list[tmp];
962 marker[col] = localRow;
967 if (marker[localRow] != localRow)
970 while (localRow > list[tmp])
972 list[localRow] = list[tmp];
973 list[tmp] = localRow;
974 marker[localRow] = localRow;
980 while (list[head] < localRow)
982 int row = list[head];
988 pv = aval[diag[row]];
992 if (fabs (mult) > droptol)
996 for (j = diag[row] + 1; j < rp[row + 1]; ++j)
999 work[col] -= (mult * aval[j]);
1002 if (marker[col] < localRow)
1004 marker[col] = localRow;
1006 while (col > list[tmp])
1008 list[col] = list[tmp];
int SubdomainGraph_dhFindOwner(SubdomainGraph_dh s, int idx, bool permuted)
static int symbolic_row_private(int localRow, int *list, int *marker, int *tmpFill, int len, int *CVAL, double *AVAL, int *o2n_col, Euclid_dh ctx, bool debug)
static int numeric_row_private(int localRow, int len, int *CVAL, double *AVAL, REAL_DH *work, int *o2n_col, Euclid_dh ctx, bool debug)
int ilut_row_private(int localRow, int *list, int *o2n_col, int *marker, int len, int *CVAL, double *AVAL, REAL_DH *work, Euclid_dh ctx, bool debug)