Addede more block interpolation routines.

This commit is contained in:
baker59 2007-02-15 19:30:14 +00:00
parent beab40f806
commit 156a6e2721
5 changed files with 4085 additions and 396 deletions

View File

@ -468,7 +468,14 @@ hypre_CSRBlockMatrixBlockNorm(int norm_type, double* data, double* out, int bloc
switch (norm_type)
{
case 6: /* sum of all elements in the block */
{
for(i = 0; i < sz; i++)
{
sum += (data[i]);
}
break;
}
case 5: /* one norm - max col sum*/
{
@ -657,7 +664,116 @@ hypre_CSRBlockMatrixBlockMultAddDiag(double* i1, double* i2, double beta,
}
return 0;
}
/*--------------------------------------------------------------------------
* hypre_CSRBlockMatrixBlockMultAddDiag2 (scales cols of il by diag of i2)
* ((o) = (i1) * diag(i2) + beta * (o))
*--------------------------------------------------------------------------*/
int
hypre_CSRBlockMatrixBlockMultAddDiag2(double* i1, double* i2, double beta,
double* o, int block_size)
{
int i, j;
if (beta == 0.0)
{
for (i = 0; i < block_size; i++)
{
for (j = 0; j < block_size; j++)
{
o[i*block_size + j] = i1[i*block_size + j] * i2[j*block_size + j];
}
}
}
else if (beta == 1.0)
{
for (i = 0; i < block_size; i++)
{
for (j = 0; j < block_size; j++)
{
o[i*block_size + j] = o[i*block_size + j] + i1[i*block_size + j] * i2[j*block_size + j];
}
}
}
else
{
for (i = 0; i < block_size; i++)
{
for (j = 0; j < block_size; j++)
{
o[i*block_size + j] = beta * o[i*block_size + j] + i1[i*block_size + j] * i2[j*block_size + j];
}
}
}
return 0;
}
/*--------------------------------------------------------------------------
* hypre_CSRBlockMatrixBlockMultAddDiag3 (scales cols of il by i2 -
whose diag elements are row sums)
* ((o) = (i1) * diag(i2) + beta * (o))
*--------------------------------------------------------------------------*/
int
hypre_CSRBlockMatrixBlockMultAddDiag3(double* i1, double* i2, double beta,
double* o, int block_size)
{
int i, j;
double *row_sum;
row_sum = hypre_CTAlloc(double, block_size);
for (i = 0; i < block_size; i++)
{
for (j = 0; j < block_size; j++)
{
row_sum[i] += i2[i*block_size + j];
}
}
if (beta == 0.0)
{
for (i = 0; i < block_size; i++)
{
for (j = 0; j < block_size; j++)
{
o[i*block_size + j] = i1[i*block_size + j] * row_sum[j];
}
}
}
else if (beta == 1.0)
{
for (i = 0; i < block_size; i++)
{
for (j = 0; j < block_size; j++)
{
o[i*block_size + j] = o[i*block_size + j] + i1[i*block_size + j] * row_sum[j];
}
}
}
else
{
for (i = 0; i < block_size; i++)
{
for (j = 0; j < block_size; j++)
{
o[i*block_size + j] = beta * o[i*block_size + j] + i1[i*block_size + j] * row_sum[j];
}
}
}
hypre_TFree(row_sum);
return 0;
}
/*--------------------------------------------------------------------------
* hypre_CSRBlockMatrixBlockMatvec
* (ov = alpha* mat * v + beta * ov)
@ -1276,6 +1392,86 @@ hypre_CSRBlockMatrixBlockInvMultDiag(double* i1, double* i2, double* o, int bloc
return (ierr);
}
/*--------------------------------------------------------------------------
* hypre_CSRBlockMatrixBlockInvMultDiag2
* (o = (i1)* diag(i2)^-1) - so this scales the cols of il by
the diag entries in i2
*--------------------------------------------------------------------------*/
int
hypre_CSRBlockMatrixBlockInvMultDiag2(double* i1, double* i2, double* o, int block_size)
{
int ierr = 0;
int i, j;
double eps, tmp;
eps = 1.0e-8;
for (i = 0; i < block_size; i++)
{
if (fabs(i2[i*block_size + i]) > eps)
{
tmp = 1 / i2[i*block_size + i];
}
else
{
tmp = 1.0;
}
for (j=0; j< block_size; j++) /* this should be re-written to access by row (not col)! */
{
o[j*block_size + i] = i1[j*block_size + i] *tmp;
}
}
return (ierr);
}
/*--------------------------------------------------------------------------
* hypre_CSRBlockMatrixBlockInvMultDiag3
* (o = (i1)* diag(i2)^-1) - so this scales the cols of il by
the i2 whose diags are row sums
*--------------------------------------------------------------------------*/
int
hypre_CSRBlockMatrixBlockInvMultDiag3(double* i1, double* i2, double* o, int block_size)
{
int ierr = 0;
int i, j;
double eps, tmp, row_sum;
eps = 1.0e-8;
for (i = 0; i < block_size; i++)
{
/* get row sum of i2, row i */
row_sum = 0.0;
for (j=0; j< block_size; j++)
{
row_sum += i2[i*block_size + j];
}
/* invert */
if (fabs(row_sum) > eps)
{
tmp = 1 / row_sum;
}
else
{
tmp = 1.0;
}
/* scale col of i1 */
for (j=0; j< block_size; j++) /* this should be re-written to access by row (not col)! */
{
o[j*block_size + i] = i1[j*block_size + i] *tmp;
}
}
return (ierr);
}
/*--------------------------------------------------------------------------

View File

@ -99,10 +99,26 @@ int hypre_CSRBlockMatrixBlockAdd(double *, double *, double*, int);
int hypre_CSRBlockMatrixBlockMultAdd(double *, double *, double, double *, int);
int hypre_CSRBlockMatrixBlockMultAddDiag(double *, double *, double, double *, int);
int
hypre_CSRBlockMatrixBlockMultAddDiag2(double* i1, double* i2, double beta,
double* o, int block_size);
int
hypre_CSRBlockMatrixBlockMultAddDiag3(double* i1, double* i2, double beta,
double* o, int block_size);
int hypre_CSRBlockMatrixBlockInvMult(double *, double *, double *, int);
int hypre_CSRBlockMatrixBlockInvMultDiag(double *, double *, double *, int);
int
hypre_CSRBlockMatrixBlockInvMultDiag2(double* i1, double* i2, double* o, int block_size);
int
hypre_CSRBlockMatrixBlockInvMultDiag3(double* i1, double* i2, double* o, int block_size);
int hypre_CSRBlockMatrixBlockMultInv(double *, double *, double *, int);
int hypre_CSRBlockMatrixBlockTranspose(double *, double *, int);

View File

@ -40,12 +40,13 @@
3 = largest element (positive or negative)
4 = 1-norm
5 = inf - norm
6 = sum of all elements
*--------------------------------------------------------------------------*/
int
hypre_BoomerAMGBlockCreateNodalA(hypre_ParCSRBlockMatrix *A,
int option,
hypre_ParCSRMatrix **AN_ptr)
int option, int diag_option,
hypre_ParCSRMatrix **AN_ptr)
{
MPI_Comm comm = hypre_ParCSRBlockMatrixComm(A);
hypre_CSRBlockMatrix *A_diag = hypre_ParCSRBlockMatrixDiag(A);
@ -107,7 +108,10 @@ hypre_BoomerAMGBlockCreateNodalA(hypre_ParCSRBlockMatrix *A,
int global_num_nodes;
int num_nodes;
int index, k;
double tmp;
double sum;
@ -169,15 +173,40 @@ hypre_BoomerAMGBlockCreateNodalA(hypre_ParCSRBlockMatrix *A,
AN_diag_data[i] = tmp;
}
#if 0
/* to compare with serial - make diag entries negative*/
for (i=0; i < num_nodes; i++)
if (diag_option ==1 )
{
AN_diag_data[AN_diag_i[i]] = - AN_diag_data[AN_diag_i[i]];
/* make the diag entry the negative of the sum of off-diag entries (NEED to get more below!)*/
/* the diagonal is the first element listed in each row - */
for (i=0; i < num_nodes; i++)
{
index = AN_diag_i[i];
sum = 0.0;
for (k = AN_diag_i[i]+1; k < AN_diag_i[i+1]; k++)
{
sum += AN_diag_data[k];
}
AN_diag_data[index] = -sum;
}
}
#endif
else if (diag_option == 2)
{
/* make all diagonal entries negative */
/* the diagonal is the first element listed in each row - */
for (i=0; i < num_nodes; i++)
{
index = AN_diag_i[i];
AN_diag_data[index] = -AN_diag_data[index];
}
}
/* copy the commpkg */
if (comm_pkg)
@ -266,6 +295,29 @@ hypre_BoomerAMGBlockCreateNodalA(hypre_ParCSRBlockMatrix *A,
hypre_CSRMatrixJ(AN_offd) = AN_offd_j;
hypre_CSRMatrixData(AN_offd) = AN_offd_data;
if (diag_option ==1 )
{
/* make the diag entry the negative of the sum of off-diag entries (here we are adding the
off_diag contribution)*/
/* the diagonal is the first element listed in each row of AN_diag_data - */
for (i=0; i < num_nodes; i++)
{
sum = 0.0;
for (k = AN_offd_i[i]; k < AN_offd_i[i+1]; k++)
{
sum += AN_offd_data[k];
}
index = AN_diag_i[i];/* location of diag entry in data */
AN_diag_data[index] -= sum; /* subtract from current value */
}
}
/* now create AN */
AN = hypre_ParCSRMatrixCreate(comm, global_num_nodes, global_num_nodes,

File diff suppressed because it is too large Load Diff

View File

@ -160,9 +160,35 @@ hypre_BoomerAMGBuildBlockInterp( hypre_ParCSRBlockMatrix *A,
int *dof_func,
int debug_flag,
double trunc_factor,
int add_weak_to_diag,
int *col_offd_S_to_A,
hypre_ParCSRBlockMatrix **P_ptr);
int
hypre_BoomerAMGBuildBlockInterpRV( hypre_ParCSRBlockMatrix *A,
int *CF_marker,
hypre_ParCSRMatrix *S,
int *num_cpts_global,
int num_functions,
int *dof_func,
int debug_flag,
double trunc_factor,
int *col_offd_S_to_A,
hypre_ParCSRBlockMatrix **P_ptr);
int
hypre_BoomerAMGBuildBlockInterpRV2( hypre_ParCSRBlockMatrix *A,
int *CF_marker,
hypre_ParCSRMatrix *S,
int *num_cpts_global,
int num_functions,
int *dof_func,
int debug_flag,
double trunc_factor,
int *col_offd_S_to_A,
hypre_ParCSRBlockMatrix **P_ptr);
int
hypre_BoomerAMGBuildBlockInterpDiag( hypre_ParCSRBlockMatrix *A,
int *CF_marker,
@ -172,6 +198,7 @@ hypre_BoomerAMGBuildBlockInterpDiag( hypre_ParCSRBlockMatrix *A,
int *dof_func,
int debug_flag,
double trunc_factor,
int add_weak_to_diag,
int *col_offd_S_to_A,
hypre_ParCSRBlockMatrix **P_ptr);
@ -253,7 +280,7 @@ hypre_ParCSRBlockMatrixSetDNumNonzeros( hypre_ParCSRBlockMatrix *matrix);
int
hypre_BoomerAMGBlockCreateNodalA(hypre_ParCSRBlockMatrix *A,
int option,
int option, int diag_option,
hypre_ParCSRMatrix **AN_ptr);
hypre_ParVector *