CFD/hypre
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Module: sstruct_mv

Browse Source 
Data Structure/Function:sstruct_matrix
Changes: added a type that for now is HYPRE_PARCSR or HYPRE_SSTRUCT
Reason: This will determine the dimension and entry location in the IJMatrix side of the structure
File: sstruct_mv.h sstruct_matrix.h

Module: sstruct_mv
Data Structure/Function: HYPRE_SStructMatrixCReate
Changes: the default of the type was set to HYPRE_SSTRUCT
         The creation of the IJmatrix was moved to the initialization phase
Reason: The ijcreation is postponed until the user decided the type of matrix and vector she/he will
        use to solve the problem
File: HYPRE_sstruct_matrix.c


Module: sstruct_mv
Data Structure/Function: HYPRE_SStruct_MAtrixSetObjectType
Changes: to inject the type into the sstruct_matrix
Reason: The user at the end decides whether to use parcsr or sstruct solver
File: HYPRE_sstruct_matrix.c

Module: sstruct_mv
Data Structure/Function: HYPRE_SStructMatrixInitialize
Changes: the ijmatrix is created here with a dimension that depends on the type.
Reason: It could have been done in the create part but it seems more natural after running the setobjectype.
File:HYPRE_sstruct_matrix.c

Module: sstruct_mv
Data Structure/Function:hypre_SStructUMatrixInitialize
Changes: When HYPRE_SSTRUCT type, you pass zero rowsizes to the IJSetRowSizes function.
Reason: To pass the rowsizes information, information about the nonzero rows is needed. HOwever, when
        the matrix is HYPRE_SSTRUCT all the rows associated with a ghost index are zero. But, we
        do not have this information stored in the sstruct_matrix, and it is too costly to have it.
        It turns out that passing zero for all rows is less costly than passing a constant size.

File: sstruct_matrix.c

Module: sstruct_mv
Data Structure/Function: hypre_SStructUMatrixSetBoxValues, UMatrixSetValues
Changes: Replace the calculation of the Global Rank for one that depends on the matrix type..
         hypre_SStructGetGlobalRank(mapentry,index,rank,type)
         instead of hypre_SStructGetGlobalRank(mapentry, index, rank)
Reason: the dimension and entries of the umatrix depend on the type of matrix
File: sstruct_matrix.c
This commit is contained in:
castilla 2002-12-13 05:17:10 +00:00
parent 0c53a29f5e
commit 8966cbc865
3 changed files with 187 additions and 33 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

50
sstruct_mv/HYPRE_sstruct_matrix.c
View File

@ -43,7 +43,6 @@ HYPRE_SStructMatrixCreate( MPI_Comm comm,
HYPRE_SStructVariable vitype, vjtype;
int part, vi, vj, i;
int size;
int ilower, iupper;
matrix = hypre_TAlloc(hypre_SStructMatrix, 1);
@ -95,10 +94,13 @@ HYPRE_SStructMatrixCreate( MPI_Comm comm,
}
}
ilower = hypre_SStructGridStartRank(grid);
iupper = ilower + hypre_SStructGridLocalSize(grid) - 1;
HYPRE_IJMatrixCreate(comm, ilower, iupper, ilower, iupper,
&hypre_SStructMatrixIJMatrix(matrix));
/* GEC0902 move the IJ creation to the initialization phase
* ilower = hypre_SStructGridGhstartRank(grid);
* iupper = ilower + hypre_SStructGridGhlocalSize(grid) - 1;
* HYPRE_IJMatrixCreate(comm, ilower, iupper, ilower, iupper,
* &hypre_SStructMatrixIJMatrix(matrix)); */
hypre_SStructMatrixIJMatrix(matrix) = NULL;
hypre_SStructMatrixParCSRMatrix(matrix) = NULL;
size = 0;
@ -120,6 +122,10 @@ HYPRE_SStructMatrixCreate( MPI_Comm comm,
hypre_SStructMatrixNSSymmetric(matrix) = 0;
hypre_SStructMatrixGlobalSize(matrix) = 0;
hypre_SStructMatrixRefCount(matrix) = 1;
/* GEC0902 setting the default of the object_type to 3333 HYPRE_SSTRUCT for now */
hypre_SStructMatrixObjectType(matrix) = 3333;
*matrix_ptr = matrix;
@ -210,6 +216,12 @@ HYPRE_SStructMatrixInitialize( HYPRE_SStructMatrix matrix )
int part, var, i;
/* GEC0902 addition of variables for ilower and iupper */
int comm;
hypre_SStructGrid *grid;
int ilower, iupper;
int matrix_type = hypre_SStructMatrixObjectType(matrix);
/* S-matrix */
for (part = 0; part < nparts; part++)
{
@ -256,6 +268,30 @@ HYPRE_SStructMatrixInitialize( HYPRE_SStructMatrix matrix )
hypre_SStructPMatrixInitialize(pmatrices[part]);
}
/* U-matrix */
/* GEC0902 knowing the kind of matrix we can create the IJMATRIX with the
* the right dimension (HYPRE_PARCSR without ghosts) */
grid = hypre_SStructGraphGrid(graph);
comm = hypre_SStructMatrixComm(matrix);
if(matrix_type == HYPRE_PARCSR)
{
ilower = hypre_SStructGridStartRank(grid);
iupper = ilower + hypre_SStructGridLocalSize(grid) - 1;
}
if(matrix_type == HYPRE_SSTRUCT)
{
ilower = hypre_SStructGridGhstartRank(grid);
iupper = ilower + hypre_SStructGridGhlocalSize(grid) - 1;
}
HYPRE_IJMatrixCreate(comm, ilower, iupper, ilower, iupper,
&hypre_SStructMatrixIJMatrix(matrix));
/* U-matrix */
hypre_SStructUMatrixInitialize(matrix);
@ -683,7 +719,11 @@ HYPRE_SStructMatrixSetObjectType( HYPRE_SStructMatrix matrix,
}
}
}
/* GEC0902 Injecting the type into the matrix data structure */
hypre_SStructMatrixObjectType(matrix) = type ;
return ierr;
}

166
sstruct_mv/sstruct_matrix.c
View File

@ -457,22 +457,39 @@ hypre_SStructUMatrixInitialize( hypre_SStructMatrix *matrix )
int *iUventries = hypre_SStructGraphIUVEntries(graph);
hypre_SStructUVEntry **Uventries = hypre_SStructGraphUVEntries(graph);
int **nvneighbors = hypre_SStructGridNVNeighbors(grid);
hypre_StructGrid *sgrid;
hypre_SStructStencil *stencil;
int *split;
int nvars;
int nrows, nnzs;
int part, var, entry, i, j;
int nrows, nnzs, snrows;
int part, var, entry, i, j, k,m,b;
int *row_sizes;
int max_row_size;
int matrix_type = hypre_SStructMatrixObjectType(matrix);
hypre_Box *gridbox;
hypre_Box *loopbox;
hypre_Box *ghostbox;
hypre_BoxArray *boxes;
int *num_ghost;
ierr = HYPRE_IJMatrixSetObjectType(ijmatrix, HYPRE_PARCSR);
nrows = hypre_SStructGridLocalSize(grid);
/* GEC1002 the ghlocalsize is used to set the number of rows */
if (matrix_type == HYPRE_PARCSR)
{
nrows = hypre_SStructGridLocalSize(grid);
}
if (matrix_type == HYPRE_SSTRUCT)
{
nrows = hypre_SStructGridGhlocalSize(grid) ;
}
/* set row sizes */
i = 0;
m = 0;
row_sizes = hypre_CTAlloc(int, nrows);
max_row_size = 0;
for (part = 0; part < nparts; part++)
@ -481,7 +498,7 @@ hypre_SStructUMatrixInitialize( hypre_SStructMatrix *matrix )
for (var = 0; var < nvars; var++)
{
sgrid = hypre_SStructPGridSGrid(pgrids[part], var);
nrows = hypre_StructGridLocalSize(sgrid);
stencil = stencils[part][var];
split = hypre_SStructMatrixSplit(matrix, part, var);
nnzs = 0;
@ -499,25 +516,83 @@ hypre_SStructUMatrixInitialize( hypre_SStructMatrix *matrix )
nnzs = 2*nnzs - 1;
}
#endif
for (j = 0; j < nrows; j++)
{
row_sizes[i] = nnzs;
max_row_size = hypre_max(max_row_size, row_sizes[i]);
i++;
/**************/
boxes = hypre_StructGridBoxes(sgrid) ;
num_ghost = hypre_StructGridNumGhost(sgrid);
for (b = 0; b < hypre_BoxArraySize(boxes); b++)
{
gridbox = hypre_BoxArrayBox(boxes, b);
ghostbox = hypre_BoxCreate();
loopbox = hypre_BoxCreate();
hypre_CopyBox(gridbox,ghostbox);
hypre_BoxExpand(ghostbox,num_ghost);
if (matrix_type == HYPRE_SSTRUCT)
{
hypre_CopyBox(ghostbox,loopbox);
}
if (matrix_type == HYPRE_PARCSR)
{
hypre_CopyBox(gridbox,loopbox);
}
for (k = hypre_BoxIMinZ(loopbox); k <= hypre_BoxIMaxZ(loopbox); k++)
{
for (j = hypre_BoxIMinY(loopbox); j <= hypre_BoxIMaxY(loopbox); j++)
{
for (i = hypre_BoxIMinX(loopbox); i <= hypre_BoxIMaxX(loopbox); i++)
{
if ( ( ( i>=hypre_BoxIMinX(gridbox) )
&& ( j>=hypre_BoxIMinY(gridbox) ) )
&& ( k>=hypre_BoxIMinZ(gridbox) ) )
{
if ( ( ( i<=hypre_BoxIMaxX(gridbox) )
&& ( j<=hypre_BoxIMaxY(gridbox) ) )
&& ( k<=hypre_BoxIMaxZ(gridbox) ) )
{
row_sizes[m] = nnzs;
max_row_size = hypre_max(max_row_size, row_sizes[m]);
}
}
m++;
}
}
}
hypre_BoxDestroy(ghostbox);
hypre_BoxDestroy(loopbox);
}
/* If there are vneighbors, use stencil size to set max_row_size */
if (nvneighbors[part][var])
{
max_row_size = hypre_max(max_row_size,
hypre_SStructStencilSize(stencil));
}
/*********************/
}
}
/* GEC0902 essentially for each UVentry we figure out how many extra columns
* we need to add to the rowsizes */
for (entry = 0; entry < nUventries; entry++)
{
i = iUventries[entry];
row_sizes[i] += hypre_SStructUVEntryNUEntries(Uventries[i]);
max_row_size = hypre_max(max_row_size, row_sizes[i]);
/* GEC0902 iUventries gives me the local rank (grid rank) whereas
* ighUventries gives me the localghostrank. Yet the Uventries are
* stored according to the grid no ghosts rank. */
i = iUventries[entry];
if (matrix_type == HYPRE_SSTRUCT)
{
row_sizes[i] += hypre_SStructUVEntryNUEntries(Uventries[i]);
max_row_size = hypre_max(max_row_size, row_sizes[i]);
}
}
/* ZTODO: Update row_sizes based on neighbor off-part couplings */
ierr += HYPRE_IJMatrixSetRowSizes (ijmatrix, (const int *) row_sizes);
@ -528,6 +603,8 @@ hypre_SStructUMatrixInitialize( hypre_SStructMatrix *matrix )
hypre_SStructMatrixTmpCoeffs(matrix) =
hypre_CTAlloc(double, max_row_size);
/* GEC1002 at this point the processor has the partitioning (creation of ij) */
ierr += HYPRE_IJMatrixInitialize(ijmatrix);
return ierr;
@ -567,6 +644,8 @@ hypre_SStructUMatrixSetValues( hypre_SStructMatrix *matrix,
double *coeffs;
int i, entry;
int proc, myproc;
/* GEC1002 the matrix type */
int matrix_type = hypre_SStructMatrixObjectType(matrix);
hypre_SStructGridFindMapEntry(grid, part, index, var, &map_entry);
if (map_entry == NULL)
@ -601,8 +680,12 @@ hypre_SStructUMatrixSetValues( hypre_SStructMatrix *matrix,
return ierr;
}
hypre_SStructMapEntryGetGlobalRank(map_entry, index, &row_coord);
/* GEC1002 get the rank using the function with the type=matrixtype*/
hypre_SStructMapEntryGetGlobalRank(map_entry, index, &row_coord, matrix_type);
col_coords = hypre_SStructMatrixTmpColCoords(matrix);
coeffs = hypre_SStructMatrixTmpCoeffs(matrix);
ncoeffs = 0;
@ -622,19 +705,31 @@ hypre_SStructUMatrixSetValues( hypre_SStructMatrix *matrix,
&map_entry);
if (map_entry != NULL)
{
hypre_SStructMapEntryGetGlobalRank(map_entry, to_index,
&col_coords[ncoeffs]);
coeffs[ncoeffs] = values[i];
ncoeffs++;
}
{
hypre_SStructMapEntryGetGlobalRank(map_entry, to_index,
&col_coords[ncoeffs],matrix_type);
coeffs[ncoeffs] = values[i];
ncoeffs++;
}
}
else
{
/* non-stencil entries */
entry -= size;
hypre_SStructGraphFindUVEntry(graph, part, index, var, &Uventry);
col_coords[ncoeffs] = hypre_SStructUVEntryRank(Uventry, entry);
if (matrix_type == HYPRE_PARCSR)
{
col_coords[ncoeffs] = hypre_SStructUVEntryRank(Uventry, entry);
}
if (matrix_type == HYPRE_SSTRUCT)
{
col_coords[ncoeffs] = hypre_SStructUVEntryGhrank(Uventry, entry);
}
coeffs[ncoeffs] = values[i];
ncoeffs++;
}
@ -701,6 +796,8 @@ hypre_SStructUMatrixSetBoxValues( hypre_SStructMatrix *matrix,
int row_base, col_base, val_base;
int e, entry, ii, jj, i, j, k;
int proc, myproc;
/* GEC1002 the matrix type */
int matrix_type = hypre_SStructMatrixObjectType(matrix);
box = hypre_BoxCreate();
@ -743,7 +840,9 @@ hypre_SStructUMatrixSetBoxValues( hypre_SStructMatrix *matrix,
continue;
}
hypre_SStructMapEntryGetStrides(map_entries[ii], rs);
/* GEC1002 introducing the strides based on the type of the matrix */
hypre_SStructMapEntryGetStrides(map_entries[ii], rs, matrix_type);
hypre_CopyIndex(ilower, hypre_BoxIMin(box));
hypre_CopyIndex(iupper, hypre_BoxIMax(box));
@ -775,7 +874,10 @@ hypre_SStructUMatrixSetBoxValues( hypre_SStructMatrix *matrix,
for (jj = 0; jj < nmap_to_entries; jj++)
{
hypre_SStructMapEntryGetStrides(map_to_entries[jj], cs);
/* GEC1002 introducing the strides based on the type of the matrix */
hypre_SStructMapEntryGetStrides(map_to_entries[jj], cs, matrix_type);
hypre_BoxMapEntryGetExtents(map_to_entries[jj],
hypre_BoxIMin(map_box),
@ -783,13 +885,21 @@ hypre_SStructUMatrixSetBoxValues( hypre_SStructMatrix *matrix,
hypre_IntersectBoxes(to_box, map_box, int_box);
hypre_CopyIndex(hypre_BoxIMin(int_box), index);
/* GEC1002 introducing the rank based on the type of the matrix */
hypre_SStructMapEntryGetGlobalRank(map_to_entries[jj],
index, &col_base);
index, &col_base,matrix_type);
hypre_IndexX(index) -= hypre_IndexX(offset);
hypre_IndexY(index) -= hypre_IndexY(offset);
hypre_IndexZ(index) -= hypre_IndexZ(offset);
/* GEC1002 introducing the rank based on the type of the matrix */
hypre_SStructMapEntryGetGlobalRank(map_entries[ii],
index, &row_base);
index, &row_base,matrix_type);
hypre_IndexX(index) -= hypre_IndexX(ilower);
hypre_IndexY(index) -= hypre_IndexY(ilower);
hypre_IndexZ(index) -= hypre_IndexZ(ilower);
@ -846,7 +956,7 @@ hypre_SStructUMatrixSetBoxValues( hypre_SStructMatrix *matrix,
hypre_BoxDestroy(to_box);
hypre_BoxDestroy(map_box);
hypre_BoxDestroy(int_box);
}
}
/*------------------------------------------
* non-stencil entries

4
sstruct_mv/sstruct_matrix.h
View File

@ -72,6 +72,9 @@ typedef struct hypre_SStructMatrix_struct
int ref_count;
/* GEC0902 adding an object type to the matrix */
int object_type;
} hypre_SStructMatrix;
/*--------------------------------------------------------------------------
@ -98,6 +101,7 @@ typedef struct hypre_SStructMatrix_struct
#define hypre_SStructMatrixComplex(mat) ((mat) -> complex)
#define hypre_SStructMatrixGlobalSize(mat) ((mat) -> global_size)
#define hypre_SStructMatrixRefCount(mat) ((mat) -> ref_count)
#define hypre_SStructMatrixObjectType(mat) ((mat) -> object_type)
/*--------------------------------------------------------------------------
* Accessor macros: hypre_SStructPMatrix