[cig-commits] r4847 - in long/3D/Gale/trunk/src/StGermain: .
Discretisation/Mesh/src
walter at geodynamics.org
walter at geodynamics.org
Wed Oct 11 13:47:00 PDT 2006
Author: walter
Date: 2006-10-11 13:46:59 -0700 (Wed, 11 Oct 2006)
New Revision: 4847
Modified:
long/3D/Gale/trunk/src/StGermain/
long/3D/Gale/trunk/src/StGermain/Discretisation/Mesh/src/MeshTopology.c
long/3D/Gale/trunk/src/StGermain/Discretisation/Mesh/src/MeshTopology.h
Log:
r2904 at earth: boo | 2006-10-11 13:42:37 -0700
r2820 at earth (orig r3808): LukeHodkinson | 2006-09-26 20:25:33 -0700
Reworking mesh topology to automate a couple of
features:
* Added a method to allow for 'cascading' of
incidence relations. This means a smaller
subset of incidence needs to be specified in
order to calculate the remaining relations.
* Added a feature for automatic calculation of
shadow elements, including communication of
their incidence relations to the appropriate
processors.
Property changes on: long/3D/Gale/trunk/src/StGermain
___________________________________________________________________
Name: svk:merge
- 1ef209d2-b310-0410-a72d-e20c9eb0015c:/cig:2903
afb6c753-b9d0-0310-b4e7-dbd8d91cdd35:/trunk/StGermain:3807
+ 1ef209d2-b310-0410-a72d-e20c9eb0015c:/cig:2904
afb6c753-b9d0-0310-b4e7-dbd8d91cdd35:/trunk/StGermain:3808
Modified: long/3D/Gale/trunk/src/StGermain/Discretisation/Mesh/src/MeshTopology.c
===================================================================
--- long/3D/Gale/trunk/src/StGermain/Discretisation/Mesh/src/MeshTopology.c 2006-10-11 20:46:57 UTC (rev 4846)
+++ long/3D/Gale/trunk/src/StGermain/Discretisation/Mesh/src/MeshTopology.c 2006-10-11 20:46:59 UTC (rev 4847)
@@ -41,6 +41,9 @@
#include "types.h"
#include "shortcuts.h"
+#include "CommTopology.h"
+#include "Decomp.h"
+#include "Decomp_Sync.h"
#include "MeshTopology.h"
@@ -84,7 +87,14 @@
void _MeshTopology_Init( MeshTopology* self ) {
self->nDims = 0;
- self->nEls = NULL;
+ self->nTDims = 0;
+
+ self->domains = NULL;
+ self->nRentals = NULL;
+ self->rentals = NULL;
+ self->nDomainEls = NULL;
+ self->shadowDepth = 0;
+
self->nIncEls = NULL;
self->incEls = NULL;
}
@@ -97,6 +107,8 @@
void _MeshTopology_Delete( void* topo ) {
MeshTopology* self = (MeshTopology*)topo;
+ MeshTopology_Destruct( self );
+
/* Delete the parent. */
_Stg_Component_Delete( self );
}
@@ -170,194 +182,773 @@
MeshTopology* self = (MeshTopology*)topo;
assert( self );
+ assert( nDims > 0 );
/* If we're changing dimensions, kill everything and begin again. */
MeshTopology_Destruct( self );
/* Set and allocate. */
self->nDims = nDims;
+ self->nTDims = nDims + 1;
if( nDims ) {
- self->nEls = Memory_Alloc_Array( unsigned, nDims, "MeshTopology::nEls" );
- self->nIncEls = Memory_Alloc_2DArray( unsigned*, nDims, nDims, "MeshTopology::nIncEls" );
- self->incEls = Memory_Alloc_2DArray( unsigned**, nDims, nDims, "MeshTopology::incEls" );
+ unsigned d_i;
+
+ self->domains = Memory_Alloc_Array( Decomp_Sync*, self->nTDims, "MeshTopology::domains" );
+ self->nDomainEls = Memory_Alloc_Array( unsigned, self->nTDims, "MeshTopology::nDomainEls" );
+ self->nIncEls = Memory_Alloc_2DArray( unsigned*, self->nTDims, self->nTDims, "MeshTopology::nIncEls" );
+ self->incEls = Memory_Alloc_2DArray( unsigned**, self->nTDims, self->nTDims, "MeshTopology::incEls" );
+ for( d_i = 0; d_i < self->nTDims; d_i++ ) {
+ unsigned d_j;
+
+ self->domains[d_i] = Decomp_Sync_New( "" );
+ self->nDomainEls[d_i] = 0;
+ Decomp_Sync_SetDecomp( self->domains[d_i], Decomp_New( "" ) );
+
+ for( d_j = 0; d_j < self->nTDims; d_j++ ) {
+ self->nIncEls[d_i][d_j] = NULL;
+ self->incEls[d_i][d_j] = NULL;
+ }
+ }
}
}
-void MeshTopology_SetNElements( void* topo, MeshTopology_ElementDim dim, unsigned nEls ) {
+void MeshTopology_SetElements( void* topo, MeshTopology_Dim dim, unsigned nEls, unsigned* els ) {
MeshTopology* self = (MeshTopology*)topo;
- unsigned d_i;
assert( self );
- assert( dim < self->nDims );
+ assert( dim < self->nTDims );
+ assert( !nEls || els );
- /* Clear existing incidence relating to or from this dimension. */
- for( d_i = 0; d_i < self->nDims; d_i++ ) {
- KillArray( self->nIncEls[dim][d_i] );
- KillArray( self->incEls[dim][d_i] );
- KillArray( self->incEls[d_i][dim] );
- }
+ /* Decompose elements. */
+ Decomp_Sync_Decompose( self->domains[dim], nEls, els );
- /* Set size. */
- self->nEls[dim] = nEls;
+ /* Set the number of domain elements. */
+ self->nDomainEls[dim] = self->domains[dim]->decomp->nLocals + self->domains[dim]->nRemotes;
}
-void MeshTopology_SetIncidence( void* topo,
- MeshTopology_ElementDim fromDim, unsigned elInd,
- MeshTopology_ElementDim toDim,
+void MeshTopology_SetIncidence( void* topo, MeshTopology_Dim fromDim, MeshTopology_Dim toDim,
unsigned* nIncEls, unsigned** incEls )
{
MeshTopology* self = (MeshTopology*)topo;
+ unsigned size;
unsigned e_i;
assert( self );
- assert( fromDim < self->nDims );
- assert( elInd < self->nEls[fromDim] );
- assert( toDim < self->nDims );
+ assert( fromDim < self->nTDims );
+ assert( toDim < self->nTDims );
assert( nIncEls );
assert( incEls );
+ assert( Decomp_Sync_GetDomainSize( self->domains[fromDim] ) );
+ assert( Decomp_Sync_GetDomainSize( self->domains[toDim] ) );
/* Clear existing incidence. */
KillArray( self->incEls[fromDim][toDim] );
KillArray( self->nIncEls[fromDim][toDim] );
/* Allocate for incoming data. */
- self->nIncEls[fromDim][toDim] = Memory_Alloc_Array( unsigned, self->nEls[fromDim],
- "MeshTopology::nIncEls[][]" );
- self->incEls[fromDim][toDim] = Memory_Alloc_2DComplex( unsigned, self->nEls[fromDim], nIncEls,
- "MeshTopology::incEls[][]" );
+ size = Decomp_Sync_GetDomainSize( self->domains[fromDim] );
+ self->nIncEls[fromDim][toDim] = Memory_Alloc_Array( unsigned, size, "MeshTopology::nIncEls[][]" );
+ self->incEls[fromDim][toDim] = Memory_Alloc_2DComplex( unsigned, size, nIncEls, "MeshTopology::incEls[][]" );
/* Copy the lot. */
- memcpy( self->nIncEls[fromDim][toDim], nIncEls, self->nEls[fromDim] * sizeof(unsigned) );
- for( e_i = 0; e_i < self->nEls[fromDim]; e_i++ ) {
+ memcpy( self->nIncEls[fromDim][toDim], nIncEls, size * sizeof(unsigned) );
+ for( e_i = 0; e_i < size; e_i++ ) {
memcpy( self->incEls[fromDim][toDim][e_i], incEls[e_i],
self->nIncEls[fromDim][toDim][e_i] * sizeof(unsigned) );
}
}
-void MeshTopology_Complete( void* topo, MeshTopology_ElementDim fromDim, MeshTopology_ElementDim toDim ) {
- static const unsigned maxEls[4][4][4] = {{{-1, -1, -1, -1}, {-1, -1, -1, -1}, {-1, -1, -1, -1}, {-1, -1, -1, -1}},
- {{2, 2, -1, -1}, {-1, 2, -1, -1}, {-1, -1, -1, -1}, {-1, -1, -1, -1}},
- {{4, 4, 4, -1}, {-1, 6, 2, -1}, {-1, -1, 8, -1}, {-1, -1, -1, -1}},
- {{6, 6, 12, 8}, {-1, 10, 4, 4}, {-1, -1, 24, 2}, {-1, -1, -1, 6}}};
- MeshTopology* self = (MeshTopology*)topo;
+void MeshTopology_Complete( void* topo ) {
+ MeshTopology* self = (MeshTopology*)topo;
+ unsigned d_i, d_j;
assert( self );
- /* Neighbours are a little different; treat them separately. */
- if( fromDim == toDim ) {
- MeshTopology_CompleteNbrs( self, fromDim, maxEls[self->nDims][fromDim][fromDim] );
+ /* Complete downwards. */
+ for( d_i = 0; d_i < self->nTDims - 1; d_i++ ) {
+ for( d_j = d_i + 2; d_j < self->nTDims; d_j++ ) {
+ if( !self->nIncEls[d_j][d_i] )
+ MeshTopology_Cascade( self, d_j, d_i );
+ }
}
- else {
- unsigned* invNIncEls;
- unsigned** invIncEls;
- unsigned* nIncEls;
- unsigned** incEls;
- unsigned max = maxEls[self->nDims][fromDim][toDim];
- unsigned e_i;
- assert( fromDim < self->nDims );
- assert( toDim < self->nDims );
- assert( self->nEls[fromDim] );
- assert( self->nEls[toDim] );
- assert( !self->nIncEls[fromDim][toDim] && !self->incEls[fromDim][toDim] );
- assert( self->nIncEls[toDim][fromDim] && self->incEls[toDim][fromDim] );
+ /* Invert missing up relations. */
+ for( d_i = 0; d_i < self->nTDims - 1; d_i++ ) {
+ for( d_j = d_i + 1; d_j < self->nTDims; d_j++ ) {
+ if( !self->nIncEls[d_i][d_j] )
+ MeshTopology_Invert( self, d_i, d_j );
+ }
+ }
- /* Shortcuts. */
- invNIncEls = self->nIncEls[toDim][fromDim];
- invIncEls = self->incEls[toDim][fromDim];
+ /* Build neighbourhoods. */
+ for( d_i = 0; d_i < self->nTDims; d_i++ ) {
+ if( !self->nIncEls[d_i][d_i] )
+ MeshTopology_Neighbourhood( self, d_i );
+ }
+}
- /* Allocate some tables. */
- nIncEls = Memory_Alloc_Array( unsigned, self->nEls[fromDim], "MeshTopology::nIncEls[][]" );
- incEls = Memory_Alloc_2DArray( unsigned, self->nEls[fromDim], max, "" );
+void MeshTopology_SetShadowDepth( void* topo, unsigned depth ) {
+ MeshTopology* self = (MeshTopology*)topo;
+ RangeSet*** shadows;
+ unsigned** incSizes;
+ unsigned*** shadowInc;
- /* Build up the tables. */
- for( e_i = 0; e_i < self->nEls[toDim]; e_i++ ) {
- unsigned inc_i;
+ assert( self );
+ assert( depth );
- /* For each element incident on this 'toDim' element, add to table. */
- for( inc_i = 0; inc_i < 4; inc_i++ ) {
- unsigned elInd = self->incEls[toDim][e_i][fromDim][inc_i];
+ /* Need to store shadow information for each processor in each topological dimension. */
+ self->shadowDepth = depth;
- incEls[elInd][nIncEls[elInd]++] = e_i;
- }
+ /* If there are no neighbouring processors, forget about it. */
+ if( !self->domains[MT_VERTEX]->commTopo->nInc )
+ return;
+
+ shadows = Memory_Alloc_2DArray_Unnamed( RangeSet**, self->nTDims, self->domains[MT_VERTEX]->commTopo->nInc );
+ MeshTopology_BuildTopShadows( self, shadows + self->nDims );
+ MeshTopology_BuildAllShadows( self, shadows );
+ MeshTopology_BuildShadowInc( self, shadows, &incSizes, &shadowInc );
+
+ /* Send to neighbours. */
+ MeshTopology_SendRecvShadows( self, shadows, incSizes, shadowInc );
+}
+
+void MeshTopology_Invert( void* topo, MeshTopology_Dim fromDim, MeshTopology_Dim toDim ) {
+ MeshTopology* self = (MeshTopology*)topo;
+ unsigned fromSize, toSize;
+ unsigned* invNIncEls;
+ unsigned** invIncEls;
+ unsigned* nIncEls;
+ unsigned** incEls;
+ unsigned e_i;
+
+ /* Sanity check. */
+ assert( self );
+ assert( fromDim < self->nTDims );
+ assert( toDim < self->nTDims );
+ assert( Decomp_Sync_GetDomainSize( self->domains[fromDim] ) );
+ assert( Decomp_Sync_GetDomainSize( self->domains[toDim] ) );
+ assert( !self->nIncEls[fromDim][toDim] && !self->incEls[fromDim][toDim] );
+ assert( self->nIncEls[toDim][fromDim] && self->incEls[toDim][fromDim] );
+
+ /* Shortcuts. */
+ fromSize = Decomp_Sync_GetDomainSize( self->domains[fromDim] );
+ toSize = Decomp_Sync_GetDomainSize( self->domains[toDim] );
+ invNIncEls = self->nIncEls[toDim][fromDim];
+ invIncEls = self->incEls[toDim][fromDim];
+
+ /* Allocate some tables. */
+ nIncEls = Memory_Alloc_Array( unsigned, fromSize, "MeshTopology::nIncEls[][]" );
+ memset( nIncEls, 0, fromSize * sizeof(unsigned) );
+
+ /* Two phase process: count the numbers then allocate and do it again. */
+ for( e_i = 0; e_i < toSize; e_i++ ) {
+ unsigned inc_i;
+
+ for( inc_i = 0; inc_i < invNIncEls[e_i]; inc_i++ ) {
+ unsigned elInd = invIncEls[e_i][inc_i];
+
+ nIncEls[elInd]++;
}
+ }
- /* Transfer to permanent storage. */
- self->nIncEls[fromDim][toDim] = nIncEls;
- self->incEls[fromDim][toDim] = Memory_Alloc_2DComplex( unsigned, self->nEls[fromDim], nIncEls,
- "MeshTopology::incEls[][]" );
- for( e_i = 0; e_i < self->nEls[fromDim]; e_i++ )
- memcpy( self->incEls[fromDim][toDim][e_i], incEls[e_i], nIncEls[e_i] * sizeof(unsigned) );
- Memory_Free( incEls );
+ /* Build up the tables. */
+ incEls = Memory_Alloc_2DComplex( unsigned, fromSize, nIncEls, "MeshTopology::incEls[][]" );
+ memset( nIncEls, 0, fromSize * sizeof(unsigned) );
+ for( e_i = 0; e_i < toSize; e_i++ ) {
+ unsigned inc_i;
+
+ for( inc_i = 0; inc_i < invNIncEls[e_i]; inc_i++ ) {
+ unsigned elInd = invIncEls[e_i][inc_i];
+
+ incEls[elInd][nIncEls[elInd]++] = e_i;
+ }
}
+
+ /* Transfer to permanent storage. */
+ self->nIncEls[fromDim][toDim] = nIncEls;
+ self->incEls[fromDim][toDim] = incEls;
}
-/*----------------------------------------------------------------------------------------------------------------------------------
-** Private Functions
-*/
+void MeshTopology_Cascade( void* topo, MeshTopology_Dim fromDim, MeshTopology_Dim toDim ) {
+ MeshTopology* self = (MeshTopology*)topo;
+ unsigned maxInc = 0;
+ unsigned* nIncEls;
+ unsigned** incEls;
+ unsigned e_i;
-void MeshTopology_CompleteNbrs( MeshTopology* self, MeshTopology_ElementDim dim, unsigned maxNbrs ) {
- MeshTopology_ElementDim toDim;
+ assert( self );
+ assert( fromDim < self->nTDims );
+ assert( toDim < self->nTDims );
+ assert( fromDim > toDim + 1 );
+
+ /* Determine maximum incidence for any element. */
+ for( e_i = 0; e_i < self->nDomainEls[fromDim]; e_i++ ) {
+ unsigned curInc = 0;
+ unsigned inc_i;
+
+ for( inc_i = 0; inc_i < self->nIncEls[fromDim][fromDim - 1][e_i]; inc_i++ ) {
+ unsigned incEl = self->incEls[fromDim][fromDim - 1][e_i][inc_i];
+
+ curInc += self->nIncEls[fromDim - 1][toDim][incEl];
+ }
+
+ if( curInc > maxInc )
+ maxInc = curInc;
+ }
+
+ /* Allocate. */
+ nIncEls = Memory_Alloc_Array_Unnamed( unsigned, self->nDomainEls[fromDim] );
+ incEls = Memory_Alloc_2DArray_Unnamed( unsigned, self->nDomainEls[fromDim], maxInc );
+ memset( nIncEls, 0, self->nDomainEls[fromDim] * sizeof(unsigned) );
+
+ /* Determine actual incidence. */
+ for( e_i = 0; e_i < self->nDomainEls[fromDim]; e_i++ ) {
+ unsigned inc_i;
+
+ for( inc_i = 0; inc_i < self->nIncEls[fromDim][fromDim - 1][e_i]; inc_i++ ) {
+ unsigned incEl = self->incEls[fromDim][fromDim - 1][e_i][inc_i];
+ unsigned inc_j;
+
+ for( inc_j = 0; inc_j < self->nIncEls[fromDim - 1][toDim][incEl]; inc_j++ ) {
+ unsigned target = self->incEls[fromDim - 1][toDim][incEl][inc_j];
+ unsigned e_j;
+
+ for( e_j = 0; e_j < nIncEls[e_i]; e_j++ ) {
+ if( incEls[e_i][e_j] == target )
+ break;
+ }
+ if( e_j == nIncEls[e_i] )
+ incEls[e_i][nIncEls[e_i]++] = target;
+ }
+ }
+ }
+
+ /* Set incidence and kill temporary arrays. */
+ MeshTopology_SetIncidence( self, fromDim, toDim, nIncEls, incEls );
+ FreeArray( nIncEls );
+ FreeArray( incEls );
+}
+
+void MeshTopology_Neighbourhood( void* topo, MeshTopology_Dim dim ) {
+ MeshTopology* self = (MeshTopology*)topo;
+ unsigned size;
+ MeshTopology_Dim toDim;
unsigned* nNbrs;
unsigned** nbrs;
unsigned e_i;
assert( self );
- assert( dim < self->nDims );
- assert( maxNbrs != -1 );
+ assert( dim < self->nTDims );
+ assert( Decomp_Sync_GetDomainSize( self->domains[dim] ) );
/* Vertex neighbours search upwards, all others search downwards. */
- toDim = (dim == vertex) ? edge : vertex;
+ size = Decomp_Sync_GetDomainSize( self->domains[dim] );
+ toDim = (dim == MT_VERTEX) ? MT_EDGE : MT_VERTEX;
- /* Allocate some space. */
- nNbrs = Memory_Alloc_Array( unsigned, self->nEls[dim], "MeshTopology::nIncEls[][]" );
- nbrs = Memory_Alloc_2DArray( unsigned, self->nEls[dim], maxNbrs, "" );
- memset( nNbrs, 0, self->nEls[dim] * sizeof(unsigned) );
+ /* Allocate some space for neighbour counts. */
+ nNbrs = Memory_Alloc_Array_Unnamed( unsigned, size );
+ /* Calculate maximum neighbours for each element. */
+ for( e_i = 0; e_i < size; e_i++ ) {
+ unsigned nNodes = self->nIncEls[dim][toDim][e_i];
+ unsigned* nodes = self->incEls[dim][toDim][e_i];
+ unsigned n_i;
+
+ nNbrs[e_i] = 0;
+ for( n_i = 0; n_i < nNodes; n_i++ ) {
+ unsigned nCurNbrs = self->nIncEls[toDim][dim][nodes[n_i]];
+
+ nNbrs[e_i] += nCurNbrs - 1;
+ }
+ }
+
+ /* Allocate for maximum neighbours and clear neighbour counts. */
+ nbrs = Memory_Alloc_2DComplex_Unnamed( unsigned, size, nNbrs );
+ memset( nNbrs, 0, size * sizeof(unsigned) );
+
/* Build neighbours for each element of dimension 'dim'. */
- for( e_i = 0; e_i < self->nEls[dim]; e_i++ ) {
+ for( e_i = 0; e_i < size; e_i++ ) {
unsigned nNodes = self->nIncEls[dim][toDim][e_i];
unsigned* nodes = self->incEls[dim][toDim][e_i];
unsigned n_i;
/* Build a set of unique element nodes' elements. */
for( n_i = 0; n_i < nNodes; n_i++ ) {
- unsigned nbr_i;
+ unsigned nCurNbrs = self->nIncEls[toDim][dim][nodes[n_i]];
+ unsigned e_j;
- /* Don't add current element. */
- if( nodes[n_i] == e_i ) continue;
+ for( e_j = 0; e_j < nCurNbrs; e_j++ ) {
+ unsigned curNbr = self->incEls[toDim][dim][nodes[n_i]][e_j];
+ unsigned nbr_i;
- /* Ensure is unique. */
- for( nbr_i = 0; nbr_i < nNbrs[e_i]; nbr_i++ )
- if( nbrs[e_i][nbr_i] == nodes[n_i] ) break;
- if( nbr_i < nNbrs[e_i] )
- nbrs[e_i][nNbrs[e_i]++] = nodes[n_i];
+ /* Don't add current element. */
+ if( curNbr == e_i ) continue;
+
+ /* Ensure is unique. */
+ for( nbr_i = 0; nbr_i < nNbrs[e_i]; nbr_i++ ) {
+ if( nbrs[e_i][nbr_i] == curNbr )
+ break;
+ }
+ if( nbr_i == nNbrs[e_i] )
+ nbrs[e_i][nNbrs[e_i]++] = curNbr;
+ }
}
}
/* Transfer. */
- self->nIncEls[dim][dim] = nNbrs;
- self->incEls[dim][dim] = Memory_Alloc_2DComplex( unsigned, self->nEls[dim], nNbrs,
- "MeshTopology::incEls[][]" );
- for( e_i = 0; e_i < self->nEls[dim]; e_i++ )
- memcpy( self->incEls[dim][dim][e_i], nbrs[e_i], nNbrs[e_i] * sizeof(unsigned) );
- Memory_Free( nbrs );
+ MeshTopology_SetIncidence( self, dim, dim, nNbrs, nbrs );
+ FreeArray( nNbrs );
+ FreeArray( nbrs );
}
+unsigned MeshTopology_GetLocalSize( void* meshTopology, MeshTopology_Dim dim ) {
+ MeshTopology* self = (MeshTopology*)meshTopology;
+
+ assert( self );
+ assert( dim < self->nTDims );
+ assert( self->domains );
+ assert( self->domains[dim] );
+ assert( self->domains[dim]->decomp );
+
+ return self->domains[dim]->decomp->nLocals;
+}
+
+unsigned MeshTopology_GetShadowSize( void* meshTopology, MeshTopology_Dim dim ) {
+ MeshTopology* self = (MeshTopology*)meshTopology;
+
+ assert( self );
+ assert( dim < self->nTDims );
+ assert( self->domains );
+ assert( self->domains[dim] );
+
+ return self->domains[dim]->nRemotes;
+}
+
+unsigned MeshTopology_GetDomainSize( void* meshTopology, MeshTopology_Dim dim ) {
+ MeshTopology* self = (MeshTopology*)meshTopology;
+
+ assert( self );
+ assert( dim < self->nTDims );
+ assert( self->nDomainEls );
+
+ return self->nDomainEls[dim];
+}
+
+
+/*----------------------------------------------------------------------------------------------------------------------------------
+** Private Functions
+*/
+
+void MeshTopology_BuildTopShadows( void* meshTopology, RangeSet*** shadows ) {
+ MeshTopology* self = (MeshTopology*)meshTopology;
+ CommTopology* commTopo;
+ unsigned nIncRanks;
+ unsigned* incRanks;
+ void** shdSets;
+ Decomp_Sync* sync;
+ unsigned p_i, s_i;
+
+ assert( self );
+ assert( shadows );
+
+ /* Extract communicator incidence. */
+ commTopo = self->domains[MT_VERTEX]->commTopo;
+ CommTopology_GetIncidence( commTopo, commTopo->rank, &nIncRanks, &incRanks );
+
+ /* Build some sets, initially index sets. */
+ shdSets = (void**)Memory_Alloc_Array_Unnamed( IndexSet*, nIncRanks );
+ for( p_i = 0; p_i < nIncRanks; p_i++ )
+ shdSets[p_i] = (void*)IndexSet_New( Decomp_Sync_GetDomainSize( self->domains[MT_VERTEX] ) );
+
+ /* Store shared vertices. */
+ sync = self->domains[0];
+ for( s_i = 0; s_i < sync->nShared; s_i++ ) {
+ for( p_i = 0; p_i < sync->nSharers[s_i]; p_i++ )
+ IndexSet_Add( shdSets[sync->sharers[s_i][p_i]], sync->shared[s_i] );
+ }
+
+ /* Build range sets of shadowed elements. */
+ for( p_i = 0; p_i < nIncRanks; p_i++ ) {
+ unsigned nEdges;
+ unsigned* edges;
+ unsigned nShdEls;
+ unsigned* shdEls;
+ unsigned e_i;
+
+ IndexSet_GetMembers( shdSets[p_i], &nEdges, &edges );
+ FreeObject( shdSets[p_i] );
+ shdSets[p_i] = (void*)IndexSet_New( MeshTopology_GetLocalSize( self, self->nDims ) );
+
+ /* Build initial set of elements. */
+ for( e_i = 0; e_i < nEdges; e_i++ ) {
+ unsigned nIncEls = self->nIncEls[MT_VERTEX][self->nDims][edges[e_i]];
+ unsigned* incEls = self->incEls[MT_VERTEX][self->nDims][edges[e_i]];
+ unsigned inc_i;
+
+ for( inc_i = 0; inc_i < nIncEls; inc_i++ ) {
+ if( incEls[inc_i] >= MeshTopology_GetLocalSize( self, self->nDims ) )
+ continue;
+ IndexSet_Add( shdSets[p_i], incEls[inc_i] );
+ }
+ }
+
+ /* Free edge set. */
+ FreeArray( edges );
+
+ /* Expand the elements to satisfy shadow depth. */
+ if( self->shadowDepth > 1 ) {
+ IndexSet* bndSet;
+ unsigned d_i;
+
+ bndSet = IndexSet_DeepCopy( shdSets[p_i] );
+
+ for( d_i = 1; d_i < self->shadowDepth; d_i++ )
+ MeshTopology_ExpandShadows( self, (IndexSet*)shdSets[p_i], bndSet );
+
+ FreeObject( bndSet );
+ }
+
+ /* Get members and convert to global indices. */
+ IndexSet_GetMembers( shdSets[p_i], &nShdEls, &shdEls );
+ FreeObject( shdSets[p_i] );
+ for( e_i = 0; e_i < nShdEls; e_i++ )
+ shdEls[e_i] = Decomp_Sync_DomainToGlobal( self->domains[self->nDims], shdEls[e_i] );
+
+ /* Create range set. */
+ shdSets[p_i] = (void*)RangeSet_New();
+ RangeSet_SetIndices( shdSets[p_i], nShdEls, shdEls );
+ FreeArray( shdEls );
+ }
+
+ /* Store result. */
+ *shadows = (RangeSet**)shdSets;
+}
+
+void MeshTopology_ExpandShadows( void* meshTopology, IndexSet* shadows, IndexSet* boundary ) {
+ MeshTopology* self = (MeshTopology*)meshTopology;
+ unsigned nLocals = MeshTopology_GetLocalSize( self, self->nDims );
+ unsigned nBndEls;
+ unsigned* bndEls;
+ unsigned b_i;
+
+ assert( self );
+ assert( shadows );
+ assert( boundary );
+
+ IndexSet_GetMembers( boundary, &nBndEls, &bndEls );
+ FreeObject( boundary );
+ boundary = IndexSet_New( nLocals );
+
+ for( b_i = 0; b_i < nBndEls; b_i++ ) {
+ unsigned nNbrs = self->nIncEls[self->nDims][self->nDims][bndEls[b_i]];
+ unsigned* nbrs = self->incEls[self->nDims][self->nDims][bndEls[b_i]];
+ unsigned n_i;
+
+ for( n_i = 0; n_i < nNbrs; n_i++ ) {
+ if( nbrs[n_i] >= nLocals || IndexSet_IsMember( shadows, nbrs[n_i] ) )
+ continue;
+ IndexSet_Add( boundary, nbrs[n_i] );
+ }
+ }
+ FreeArray( bndEls );
+
+ IndexSet_GetMembers( boundary, &nBndEls, &bndEls );
+ for( b_i = 0; b_i < nBndEls; b_i++ )
+ IndexSet_Add( shadows, bndEls[b_i] );
+ FreeArray( bndEls );
+}
+
+void MeshTopology_BuildAllShadows( void* meshTopology, RangeSet*** shadows ) {
+ MeshTopology* self = (MeshTopology*)meshTopology;
+ CommTopology* commTopo;
+ unsigned nIncRanks;
+ unsigned* incRanks;
+ unsigned p_i;
+
+ assert( self );
+ assert( shadows );
+
+ /* Extract communicator incidence. */
+ commTopo = self->domains[MT_VERTEX]->commTopo;
+ CommTopology_GetIncidence( commTopo, commTopo->rank, &nIncRanks, &incRanks );
+
+ for( p_i = 0; p_i < nIncRanks; p_i++ ) {
+ unsigned nElInds;
+ unsigned* elInds;
+ unsigned d_i, e_i;
+
+ /* Convert range set of top level to indices. */
+ RangeSet_GetIndices( shadows[self->nDims][p_i], &nElInds, &elInds );
+ for( e_i = 0; e_i < nElInds; e_i++ )
+ elInds[e_i] = Decomp_Sync_GlobalToDomain( self->domains[self->nDims], elInds[e_i] );
+
+ /* */
+ for( d_i = 0; d_i < self->nDims; d_i++ ) {
+ IndexSet* iSet;
+ unsigned nInds;
+ unsigned* inds;
+
+ iSet = IndexSet_New( MeshTopology_GetLocalSize( self, d_i ) );
+
+ for( e_i = 0; e_i < nElInds; e_i++ ) {
+ unsigned nIncEls;
+ unsigned* incEls;
+ unsigned inc_i;
+
+ nIncEls = self->nIncEls[self->nDims][d_i][elInds[e_i]];
+ incEls = self->incEls[self->nDims][d_i][elInds[e_i]];
+ for( inc_i = 0; inc_i < nIncEls; inc_i++ ) {
+ if( incEls[inc_i] >= MeshTopology_GetLocalSize( self, d_i ) )
+ continue;
+ IndexSet_Add( iSet, incEls[inc_i] );
+ }
+ }
+
+ /* Get members and convert to global indices. */
+ IndexSet_GetMembers( iSet, &nInds, &inds );
+ FreeObject( iSet );
+ for( e_i = 0; e_i < nInds; e_i++ )
+ inds[e_i] = Decomp_Sync_DomainToGlobal( self->domains[d_i], inds[e_i] );
+
+ /* Create shadow set. */
+ shadows[d_i][p_i] = RangeSet_New();
+ RangeSet_SetIndices( shadows[d_i][p_i], nInds, inds );
+ FreeArray( inds );
+ }
+ }
+}
+
+void MeshTopology_BuildShadowInc( void* meshTopology, RangeSet*** shadows,
+ unsigned*** incDataSizes, unsigned**** incData )
+{
+ MeshTopology* self = (MeshTopology*)meshTopology;
+ CommTopology* commTopo;
+ unsigned nIncRanks;
+ unsigned* incRanks;
+ unsigned p_i;
+
+ assert( self );
+ assert( shadows );
+ assert( incDataSizes );
+ assert( incData );
+
+ /* Extract communicator incidence. */
+ commTopo = self->domains[MT_VERTEX]->commTopo;
+ CommTopology_GetIncidence( commTopo, commTopo->rank, &nIncRanks, &incRanks );
+
+ /* Allocate for the incidence data sizes. */
+ *incDataSizes = Memory_Alloc_2DArray_Unnamed( unsigned, self->nTDims, nIncRanks );
+ *incData = Memory_Alloc_2DArray_Unnamed( unsigned*, self->nTDims, nIncRanks );
+
+ for( p_i = 0; p_i < nIncRanks; p_i++ ) {
+ unsigned d_i;
+
+ for( d_i = 1; d_i < self->nTDims; d_i++ ) {
+ unsigned nEls;
+ unsigned* els;
+ unsigned offs = 0;
+ unsigned e_i;
+
+ /* Extract the indices. */
+ RangeSet_GetIndices( shadows[d_i][p_i], &nEls, &els );
+ for( e_i = 0; e_i < nEls; e_i++ )
+ els[e_i] = Decomp_Sync_GlobalToDomain( self->domains[d_i], els[e_i] );
+
+ /* Sum the data size for this dimension. */
+ (*incDataSizes)[d_i][p_i] = 1;
+ for( e_i = 0; e_i < nEls; e_i++ )
+ (*incDataSizes)[d_i][p_i] += self->nIncEls[d_i][d_i - 1][els[e_i]] + 2;
+
+ /* Allocate for incidence data. */
+ (*incData)[d_i][p_i] = Memory_Alloc_Array_Unnamed( unsigned, (*incDataSizes)[d_i][p_i] );
+
+ /* Store data. */
+ (*incData)[d_i][p_i][offs++] = nEls;
+ for( e_i = 0; e_i < nEls; e_i++ ) {
+ unsigned nIncEls = self->nIncEls[d_i][d_i - 1][els[e_i]];
+ unsigned* incEls = self->incEls[d_i][d_i - 1][els[e_i]];
+ unsigned gElInd = Decomp_Sync_DomainToGlobal( self->domains[d_i], els[e_i] );
+ unsigned inc_i;
+
+ (*incData)[d_i][p_i][offs++] = gElInd;
+ (*incData)[d_i][p_i][offs++] = nIncEls;
+ for( inc_i = 0; inc_i < nIncEls; inc_i++ )
+ (*incData)[d_i][p_i][offs++] = Decomp_Sync_DomainToGlobal( self->domains[d_i - 1], incEls[inc_i] );
+ }
+
+ /* Free the indices. */
+ FreeArray( els );
+ }
+ }
+}
+
+void MeshTopology_SendRecvShadows( void* meshTopology, RangeSet*** shadows,
+ unsigned** incDataSizes, unsigned*** incData )
+{
+ MeshTopology* self = (MeshTopology*)meshTopology;
+ CommTopology* commTopo;
+ unsigned nIncRanks;
+ unsigned* incRanks;
+ unsigned* nSrcBytes;
+ Stg_Byte** srcBytes;
+ unsigned* nDstBytes;
+ Stg_Byte** dstBytes;
+ unsigned p_i, d_i;
+
+ assert( self );
+ assert( shadows );
+ assert( incDataSizes );
+ assert( incData );
+
+ /* Extract communicator incidence. */
+ commTopo = self->domains[MT_VERTEX]->commTopo;
+ CommTopology_GetIncidence( commTopo, commTopo->rank, &nIncRanks, &incRanks );
+
+ /* Convert downward incidence relations to global indices. */
+ for( d_i = 1; d_i < self->nTDims; d_i++ ) {
+ unsigned e_i;
+
+ for( e_i = 0; e_i < self->nDomainEls[d_i]; e_i++ ) {
+ unsigned nIncEls = self->nIncEls[d_i][d_i - 1][e_i];
+ unsigned* incEls = self->incEls[d_i][d_i - 1][e_i];
+ unsigned inc_i;
+
+ for( inc_i = 0; inc_i < nIncEls; inc_i++ )
+ incEls[inc_i] = Decomp_Sync_DomainToGlobal( self->domains[d_i - 1], incEls[inc_i] );
+ }
+ }
+
+ /* Pickle all the range sets and pipe them off to our neighbours. */
+ nSrcBytes = Memory_Alloc_Array_Unnamed( unsigned, nIncRanks );
+ srcBytes = Memory_Alloc_Array_Unnamed( Stg_Byte*, nIncRanks );
+ for( d_i = 0; d_i < self->nTDims; d_i++ ) {
+ unsigned nRemotes;
+ unsigned* remotes;
+
+ for( p_i = 0; p_i < nIncRanks; p_i++ )
+ RangeSet_Pickle( shadows[d_i][p_i], nSrcBytes + p_i, srcBytes + p_i );
+ CommTopology_Alltoall( commTopo, nSrcBytes, (void**)srcBytes,
+ &nDstBytes, (void***)&dstBytes, sizeof(Stg_Byte) );
+ for( p_i = 0; p_i < nIncRanks; p_i++ ) {
+ FreeArray( srcBytes[p_i] );
+ RangeSet_Unpickle( shadows[d_i][p_i], nDstBytes[p_i], dstBytes[p_i] );
+ }
+ FreeArray( nDstBytes );
+ FreeArray( dstBytes );
+
+ /* Take the union of all other procs' shadow elements. */
+ for( p_i = 1; p_i < nIncRanks; p_i++ ) {
+ RangeSet_Union( shadows[d_i][0], shadows[d_i][p_i] );
+ FreeObject( shadows[d_i][p_i] );
+ }
+
+ /* Set the remote elements. */
+ RangeSet_GetIndices( shadows[d_i][0], &nRemotes, &remotes );
+ Decomp_Sync_SetRemotes( self->domains[d_i], nRemotes, remotes );
+ if( d_i == MT_VERTEX )
+ commTopo = self->domains[MT_VERTEX]->commTopo;
+ self->nDomainEls[d_i] = self->domains[d_i]->decomp->nLocals + nRemotes;
+ FreeArray( remotes );
+
+ /* Free range set. */
+ FreeObject( shadows[d_i][0] );
+ }
+ FreeArray( nSrcBytes );
+ FreeArray( srcBytes );
+ FreeArray( shadows );
+
+ /* Send and recieve all the incidence relations. */
+ for( d_i = 1; d_i < self->nTDims; d_i++ ) {
+ unsigned* dataSize;
+ unsigned** data;
+ unsigned** oldInc;
+ unsigned e_i;
+
+ CommTopology_Alltoall( commTopo, incDataSizes[d_i], (void**)incData[d_i],
+ &dataSize, (void***)&data, sizeof(unsigned) );
+
+ /* Resize old incidence relations to include new remote sizes. */
+ self->nIncEls[d_i][d_i - 1] = Memory_Realloc_Array( self->nIncEls[d_i][d_i - 1],
+ unsigned, self->nDomainEls[d_i] );
+
+ /* Insert each other procs' incidence relations. */
+ for( p_i = 0; p_i < nIncRanks; p_i++ ) {
+ unsigned* curData = data[p_i];
+ unsigned nEls = curData[0];
+
+ curData++;
+ for( e_i = 0; e_i < nEls; e_i++ ) {
+ unsigned dInd = Decomp_Sync_GlobalToDomain( self->domains[d_i], curData[0] );
+
+ self->nIncEls[d_i][d_i - 1][dInd] = curData[1];
+ curData += curData[1] + 2;
+ }
+ }
+
+ /* Resize incidence relations. */
+ oldInc = self->incEls[d_i][d_i - 1];
+ self->incEls[d_i][d_i - 1] = Memory_Alloc_2DComplex( unsigned, self->nDomainEls[d_i], self->nIncEls[d_i][d_i - 1],
+ "MeshTopology::incEls[][]" );
+ for( e_i = 0; e_i < MeshTopology_GetLocalSize( self, d_i ); e_i++ ) {
+ memcpy( self->incEls[d_i][d_i - 1][e_i], oldInc[e_i],
+ self->nIncEls[d_i][d_i - 1][e_i] * sizeof(unsigned) );
+ }
+ FreeArray( oldInc );
+
+ /* Copy new relations. */
+ for( p_i = 0; p_i < nIncRanks; p_i++ ) {
+ unsigned* curData = data[p_i];
+ unsigned nEls = curData[0];
+
+ curData++;
+ for( e_i = 0; e_i < nEls; e_i++ ) {
+ unsigned dInd = Decomp_Sync_GlobalToDomain( self->domains[d_i], curData[0] );
+ unsigned inc_i;
+
+ for( inc_i = 0; inc_i < curData[1]; inc_i++ )
+ self->incEls[d_i][d_i - 1][dInd][inc_i] = curData[2 + inc_i];
+ curData += curData[1] + 2;
+ }
+ }
+ }
+
+ /* Convert downward incidence relations to back to domain indices. */
+ for( d_i = 1; d_i < self->nTDims; d_i++ ) {
+ unsigned e_i;
+
+ for( e_i = 0; e_i < self->nDomainEls[d_i]; e_i++ ) {
+ unsigned nIncEls = self->nIncEls[d_i][d_i - 1][e_i];
+ unsigned* incEls = self->incEls[d_i][d_i - 1][e_i];
+ unsigned inc_i;
+
+ for( inc_i = 0; inc_i < nIncEls; inc_i++ )
+ incEls[inc_i] = Decomp_Sync_GlobalToDomain( self->domains[d_i - 1], incEls[inc_i] );
+ }
+ }
+}
+
void MeshTopology_Destruct( MeshTopology* self ) {
unsigned d_i;
assert( self );
- for( d_i = 0; d_i < self->nDims; d_i++ ) {
+ KillArray( self->nRentals );
+ KillArray( self->rentals );
+ KillArray( self->nDomainEls );
+ self->shadowDepth = 0;
+
+ for( d_i = 0; d_i < self->nTDims; d_i++ ) {
unsigned d_j;
- for( d_j = 0; d_j < self->nDims; d_j++ ) {
+ FreeObject( self->domains[d_i]->decomp );
+ FreeObject( self->domains[d_i] );
+
+ for( d_j = 0; d_j < self->nTDims; d_j++ ) {
KillArray( self->incEls[d_i][d_j] );
KillArray( self->nIncEls[d_i][d_j] );
}
}
KillArray( self->incEls );
KillArray( self->nIncEls );
- KillArray( self->nEls );
+ KillArray( self->domains );
}
Modified: long/3D/Gale/trunk/src/StGermain/Discretisation/Mesh/src/MeshTopology.h
===================================================================
--- long/3D/Gale/trunk/src/StGermain/Discretisation/Mesh/src/MeshTopology.h 2006-10-11 20:46:57 UTC (rev 4846)
+++ long/3D/Gale/trunk/src/StGermain/Discretisation/Mesh/src/MeshTopology.h 2006-10-11 20:46:59 UTC (rev 4847)
@@ -46,29 +46,43 @@
/** Virtual function types */
- /** Mesh class contents */
- #define __MeshTopology \
- /* General info */ \
- __Stg_Component \
- \
- /* Virtual info */ \
- \
- /* MeshTopology info */ \
- unsigned nDims; \
- unsigned* nEls; \
- unsigned*** nIncEls; \
- unsigned**** incEls; \
+ /** MeshTopology class contents */
+ typedef enum {
+ MT_VERTEX = 0,
+ MT_EDGE = 1,
+ MT_FACE = 2,
+ MT_VOLUME = 3
+ } MeshTopology_Dim;
+ #define __MeshTopology \
+ /* General info */ \
+ __Stg_Component \
+ \
+ /* Virtual info */ \
+ \
+ /* MeshTopology info */ \
+ unsigned nDims; \
+ unsigned nTDims; \
+ \
+ Decomp_Sync** domains; \
+ unsigned* nRentals; \
+ unsigned** rentals; \
+ unsigned* nDomainEls; \
+ unsigned shadowDepth; \
+ \
+ unsigned*** nIncEls; \
+ unsigned**** incEls;
+
struct MeshTopology { __MeshTopology };
/*--------------------------------------------------------------------------------------------------------------------------
** Constructors
*/
- #define MESHTOPOLOGY_DEFARGS \
+ #define MESHTOPOLOGY_DEFARGS \
STG_COMPONENT_DEFARGS
- #define MESHTOPOLOGY_PASSARGS \
+ #define MESHTOPOLOGY_PASSARGS \
STG_COMPONENT_PASSARGS
MeshTopology* MeshTopology_New( Name name );
@@ -88,8 +102,8 @@
(Mesh*)Stg_Class_Copy( self, NULL, True, NULL, NULL )
void* _MeshTopology_Copy( void* topo, void* dest, Bool deep, Name nameExt, PtrMap* ptrMap );
+ void _MeshTopology_Construct( void* topo, Stg_ComponentFactory* cf );
void _MeshTopology_Build( void* topo, void* data );
- void _MeshTopology_Construct( void* topo, Stg_ComponentFactory* cf );
void _MeshTopology_Initialise( void* topo, void* data );
void _MeshTopology_Execute( void* topo, void* data );
void _MeshTopology_Destroy( void* topo, void* data );
@@ -99,18 +113,30 @@
*/
void MeshTopology_SetNDims( void* topo, unsigned nDims );
- void MeshTopology_SetNElements( void* topo, MeshTopology_ElementDim dim, unsigned nEls );
- void MeshTopology_SetIncidence( void* topo,
- MeshTopology_ElementDim fromDim, unsigned elInd,
- MeshTopology_ElementDim toDim,
+ void MeshTopology_SetElements( void* topo, MeshTopology_Dim dim, unsigned nEls, unsigned* els );
+ void MeshTopology_SetIncidence( void* topo, MeshTopology_Dim fromDim, MeshTopology_Dim toDim,
unsigned* nIncEls, unsigned** incEls );
- void MeshTopology_Complete( void* topo, MeshTopology_ElementDim fromDim, MeshTopology_ElementDim toDim );
+ void MeshTopology_Complete( void* topo );
+ void MeshTopology_SetShadowDepth( void* topo, unsigned depth );
+ void MeshTopology_Invert( void* topo, MeshTopology_Dim fromDim, MeshTopology_Dim toDim );
+ void MeshTopology_Cascade( void* topo, MeshTopology_Dim fromDim, MeshTopology_Dim toDim );
+ void MeshTopology_Neighbourhood( void* topo, MeshTopology_Dim dim );
+ unsigned MeshTopology_GetLocalSize( void* meshTopology, MeshTopology_Dim dim );
+ unsigned MeshTopology_GetShadowSize( void* meshTopology, MeshTopology_Dim dim );
+ unsigned MeshTopology_GetDomainSize( void* meshTopology, MeshTopology_Dim dim );
+
/*--------------------------------------------------------------------------------------------------------------------------
** Private Member functions
*/
- void MeshTopology_CompleteNbrs( MeshTopology* self, MeshTopology_ElementDim dim, unsigned maxNbrs );
+ void MeshTopology_BuildTopShadows( void* meshTopology, RangeSet*** shadows );
+ void MeshTopology_ExpandShadows( void* meshTopology, IndexSet* shadows, IndexSet* boundary );
+ void MeshTopology_BuildAllShadows( void* meshTopology, RangeSet*** shadows );
+ void MeshTopology_BuildShadowInc( void* meshTopology, RangeSet*** shadows,
+ unsigned*** incSizes, unsigned**** shadowInc );
+ void MeshTopology_SendRecvShadows( void* meshTopology, RangeSet*** shadows,
+ unsigned** incSizes, unsigned*** shadowInc );
void MeshTopology_Destruct( MeshTopology* self );
#endif /* __Discretisaton_Mesh_MeshTopology_h__ */
More information about the cig-commits
mailing list