diff --git a/sstruct_mv/sstruct_matrix.c b/sstruct_mv/sstruct_matrix.c index 220ec5f18..8c0a16ad7 100644 --- a/sstruct_mv/sstruct_matrix.c +++ b/sstruct_mv/sstruct_matrix.c @@ -290,6 +290,8 @@ hypre_SStructPMatrixSetValues( hypre_SStructPMatrix *pmatrix, int *smap = hypre_SStructPMatrixSMap(pmatrix, var); int *vars = hypre_SStructStencilVars(stencil); hypre_StructMatrix *smatrix; + hypre_BoxArray *grid_boxes; + hypre_Box *box; int *sentries; int i; @@ -305,13 +307,12 @@ hypre_SStructPMatrixSetValues( hypre_SStructPMatrix *pmatrix, hypre_StructMatrixSetValues(smatrix, index, nentries, sentries, values, action, -1, 0); - /* set values outside the grid in ghost zones (for AddTo and Get) */ + /* set (AddTo/Get) or clear (Set) values outside the grid in ghost zones */ if (action != 0) { + /* AddTo/Get */ hypre_SStructPGrid *pgrid = hypre_SStructPMatrixPGrid(pmatrix); hypre_Index varoffset; - hypre_BoxArray *grid_boxes; - hypre_Box *box; int done = 0; grid_boxes = hypre_StructGridBoxes(hypre_StructMatrixGrid(smatrix)); @@ -358,6 +359,25 @@ hypre_SStructPMatrixSetValues( hypre_SStructPMatrix *pmatrix, } } } + else + { + /* Set */ + grid_boxes = hypre_StructGridBoxes(hypre_StructMatrixGrid(smatrix)); + + hypre_ForBoxI(i, grid_boxes) + { + box = hypre_BoxArrayBox(grid_boxes, i); + if ((hypre_IndexX(index) < hypre_BoxIMinX(box)) || + (hypre_IndexX(index) > hypre_BoxIMaxX(box)) || + (hypre_IndexY(index) < hypre_BoxIMinY(box)) || + (hypre_IndexY(index) > hypre_BoxIMaxY(box)) || + (hypre_IndexZ(index) < hypre_BoxIMinZ(box)) || + (hypre_IndexZ(index) > hypre_BoxIMaxZ(box)) ) + { + hypre_StructMatrixClearValues(smatrix, index, nentries, sentries, i, 1); + } + } + } return hypre_error_flag; } @@ -383,10 +403,11 @@ hypre_SStructPMatrixSetBoxValues( hypre_SStructPMatrix *pmatrix, int *smap = hypre_SStructPMatrixSMap(pmatrix, var); int *vars = hypre_SStructStencilVars(stencil); hypre_StructMatrix *smatrix; + hypre_BoxArray *grid_boxes; hypre_Box *box; hypre_Box *value_box; int *sentries; - int i; + int i, j; smatrix = hypre_SStructPMatrixSMatrix(pmatrix, var, vars[entries[0]]); @@ -405,15 +426,14 @@ hypre_SStructPMatrixSetBoxValues( hypre_SStructPMatrix *pmatrix, hypre_StructMatrixSetBoxValues(smatrix, box, value_box, nentries, sentries, values, action, -1, 0); - /* set values outside the grid in ghost zones (for AddTo and Get) */ + /* set (AddTo/Get) or clear (Set) values outside the grid in ghost zones */ if (action != 0) { + /* AddTo/Get */ hypre_SStructPGrid *pgrid = hypre_SStructPMatrixPGrid(pmatrix); hypre_Index varoffset; - hypre_BoxArray *grid_boxes; hypre_BoxArray *left_boxes, *done_boxes, *temp_boxes; hypre_Box *left_box, *done_box, *int_box; - int j; hypre_SStructVariableGetOffset(hypre_SStructPGridVarType(pgrid, var), hypre_SStructPGridNDim(pgrid), varoffset); @@ -458,6 +478,30 @@ hypre_SStructPMatrixSetBoxValues( hypre_SStructPMatrix *pmatrix, hypre_BoxArrayDestroy(done_boxes); hypre_BoxArrayDestroy(temp_boxes); } + else + { + /* Set */ + hypre_BoxArray *diff_boxes; + hypre_Box *grid_box, *diff_box; + + grid_boxes = hypre_StructGridBoxes(hypre_StructMatrixGrid(smatrix)); + diff_boxes = hypre_BoxArrayCreate(0); + + hypre_ForBoxI(i, grid_boxes) + { + grid_box = hypre_BoxArrayBox(grid_boxes, i); + hypre_BoxArraySetSize(diff_boxes, 0); + hypre_SubtractBoxes(box, grid_box, diff_boxes); + + hypre_ForBoxI(j, diff_boxes) + { + diff_box = hypre_BoxArrayBox(diff_boxes, j); + hypre_StructMatrixClearBoxValues(smatrix, diff_box, nentries, sentries, + i, 1); + } + } + hypre_BoxArrayDestroy(diff_boxes); + } hypre_BoxDestroy(box); diff --git a/sstruct_mv/sstruct_vector.c b/sstruct_mv/sstruct_vector.c index ff91fd1b9..433d049ed 100644 --- a/sstruct_mv/sstruct_vector.c +++ b/sstruct_mv/sstruct_vector.c @@ -179,19 +179,20 @@ hypre_SStructPVectorSetValues( hypre_SStructPVector *pvector, int action ) { hypre_StructVector *svector = hypre_SStructPVectorSVector(pvector, var); + hypre_BoxArray *grid_boxes; + hypre_Box *box; + int i; /* set values inside the grid */ hypre_StructVectorSetValues(svector, index, value, action, -1, 0); - /* set values outside the grid in ghost zones (for AddTo and Get) */ + /* set (AddTo/Get) or clear (Set) values outside the grid in ghost zones */ if (action != 0) { - hypre_SStructPGrid *pgrid = hypre_SStructPVectorPGrid(pvector); - hypre_Index varoffset; - hypre_BoxArray *grid_boxes; - hypre_Box *box; - int i; - int done = 0; + /* AddTo/Get */ + hypre_SStructPGrid *pgrid = hypre_SStructPVectorPGrid(pvector); + hypre_Index varoffset; + int done = 0; grid_boxes = hypre_StructGridBoxes(hypre_StructVectorGrid(svector)); @@ -236,6 +237,25 @@ hypre_SStructPVectorSetValues( hypre_SStructPVector *pvector, } } } + else + { + /* Set */ + grid_boxes = hypre_StructGridBoxes(hypre_StructVectorGrid(svector)); + + hypre_ForBoxI(i, grid_boxes) + { + box = hypre_BoxArrayBox(grid_boxes, i); + if ((hypre_IndexX(index) < hypre_BoxIMinX(box)) || + (hypre_IndexX(index) > hypre_BoxIMaxX(box)) || + (hypre_IndexY(index) < hypre_BoxIMinY(box)) || + (hypre_IndexY(index) > hypre_BoxIMaxY(box)) || + (hypre_IndexZ(index) < hypre_BoxIMinZ(box)) || + (hypre_IndexZ(index) > hypre_BoxIMaxZ(box)) ) + { + hypre_StructVectorClearValues(svector, index, i, 1); + } + } + } return hypre_error_flag; } @@ -255,8 +275,10 @@ hypre_SStructPVectorSetBoxValues( hypre_SStructPVector *pvector, int action ) { hypre_StructVector *svector = hypre_SStructPVectorSVector(pvector, var); + hypre_BoxArray *grid_boxes; hypre_Box *box; hypre_Box *value_box; + int i, j; box = hypre_BoxCreate(); hypre_CopyIndex(ilower, hypre_BoxIMin(box)); @@ -266,15 +288,14 @@ hypre_SStructPVectorSetBoxValues( hypre_SStructPVector *pvector, /* set values inside the grid */ hypre_StructVectorSetBoxValues(svector, box, value_box, values, action, -1, 0); - /* set values outside the grid in ghost zones (for AddTo and Get) */ + /* set (AddTo/Get) or clear (Set) values outside the grid in ghost zones */ if (action != 0) { + /* AddTo/Get */ hypre_SStructPGrid *pgrid = hypre_SStructPVectorPGrid(pvector); hypre_Index varoffset; - hypre_BoxArray *grid_boxes; hypre_BoxArray *left_boxes, *done_boxes, *temp_boxes; hypre_Box *left_box, *done_box, *int_box; - int i, j; hypre_SStructVariableGetOffset(hypre_SStructPGridVarType(pgrid, var), hypre_SStructPGridNDim(pgrid), varoffset); @@ -318,6 +339,29 @@ hypre_SStructPVectorSetBoxValues( hypre_SStructPVector *pvector, hypre_BoxArrayDestroy(done_boxes); hypre_BoxArrayDestroy(temp_boxes); } + else + { + /* Set */ + hypre_BoxArray *diff_boxes; + hypre_Box *grid_box, *diff_box; + + grid_boxes = hypre_StructGridBoxes(hypre_StructVectorGrid(svector)); + diff_boxes = hypre_BoxArrayCreate(0); + + hypre_ForBoxI(i, grid_boxes) + { + grid_box = hypre_BoxArrayBox(grid_boxes, i); + hypre_BoxArraySetSize(diff_boxes, 0); + hypre_SubtractBoxes(box, grid_box, diff_boxes); + + hypre_ForBoxI(j, diff_boxes) + { + diff_box = hypre_BoxArrayBox(diff_boxes, j); + hypre_StructVectorClearBoxValues(svector, diff_box, i, 1); + } + } + hypre_BoxArrayDestroy(diff_boxes); + } hypre_BoxDestroy(box); diff --git a/struct_mv/_hypre_struct_mv.h b/struct_mv/_hypre_struct_mv.h index 5fb6031d2..b98df9bc5 100644 --- a/struct_mv/_hypre_struct_mv.h +++ b/struct_mv/_hypre_struct_mv.h @@ -2349,6 +2349,8 @@ int hypre_StructMatrixInitialize ( hypre_StructMatrix *matrix ); int hypre_StructMatrixSetValues ( hypre_StructMatrix *matrix , hypre_Index grid_index , int num_stencil_indices , int *stencil_indices , double *values , int action , int boxnum , int outside ); int hypre_StructMatrixSetBoxValues ( hypre_StructMatrix *matrix , hypre_Box *set_box , hypre_Box *value_box , int num_stencil_indices , int *stencil_indices , double *values , int action , int boxnum , int outside ); int hypre_StructMatrixSetConstantValues ( hypre_StructMatrix *matrix , int num_stencil_indices , int *stencil_indices , double *values , int action ); +int hypre_StructMatrixClearValues ( hypre_StructMatrix *matrix , hypre_Index grid_index , int num_stencil_indices , int *stencil_indices , int boxnum , int outside ); +int hypre_StructMatrixClearBoxValues ( hypre_StructMatrix *matrix , hypre_Box *clear_box , int num_stencil_indices , int *stencil_indices , int boxnum , int outside ); int hypre_StructMatrixAssemble ( hypre_StructMatrix *matrix ); int hypre_StructMatrixSetNumGhost ( hypre_StructMatrix *matrix , int *num_ghost ); int hypre_StructMatrixSetConstantCoefficient ( hypre_StructMatrix *matrix , int constant_coefficient ); @@ -2393,6 +2395,8 @@ int hypre_StructVectorInitializeData ( hypre_StructVector *vector , double *data int hypre_StructVectorInitialize ( hypre_StructVector *vector ); int hypre_StructVectorSetValues ( hypre_StructVector *vector , hypre_Index grid_index , double *values , int action , int boxnum , int outside ); int hypre_StructVectorSetBoxValues ( hypre_StructVector *vector , hypre_Box *set_box , hypre_Box *value_box , double *values , int action , int boxnum , int outside ); +int hypre_StructVectorClearValues ( hypre_StructVector *vector , hypre_Index grid_index , int boxnum , int outside ); +int hypre_StructVectorClearBoxValues ( hypre_StructVector *vector , hypre_Box *clear_box , int boxnum , int outside ); int hypre_StructVectorSetNumGhost ( hypre_StructVector *vector , int *num_ghost ); int hypre_StructVectorAssemble ( hypre_StructVector *vector ); int hypre_StructVectorCopy ( hypre_StructVector *x , hypre_StructVector *y ); diff --git a/struct_mv/struct_matrix.c b/struct_mv/struct_matrix.c index 066bcaf18..07d3d1c22 100644 --- a/struct_mv/struct_matrix.c +++ b/struct_mv/struct_matrix.c @@ -1057,6 +1057,177 @@ hypre_StructMatrixSetConstantValues( hypre_StructMatrix *matrix, return hypre_error_flag; } +/*-------------------------------------------------------------------------- + * (outside > 0): clear values possibly outside of the grid extents + * (outside = 0): clear values only inside the grid extents + *--------------------------------------------------------------------------*/ + +int +hypre_StructMatrixClearValues( hypre_StructMatrix *matrix, + hypre_Index grid_index, + int num_stencil_indices, + int *stencil_indices, + int boxnum, + int outside ) +{ + hypre_BoxArray *grid_boxes; + hypre_Box *grid_box; + + double *matp; + + int i, s, istart, istop; + + /*----------------------------------------------------------------------- + * Initialize some things + *-----------------------------------------------------------------------*/ + + if (outside > 0) + { + grid_boxes = hypre_StructMatrixDataSpace(matrix); + } + else + { + grid_boxes = hypre_StructGridBoxes(hypre_StructMatrixGrid(matrix)); + } + + if (boxnum < 0) + { + istart = 0; + istop = hypre_BoxArraySize(grid_boxes); + } + else + { + istart = boxnum; + istop = istart + 1; + } + + /*----------------------------------------------------------------------- + * Clear the matrix coefficients + *-----------------------------------------------------------------------*/ + + for (i = istart; i < istop; i++) + { + grid_box = hypre_BoxArrayBox(grid_boxes, i); + + if ((hypre_IndexX(grid_index) >= hypre_BoxIMinX(grid_box)) && + (hypre_IndexX(grid_index) <= hypre_BoxIMaxX(grid_box)) && + (hypre_IndexY(grid_index) >= hypre_BoxIMinY(grid_box)) && + (hypre_IndexY(grid_index) <= hypre_BoxIMaxY(grid_box)) && + (hypre_IndexZ(grid_index) >= hypre_BoxIMinZ(grid_box)) && + (hypre_IndexZ(grid_index) <= hypre_BoxIMaxZ(grid_box)) ) + { + for (s = 0; s < num_stencil_indices; s++) + { + matp = hypre_StructMatrixBoxDataValue(matrix, i, + stencil_indices[s], + grid_index); + *matp = 0.0; + } + } + } + + return hypre_error_flag; +} + +/*-------------------------------------------------------------------------- + * (outside > 0): clear values possibly outside of the grid extents + * (outside = 0): clear values only inside the grid extents + *--------------------------------------------------------------------------*/ + +int +hypre_StructMatrixClearBoxValues( hypre_StructMatrix *matrix, + hypre_Box *clear_box, + int num_stencil_indices, + int *stencil_indices, + int boxnum, + int outside ) +{ + hypre_BoxArray *grid_boxes; + hypre_Box *grid_box; + hypre_Box *int_box; + + hypre_BoxArray *data_space; + hypre_Box *data_box; + hypre_IndexRef data_start; + hypre_Index data_stride; + int datai; + double *datap; + + hypre_Index loop_size; + + int i, s, istart, istop; + int loopi, loopj, loopk; + + /*----------------------------------------------------------------------- + * Initialize some things + *-----------------------------------------------------------------------*/ + + if (outside > 0) + { + grid_boxes = hypre_StructMatrixDataSpace(matrix); + } + else + { + grid_boxes = hypre_StructGridBoxes(hypre_StructMatrixGrid(matrix)); + } + data_space = hypre_StructMatrixDataSpace(matrix); + + if (boxnum < 0) + { + istart = 0; + istop = hypre_BoxArraySize(grid_boxes); + } + else + { + istart = boxnum; + istop = istart + 1; + } + + /*----------------------------------------------------------------------- + * Clear the matrix coefficients + *-----------------------------------------------------------------------*/ + + hypre_SetIndex(data_stride, 1, 1, 1); + + int_box = hypre_BoxCreate(); + + for (i = istart; i < istop; i++) + { + grid_box = hypre_BoxArrayBox(grid_boxes, i); + data_box = hypre_BoxArrayBox(data_space, i); + + hypre_IntersectBoxes(clear_box, grid_box, int_box); + + /* if there was an intersection */ + if (int_box) + { + data_start = hypre_BoxIMin(int_box); + + for (s = 0; s < num_stencil_indices; s++) + { + datap = hypre_StructMatrixBoxData(matrix, i, + stencil_indices[s]); + + hypre_BoxGetSize(int_box, loop_size); + + hypre_BoxLoop1Begin(loop_size, + data_box,data_start,data_stride,datai); +#define HYPRE_BOX_SMP_PRIVATE loopk,loopi,loopj,datai +#include "hypre_box_smp_forloop.h" + hypre_BoxLoop1For(loopi, loopj, loopk, datai) + { + datap[datai] = 0.0; + } + hypre_BoxLoop1End(datai); + } + } + } + + hypre_BoxDestroy(int_box); + + return hypre_error_flag; +} + /*-------------------------------------------------------------------------- *--------------------------------------------------------------------------*/ diff --git a/struct_mv/struct_vector.c b/struct_mv/struct_vector.c index 95bddaf52..611db2467 100644 --- a/struct_mv/struct_vector.c +++ b/struct_mv/struct_vector.c @@ -426,6 +426,156 @@ hypre_StructVectorSetBoxValues( hypre_StructVector *vector, return hypre_error_flag; } +/*-------------------------------------------------------------------------- + * (outside > 0): clear values possibly outside of the grid extents + * (outside = 0): clear values only inside the grid extents + *--------------------------------------------------------------------------*/ + +int +hypre_StructVectorClearValues( hypre_StructVector *vector, + hypre_Index grid_index, + int boxnum, + int outside ) +{ + hypre_BoxArray *grid_boxes; + hypre_Box *grid_box; + + double *vecp; + + int i, istart, istop; + + if (outside > 0) + { + grid_boxes = hypre_StructVectorDataSpace(vector); + } + else + { + grid_boxes = hypre_StructGridBoxes(hypre_StructVectorGrid(vector)); + } + + if (boxnum < 0) + { + istart = 0; + istop = hypre_BoxArraySize(grid_boxes); + } + else + { + istart = boxnum; + istop = istart + 1; + } + + for (i = istart; i < istop; i++) + { + grid_box = hypre_BoxArrayBox(grid_boxes, i); + + if ((hypre_IndexX(grid_index) >= hypre_BoxIMinX(grid_box)) && + (hypre_IndexX(grid_index) <= hypre_BoxIMaxX(grid_box)) && + (hypre_IndexY(grid_index) >= hypre_BoxIMinY(grid_box)) && + (hypre_IndexY(grid_index) <= hypre_BoxIMaxY(grid_box)) && + (hypre_IndexZ(grid_index) >= hypre_BoxIMinZ(grid_box)) && + (hypre_IndexZ(grid_index) <= hypre_BoxIMaxZ(grid_box)) ) + { + vecp = hypre_StructVectorBoxDataValue(vector, i, grid_index); + *vecp = 0.0; + } + } + + return hypre_error_flag; +} + +/*-------------------------------------------------------------------------- + * (outside > 0): clear values possibly outside of the grid extents + * (outside = 0): clear values only inside the grid extents + *--------------------------------------------------------------------------*/ + +int +hypre_StructVectorClearBoxValues( hypre_StructVector *vector, + hypre_Box *clear_box, + int boxnum, + int outside ) +{ + hypre_BoxArray *grid_boxes; + hypre_Box *grid_box; + hypre_Box *int_box; + + hypre_BoxArray *data_space; + hypre_Box *data_box; + hypre_IndexRef data_start; + hypre_Index data_stride; + int datai; + double *datap; + + hypre_Index loop_size; + + int i, istart, istop; + int loopi, loopj, loopk; + + /*----------------------------------------------------------------------- + * Initialize some things + *-----------------------------------------------------------------------*/ + + if (outside > 0) + { + grid_boxes = hypre_StructVectorDataSpace(vector); + } + else + { + grid_boxes = hypre_StructGridBoxes(hypre_StructVectorGrid(vector)); + } + data_space = hypre_StructVectorDataSpace(vector); + + if (boxnum < 0) + { + istart = 0; + istop = hypre_BoxArraySize(grid_boxes); + } + else + { + istart = boxnum; + istop = istart + 1; + } + + /*----------------------------------------------------------------------- + * Set the vector coefficients + *-----------------------------------------------------------------------*/ + + hypre_SetIndex(data_stride, 1, 1, 1); + + int_box = hypre_BoxCreate(); + + for (i = istart; i < istop; i++) + { + grid_box = hypre_BoxArrayBox(grid_boxes, i); + data_box = hypre_BoxArrayBox(data_space, i); + + hypre_IntersectBoxes(clear_box, grid_box, int_box); + + /* if there was an intersection */ + if (int_box) + { + data_start = hypre_BoxIMin(int_box); + + datap = hypre_StructVectorBoxData(vector, i); + + hypre_BoxGetSize(int_box, loop_size); + + hypre_BoxLoop1Begin(loop_size, + data_box,data_start,data_stride,datai); +#define HYPRE_BOX_SMP_PRIVATE loopk,loopi,loopj,datai +#include "hypre_box_smp_forloop.h" + hypre_BoxLoop1For(loopi, loopj, loopk, datai) + { + datap[datai] = 0.0; + } + hypre_BoxLoop1End(datai); + } + } + + hypre_BoxDestroy(int_box); + + return hypre_error_flag; +} + /*-------------------------------------------------------------------------- *--------------------------------------------------------------------------*/