PetscFunctionBegin; PetscValidHeaderSpecific(mesh,DA_COOKIE,1); ALE_LOG_EVENT_BEGIN ierr = PetscObjectGetComm((PetscObject) mesh, &comm);CHKERRQ(ierr); ierr = MPI_Comm_rank(comm, &rank);CHKERRQ(ierr); ierr = MeshGetTopology(mesh, &topology);CHKERRQ(ierr); ierr = MeshGetOrientation(mesh, &orientation);CHKERRQ(ierr); dim = topology->diameter(); topology->setStratificationPolicy(ALE::Sieve::stratificationPolicyOnLocking); /* Partition the topology and orientation */ ierr = ComputePreSievePartition(orientation, topology->leaves(), "Orientation");CHKERRQ(ierr); ierr = PartitionPreSieve(orientation, "Orientation", 1);CHKERRQ(ierr); ierr = ComputeSievePartition(topology, "Topology");CHKERRQ(ierr); ierr = PartitionPreSieve(topology, "Topology", 1, &partitionTypes);CHKERRQ(ierr); topology->getLock(); topology->releaseLock(); topology->setStratificationPolicy(ALE::Sieve::stratificationPolicyOnMutation); /* Add the trivial vertex orientation */ ALE::Obj<ALE::Point_set> roots = topology->depthStratum(0); for(ALE::Point_set::iterator vertex_itor = roots->begin(); vertex_itor != roots->end(); vertex_itor++) { ALE::Point v = *vertex_itor; orientation->addCone(v, v); } /* Create element bundle */ elementBundle = ALE::IndexBundle(topology); if (debug) { elementBundle->setVerbosity(11); } elementBundle->setFiberDimensionByHeight(0, 1); elementBundle->computeOverlapIndices(); elementBundle->computeGlobalIndices(); elementBundle->getLock(); // lock the bundle so that the overlap indices do not change ierr = MeshSetElementBundle(mesh, elementBundle);CHKERRQ(ierr); ierr = createParallelCoordinates(mesh, dim, partitionTypes);CHKERRQ(ierr); ALE_LOG_EVENT_END PetscFunctionReturn(0); }
#undef __FUNCT__ #define __FUNCT__ "MeshUnify" /*@ MeshUnify - */ PetscErrorCode MeshUnify(Mesh mesh, Mesh *serialMesh) { ALE::Obj<ALE::Sieve> topology; ALE::Obj<ALE::PreSieve> orientation; ALE::Obj<ALE::IndexBundle> coordBundle; ALE::Obj<ALE::PreSieve> partitionTypes; ALE::Obj<ALE::Sieve> boundary; Vec serialCoordinates, coordinates; MPI_Comm comm; PetscMPIInt rank; PetscInt dim; PetscErrorCode ierr;
PetscFunctionBegin; PetscValidHeaderSpecific(mesh,DA_COOKIE,1); ierr = PetscObjectGetComm((PetscObject) mesh, &comm);CHKERRQ(ierr); ierr = MPI_Comm_rank(comm, &rank);CHKERRQ(ierr); ierr = MeshGetCoordinateBundle(mesh, &coordBundle);CHKERRQ(ierr); ierr = MeshGetCoordinates(mesh, &coordinates);CHKERRQ(ierr); /* Unify the topology and orientation */ ierr = MeshGetTopology(mesh, &topology);CHKERRQ(ierr); ierr = MeshGetOrientation(mesh, &orientation);CHKERRQ(ierr); ALE::Obj<ALE::Sieve> serialTopology = ALE::Sieve(comm); ALE::Obj<ALE::PreSieve> serialOrientation = ALE::PreSieve(comm); ALE::Obj<ALE::Point_set> base = topology->base(); ALE::Obj<ALE::Point_set> orientationBase = orientation->base(); ALE::Point partition(-1, 0);
ierr = MeshCreate(comm, serialMesh);CHKERRQ(ierr); serialOrientation->addCone(orientation->space(), partition); for(ALE::Point_set::iterator b_itor = orientationBase->begin(); b_itor != orientationBase->end(); b_itor++) { serialOrientation->addCone(orientation->cone(*b_itor), *b_itor); } ierr = PartitionPreSieve(serialOrientation, "Serial Orientation", 1);CHKERRQ(ierr); serialTopology->addCone(topology->space(), partition); for(ALE::Point_set::iterator b_itor = base->begin(); b_itor != base->end(); b_itor++) { serialTopology->addCone(topology->cone(*b_itor), *b_itor); } ierr = PartitionPreSieve(serialTopology, "Serial Topology", 1, &partitionTypes);CHKERRQ(ierr); ierr = MeshSetTopology(*serialMesh, serialTopology);CHKERRQ(ierr); ierr = MeshSetOrientation(*serialMesh, serialOrientation);CHKERRQ(ierr); /* Unify boundary */ ierr = MeshGetBoundary(mesh, &boundary);CHKERRQ(ierr); ALE::Obj<ALE::Sieve> serialBoundary = ALE::Sieve(comm); ALE::Obj<ALE::Point_set> boundaryBase = boundary->base();
serialBoundary->addCone(boundary->space(), partition); for(ALE::Point_set::iterator b_itor = boundaryBase->begin(); b_itor != boundaryBase->end(); b_itor++) { serialBoundary->addCone(boundary->cone(*b_itor), *b_itor); } ierr = PartitionPreSieve(serialBoundary, "Serial Boundary", 1);CHKERRQ(ierr); ierr = MeshSetBoundary(*serialMesh, serialBoundary);CHKERRQ(ierr); /* Create vertex bundle */ ALE::Obj<ALE::IndexBundle> serialVertexBundle = ALE::IndexBundle(serialTopology); serialVertexBundle->setFiberDimensionByDepth(0, 1); serialVertexBundle->computeOverlapIndices(); serialVertexBundle->computeGlobalIndices(); serialVertexBundle->getLock(); // lock the bundle so that the overlap indices do not change ierr = MeshSetVertexBundle(*serialMesh, serialVertexBundle);CHKERRQ(ierr); /* Create element bundle */ ALE::Obj<ALE::IndexBundle> serialElementBundle = ALE::IndexBundle(serialTopology); serialElementBundle->setFiberDimensionByHeight(0, 1); serialElementBundle->computeOverlapIndices(); serialElementBundle->computeGlobalIndices(); serialElementBundle->getLock(); // lock the bundle so that the overlap indices do not change ierr = MeshSetElementBundle(*serialMesh, serialElementBundle);CHKERRQ(ierr); /* Create coordinate bundle and storage */ ALE::Obj<ALE::IndexBundle> serialCoordBundle = ALE::IndexBundle(serialTopology); dim = topology->diameter(); serialCoordBundle->setFiberDimensionByDepth(0, dim); serialCoordBundle->computeOverlapIndices(); serialCoordBundle->computeGlobalIndices(); serialCoordBundle->getLock(); // lock the bundle so that the overlap indices do not change ierr = MeshCreateVector(*serialMesh, serialCoordBundle, debug, &serialCoordinates);CHKERRQ(ierr); ierr = VecSetBlockSize(serialCoordinates, dim);CHKERRQ(ierr); ierr = MeshSetCoordinateBundle(*serialMesh, serialCoordBundle);CHKERRQ(ierr); /* Setup mapping to unified storage */ VecScatter coordScatter;
ierr = MeshCreateMapping(mesh, coordBundle, partitionTypes, serialCoordBundle, &coordScatter);CHKERRQ(ierr); ierr = VecScatterBegin(coordinates, serialCoordinates, INSERT_VALUES, SCATTER_FORWARD, coordScatter);CHKERRQ(ierr); ierr = VecScatterEnd(coordinates, serialCoordinates, INSERT_VALUES, SCATTER_FORWARD, coordScatter);CHKERRQ(ierr); ierr = MeshSetCoordinates(*serialMesh, serialCoordinates);CHKERRQ(ierr); ierr = VecScatterDestroy(coordScatter);CHKERRQ(ierr); PetscFunctionReturn(0); }
ALE::Obj<ALE::Point_set> getLocal(MPI_Comm comm, ALE::Obj<ALE::Stack> spaceFootprint, ALE::Obj<ALE::Point_set> points) { ALE::Obj<ALE::Point_set> localPoints(new ALE::Point_set); ALE::Point proc(0, spaceFootprint->getCommRank());
for(ALE::Point_set::iterator p_itor = points->begin(); p_itor != points->end(); p_itor++) { if (*spaceFootprint->cone(*p_itor)->begin() != proc) continue; localPoints->insert(*p_itor); } return localPoints; }
PetscErrorCode MeshComputeOverlap(Mesh mesh) { ALE::Obj<ALE::Sieve> topology; ALE::Obj<ALE::Stack> spaceFootprint; PetscErrorCode ierr;
PetscFunctionBegin; ierr = MeshGetSpaceFootprint(mesh, &spaceFootprint);CHKERRQ(ierr); if (!spaceFootprint) { ierr = MeshGetTopology(mesh, &topology);CHKERRQ(ierr); spaceFootprint = topology->spaceFootprint(ALE::PreSieve::completionTypePoint, ALE::PreSieve::footprintTypeSupport, NULL); ierr = MeshSetSpaceFootprint(mesh, spaceFootprint);CHKERRQ(ierr); } PetscFunctionReturn(0); }
PetscErrorCode MeshGetDimLocalSize(Mesh mesh, int dim, PetscInt *size) { MPI_Comm comm; ALE::Obj<ALE::Sieve> topology; ALE::Obj<ALE::Stack> spaceFootprint; PetscErrorCode ierr;
PetscFunctionBegin; PetscValidHeaderSpecific(mesh, DA_COOKIE, 1); PetscValidIntPointer(size, 3); ierr = PetscObjectGetComm((PetscObject) mesh, &comm);CHKERRQ(ierr); ierr = MeshComputeOverlap(mesh);CHKERRQ(ierr); ierr = MeshGetTopology(mesh, &topology);CHKERRQ(ierr); ierr = MeshGetSpaceFootprint(mesh, &spaceFootprint);CHKERRQ(ierr); *size = getLocal(comm, *spaceFootprint, topology->depthStratum(dim))->size(); PetscFunctionReturn(0); }
PetscErrorCode MeshGetDimLocalRanges(Mesh mesh, int dim, PetscInt starts[]) { MPI_Comm comm; PetscInt localSize; PetscErrorCode ierr;
PetscFunctionBegin; PetscValidHeaderSpecific(mesh, DA_COOKIE, 1); PetscValidIntPointer(starts, 3); ierr = PetscObjectGetComm((PetscObject) mesh, &comm);CHKERRQ(ierr); ierr = MeshGetDimLocalSize(mesh, dim, &localSize);CHKERRQ(ierr); ierr = MPI_Allgather(&localSize, 1, MPI_INT, starts, 1, MPI_INT, comm);CHKERRQ(ierr); PetscFunctionReturn(0); }
PetscErrorCode MeshGetDimGlobalSize(Mesh mesh, int dim, PetscInt *size) { MPI_Comm comm; PetscInt localSize, globalSize; PetscErrorCode ierr;
PetscFunctionBegin; PetscValidHeaderSpecific(mesh, DA_COOKIE, 1); PetscValidIntPointer(size, 3); ierr = PetscObjectGetComm((PetscObject) mesh, &comm);CHKERRQ(ierr); ierr = MeshGetDimLocalSize(mesh, dim, &localSize);CHKERRQ(ierr); ierr = MPI_Allreduce(&localSize, &globalSize, 1, MPI_INT, MPI_SUM, comm);CHKERRQ(ierr); *size = globalSize; PetscFunctionReturn(0); }
#endif
Level:none
Location:src/dm/mesh/sieve/Simplicializer.cxx
Index of all DA routines
Table of Contents for all manual pages
Index of all manual pages