#ifdef PETSC_RCS_HEADER static char vcid[] = "$Id: meshBoundary.c,v 1.19 2000/10/17 13:48:57 knepley Exp $"; #endif /* Defines the interface to the mesh boundary functions */ #include "src/mesh/meshimpl.h" /*I "mesh.h" I*/ /*------------------------------------------- Mesh Boundary Creation ------------------------------------------------*/ #undef __FUNCT__ #define __FUNCT__ "MeshBoundary1DCreateSimple" /*@ MeshBoundary1DCreateSimple - This function creates the boundary of a simple linear mesh for the mesh generator. Input Parameter: . geomCtx - The problem geometry information Output Parameters in bdCtx: + numBd - The number of closed boundaries . numVertices - The number of vertices . vertices - The vertex coordinates . markers - The vertex markers . numSegments - The number of segments . segments - The segment endpoints . segMarkers - The segment markers . numHoles - The number of holes (particles) - holes - The hole coordinates Level: beginner .seealso MeshBoundary1DCreate(), MeshBoundary1DDestroy() .keywords mesh, boundary @*/ int MeshBoundary1DCreateSimple(MPI_Comm comm, MeshGeometryContext *geomCtx, MeshBoundary2D *bdCtx) { double length = geomCtx->size[0]; double xStart = geomCtx->start[0]; PetscTruth periodicMesh = geomCtx->isPeriodic[0]; double *vertices; int *markers; int rank; int ierr; PetscFunctionBegin; if (periodicMesh == PETSC_TRUE) { bdCtx->numBd = 0; bdCtx->numVertices = 0; } else { bdCtx->numBd = 2; bdCtx->numVertices = 2; } ierr = MPI_Comm_rank(comm, &rank); CHKERRQ(ierr); if (rank != 0) PetscFunctionReturn(0); /* Setup nodes */ ierr = PetscMalloc((bdCtx->numVertices)*2 * sizeof(double), &vertices); CHKERRQ(ierr); ierr = PetscMalloc((bdCtx->numVertices) * sizeof(int), &markers); CHKERRQ(ierr); if (periodicMesh == PETSC_FALSE) { vertices[0] = xStart; vertices[1] = xStart + length; markers[0] = BOTTOM_BD; markers[1] = TOP_BD; } bdCtx->vertices = vertices; bdCtx->markers = markers; PetscFunctionReturn(0); } #undef __FUNCT__ #define __FUNCT__ "MeshBoundary1DDestroy" /*@ MeshBoundary1DDestroy - This function frees the memory associated with the boundary of the initial mesh. Input Parameters in bdCtx: + vertices - The vertex coordinates . markers - The vertex markers . segments - The segment endpoints . segMarkers - The segment markers - holes - The hole coordinates Level: beginner .seealso MeshBoundary1DCreateSimple, MeshBoundary1DCreate .keywords mesh, boundary @*/ int MeshBoundary1DDestroy(MPI_Comm comm, MeshBoundary2D *bdCtx) { int rank; int ierr; PetscFunctionBegin; ierr = MPI_Comm_rank(comm, &rank); CHKERRQ(ierr); if (rank == 0) { ierr = PetscFree(bdCtx->vertices); CHKERRQ(ierr); ierr = PetscFree(bdCtx->markers); CHKERRQ(ierr); } PetscFunctionReturn(0); } #undef __FUNCT__ #define __FUNCT__ "MeshBoundary2DCreateSimple" /*@ MeshBoundary2DCreateSimple - This function creates the boundary of a simple square mesh for the mesh generator. Input Parameter: . geomCtx - The problem geometry information Output Parameters in bdCtx: + numBd - The number of closed boundaries . numVertices - The number of vertices . vertices - The vertex coordinates . markers - The vertex markers . numSegments - The number of segments . segments - The segment endpoints . segMarkers - The segment markers . numHoles - The number of holes (particles) - holes - The hole coordinates Level: beginner .seealso MeshBoundary2DCreate(), MeshBoundary2DDestroy() .keywords mesh, boundary @*/ int MeshBoundary2DCreateSimple(MPI_Comm comm, MeshGeometryContext *geomCtx, MeshBoundary2D *bdCtx) { double length = geomCtx->size[0]; double width = geomCtx->size[1]; double xStart = geomCtx->start[0]; double yStart = geomCtx->start[1]; PetscTruth periodicMesh = geomCtx->isPeriodic[0]; double *vertices; int *markers; int *segments; int *segMarkers; int rank; int node, seg; int ierr; PetscFunctionBegin; if (periodicMesh == PETSC_TRUE) { bdCtx->numBd = 2; bdCtx->numVertices = 12; bdCtx->numSegments = 20; bdCtx->numHoles = 0; } else { bdCtx->numBd = 1; bdCtx->numVertices = 9; bdCtx->numSegments = 12; bdCtx->numHoles = 0; } ierr = MPI_Comm_rank(comm, &rank); CHKERRQ(ierr); if (rank != 0) PetscFunctionReturn(0); /* Setup nodes */ ierr = PetscMalloc((bdCtx->numVertices)*2 * sizeof(double), &vertices); CHKERRQ(ierr); ierr = PetscMalloc((bdCtx->numVertices) * sizeof(int), &markers); CHKERRQ(ierr); if (periodicMesh == PETSC_TRUE) { vertices[0] = xStart; vertices[1] = yStart; vertices[2] = xStart + length/4.0; vertices[3] = yStart; vertices[4] = xStart + length/2.0; vertices[5] = yStart; vertices[6] = xStart + 3.0*length/4.0; vertices[7] = yStart; vertices[8] = xStart + length; vertices[9] = yStart + width; vertices[10] = xStart + 3.0*length/4.0; vertices[11] = yStart + width; vertices[12] = xStart + length/2.0; vertices[13] = yStart + width; vertices[14] = xStart + length/4.0; vertices[15] = yStart + width; vertices[16] = xStart; vertices[17] = yStart + width/2.0; vertices[18] = xStart + length/4.0; vertices[19] = yStart + width/2.0; vertices[20] = xStart + length/2.0; vertices[21] = yStart + width/2.0; vertices[22] = xStart + 3.0*length/4.0; vertices[23] = yStart + width/2.0; for(node = 0; node < 4; node++) markers[node] = BOTTOM_BD; for(node = 4; node < 8; node++) markers[node] = TOP_BD; for(node = 8; node < 12; node++) markers[node] = 0; } else { vertices[0] = xStart; vertices[1] = yStart; vertices[2] = xStart + length/2.0; vertices[3] = yStart; vertices[4] = xStart + length; vertices[5] = yStart; vertices[6] = xStart + length; vertices[7] = yStart + width/2.0; vertices[8] = xStart + length; vertices[9] = yStart + width; vertices[10] = xStart + length/2.0; vertices[11] = yStart + width; vertices[12] = xStart; vertices[13] = yStart + width; vertices[14] = xStart; vertices[15] = yStart + width/2.0; vertices[16] = xStart + length/2.0; vertices[17] = yStart + width/2.0; for(node = 0; node < 8; node++) markers[node] = OUTER_BD; markers[8] = 0; } /* Setup segments */ ierr = PetscMalloc((bdCtx->numSegments)*2 * sizeof(int), &segments); CHKERRQ(ierr); ierr = PetscMalloc((bdCtx->numSegments) * sizeof(int), &segMarkers); CHKERRQ(ierr); if (periodicMesh == PETSC_TRUE) { segments[0] = 0; segments[1] = 1; segments[2] = 1; segments[3] = 2; segments[4] = 2; segments[5] = 3; segments[6] = 3; segments[7] = 0; segments[8] = 4; segments[9] = 5; segments[10] = 5; segments[11] = 6; segments[12] = 6; segments[13] = 7; segments[14] = 7; segments[15] = 4; segments[16] = 8; segments[17] = 9; segments[18] = 9; segments[19] = 10; segments[20] = 10; segments[21] = 11; segments[22] = 11; segments[23] = 8; segments[24] = 0; segments[25] = 8; segments[26] = 4; segments[27] = 8; segments[28] = 1; segments[29] = 9; segments[30] = 7; segments[31] = 9; segments[32] = 2; segments[33] = 10; segments[34] = 6; segments[35] = 10; segments[36] = 3; segments[37] = 11; segments[38] = 5; segments[39] = 11; for(seg = 0; seg < 4; seg++) segMarkers[seg] = BOTTOM_BD; for(seg = 4; seg < 8; seg++) segMarkers[seg] = TOP_BD; for(seg = 8; seg < 20; seg++) segMarkers[seg] = 0; } else { segments[0] = 0; segments[1] = 1; segments[2] = 1; segments[3] = 2; segments[4] = 2; segments[5] = 3; segments[6] = 3; segments[7] = 4; segments[8] = 4; segments[9] = 5; segments[10] = 5; segments[11] = 6; segments[12] = 6; segments[13] = 7; segments[14] = 7; segments[15] = 0; segments[16] = 1; segments[17] = 8; segments[18] = 3; segments[19] = 8; segments[20] = 5; segments[21] = 8; segments[22] = 7; segments[23] = 8; for(seg = 0; seg < 8; seg++) segMarkers[seg] = OUTER_BD; for(seg = 8; seg < 12; seg++) segMarkers[seg] = 0; } bdCtx->vertices = vertices; bdCtx->markers = markers; bdCtx->segments = segments; bdCtx->segMarkers = segMarkers; bdCtx->holes = PETSC_NULL; PetscFunctionReturn(0); } #undef __FUNCT__ #define __FUNCT__ "MeshBoundary2DDestroy" /*@ MeshBoundary2DDestroy - This function frees the memory associated with the boundary of the initial mesh. Input Parameters in bdCtx: + vertices - The vertex coordinates . markers - The vertex markers . segments - The segment endpoints . segMarkers - The segment markers - holes - The hole coordinates Level: beginner .seealso MeshBoundary2DCreateSimple, MeshBoundary2DCreate .keywords mesh, boundary @*/ int MeshBoundary2DDestroy(MPI_Comm comm, MeshBoundary2D *bdCtx) { int rank; int ierr; PetscFunctionBegin; ierr = MPI_Comm_rank(comm, &rank); CHKERRQ(ierr); if (rank == 0) { ierr = PetscFree(bdCtx->vertices); CHKERRQ(ierr); ierr = PetscFree(bdCtx->markers); CHKERRQ(ierr); ierr = PetscFree(bdCtx->segments); CHKERRQ(ierr); ierr = PetscFree(bdCtx->segMarkers); CHKERRQ(ierr); if (bdCtx->numHoles > 0) { ierr = PetscFree(bdCtx->holes); CHKERRQ(ierr); } } PetscFunctionReturn(0); } /*---------------------------------------- Mesh Boundary Visualization ----------------------------------------------*/ #undef __FUNCT__ #define __FUNCT__ "MeshBoundary2DView_File" static int MeshBoundary2DView_File(Mesh mesh, MeshBoundary2D *bdCtx, PetscViewer viewer) { FILE *fd; int i; int ierr; PetscFunctionBegin; ierr = PetscViewerASCIIGetPointer(viewer, &fd); CHKERRQ(ierr); PetscFPrintf(PETSC_COMM_WORLD, fd, "Mesh Boundary Object:\n"); if (bdCtx->numBd == 1) { PetscFPrintf(PETSC_COMM_WORLD, fd, " %d closed boundary\n", bdCtx->numBd); } else { PetscFPrintf(PETSC_COMM_WORLD, fd, " %d closed boundaries\n", bdCtx->numBd); } PetscFPrintf(PETSC_COMM_WORLD, fd, " Boundary vertices: %d segments: %d holes: %d\n", bdCtx->numVertices, bdCtx->numSegments, bdCtx->numHoles); for(i = 0; i < bdCtx->numVertices; i++) { PetscFPrintf(PETSC_COMM_WORLD, fd, " %d %g %g %d\n", i, bdCtx->vertices[i*2], bdCtx->vertices[i*2+1], bdCtx->markers[i]); } for(i = 0; i < bdCtx->numSegments; i++) { PetscFPrintf(PETSC_COMM_WORLD, fd, " %d %d %d %d\n", i, bdCtx->segments[i*2], bdCtx->segments[i*2+1], bdCtx->segMarkers[i]); } for(i = 0; i < bdCtx->numHoles; i++) { PetscFPrintf(PETSC_COMM_WORLD, fd, " %d %g %g\n", i, bdCtx->holes[i*2], bdCtx->holes[i*2+1]); } PetscFunctionReturn(0); } #undef __FUNCT__ #define __FUNCT__ "MeshBoundary2DView_Draw_Mesh" static int MeshBoundary2DView_Draw_Mesh(Mesh mesh, MeshBoundary2D *bdCtx, PetscDraw draw) { double *vertices = bdCtx->vertices; int *segments = bdCtx->segments; int *segMarkers = bdCtx->segMarkers; int color; int seg; int ierr; PetscFunctionBegin; for(seg = 0; seg < bdCtx->numSegments; seg++) { if (segMarkers[seg]) { color = PETSC_DRAW_RED; } else { color = PETSC_DRAW_BLACK; } ierr = MeshDrawLine(mesh, draw, vertices[segments[seg*2]*2], vertices[segments[seg*2]*2+1], vertices[segments[seg*2+1]*2], vertices[segments[seg*2+1]*2+1], color); CHKERRQ(ierr); } ierr = PetscDrawSynchronizedFlush(draw); CHKERRQ(ierr); PetscFunctionReturn(0); } #undef __FUNCT__ #define __FUNCT__ "MeshBoundary2DView_Draw" static int MeshBoundary2DView_Draw(Mesh mesh, MeshBoundary2D *bdCtx, PetscViewer v) { PetscReal startX, endX; PetscReal startY, endY; PetscReal sizeX, sizeY; double *vertices; PetscDraw draw; PetscTruth isnull; PetscDrawButton button; char statusLine[256]; int pause; PetscReal xc, yc, scale; int xsize, ysize; int vert; int ierr; PetscFunctionBegin; ierr = PetscViewerDrawGetDraw(v, 0, &draw); CHKERRQ(ierr); ierr = PetscDrawIsNull(draw, &isnull); CHKERRQ(ierr); if (isnull) PetscFunctionReturn(0); ierr = MeshGetBoundingBox(mesh, &startX, &startY, PETSC_NULL, &endX, &endY, PETSC_NULL); CHKERRQ(ierr); vertices = bdCtx->vertices; for(vert = 0; vert < bdCtx->numVertices; vert++) { if (startX > vertices[vert*2]) startX = vertices[vert*2]; if (startY > vertices[vert*2+1]) startY = vertices[vert*2+1]; if (endX < vertices[vert*2]) endX = vertices[vert*2]; if (endY < vertices[vert*2+1]) endY = vertices[vert*2+1]; } sizeX = endX - startX; sizeY = endY - startY; if (sizeX > sizeY) { ysize = 300; xsize = ysize * (int) (sizeX/sizeY); } else { xsize = 300; ysize = xsize * (int) (sizeY/sizeX); } ierr = PetscDrawResizeWindow(draw, xsize, ysize); CHKERRQ(ierr); ierr = PetscDrawSynchronizedClear(draw); CHKERRQ(ierr); ierr = PetscDrawSetCoordinates(draw, startX-0.02*sizeX, startY-0.02*sizeY, endX+0.02*sizeX, endY+0.02*sizeY); CHKERRQ(ierr); ierr = MeshBoundary2DView_Draw_Mesh(mesh, bdCtx, draw); CHKERRQ(ierr); ierr = PetscDrawGetPause(draw, &pause); CHKERRQ(ierr); if (pause >= 0) { PetscSleep(pause); PetscFunctionReturn(0); } /* Allow the mesh to zoom or shrink */ ierr = PetscDrawCheckResizedWindow(draw); CHKERRQ(ierr); ierr = PetscDrawSynchronizedGetMouseButton(draw, &button, &xc, &yc, 0, 0); CHKERRQ(ierr); while ((button != BUTTON_RIGHT) && (button != BUTTON_CENTER_SHIFT)) { ierr = PetscDrawSynchronizedClear(draw); CHKERRQ(ierr); statusLine[0] = 0; switch (button) { case BUTTON_LEFT: scale = 0.5; break; case BUTTON_CENTER: scale = 2.0; break; default: scale = 1.0; break; } if (scale != 1.0) { startX = scale*(startX + sizeX - xc) + xc - sizeX*scale; endX = scale*(endX - sizeX - xc) + xc + sizeX*scale; startY = scale*(startY + sizeY - yc) + yc - sizeY*scale; endY = scale*(endY - sizeY - yc) + yc + sizeY*scale; sizeX *= scale; sizeY *= scale; } ierr = PetscDrawSetCoordinates(draw, startX-0.02*sizeX, startY-0.02*sizeY, endX+0.02*sizeX, endY+0.02*sizeY); CHKERRQ(ierr); ierr = PetscDrawString(draw, startX - 0.02*sizeX, startY - 0.02*sizeY, PETSC_DRAW_BLACK, statusLine); CHKERRQ(ierr); ierr = MeshBoundary2DView_Draw_Mesh(mesh, bdCtx, draw); CHKERRQ(ierr); ierr = PetscDrawCheckResizedWindow(draw); CHKERRQ(ierr); ierr = PetscDrawSynchronizedGetMouseButton(draw, &button, &xc, &yc, 0, 0); CHKERRQ(ierr); } PetscFunctionReturn(0); } #undef __FUNCT__ #define __FUNCT__ "MeshBoundary2DView" int MeshBoundary2DView(Mesh mesh, MeshBoundary2D *bdCtx, PetscViewer viewer) { PetscTruth isascii, isdraw; int ierr; PetscFunctionBegin; ierr = PetscTypeCompare((PetscObject) viewer, PETSC_VIEWER_ASCII, &isascii); CHKERRQ(ierr); ierr = PetscTypeCompare((PetscObject) viewer, PETSC_VIEWER_DRAW, &isdraw); CHKERRQ(ierr); if (isascii == PETSC_TRUE) { ierr = MeshBoundary2DView_File(mesh, bdCtx, viewer); CHKERRQ(ierr); } else if (isdraw == PETSC_TRUE) { ierr = MeshBoundary2DView_Draw(mesh, bdCtx, viewer); CHKERRQ(ierr); } PetscFunctionReturn(0); } /*------------------------------------------- Mesh Boundary Functions ------------------------------------------------*/ #undef __FUNCT__ #define __FUNCT__ "MeshBoundary2DSetBoundingBox" int MeshBoundary2DSetBoundingBox(Mesh mesh, MeshBoundary2D *bdCtx) { MPI_Comm comm; double startX, startY; double endX, endY; double *vertices; int rank; int p; int ierr; PetscFunctionBegin; ierr = PetscObjectGetComm((PetscObject) mesh, &comm); CHKERRQ(ierr); ierr = MPI_Comm_rank(comm, &rank); CHKERRQ(ierr); if (rank == 0) { vertices = bdCtx->vertices; startX = endX = vertices[0]; startY = endY = vertices[1]; for(p = 0; p < bdCtx->numVertices; p++) { /* Ignore moving boundaries (to accomodate split particles)? This should be done better */ if (bdCtx->markers[p] < 0) continue; if (startX > vertices[p*2]) startX = vertices[p*2]; if (startY > vertices[p*2+1]) startY = vertices[p*2+1]; if (endX < vertices[p*2]) endX = vertices[p*2]; if (endY < vertices[p*2+1]) endY = vertices[p*2+1]; } } ierr = MPI_Bcast(&startX, 1, MPI_DOUBLE, 0, comm); CHKERRQ(ierr); ierr = MPI_Bcast(&endX, 1, MPI_DOUBLE, 0, comm); CHKERRQ(ierr); ierr = MPI_Bcast(&startY, 1, MPI_DOUBLE, 0, comm); CHKERRQ(ierr); ierr = MPI_Bcast(&endY, 1, MPI_DOUBLE, 0, comm); CHKERRQ(ierr); ierr = MeshSetBoundingBox(mesh, startX, startY, 0.0, endX, endY, 0.0); CHKERRQ(ierr); PetscFunctionReturn(0); }