hypre_SStructSendInfoData * hypre_SStructSendInfo( hypre_StructGrid *fgrid, hypre_BoxManager *cboxman, hypre_Index rfactor ) { hypre_SStructSendInfoData *sendinfo_data; MPI_Comm comm= hypre_SStructVectorComm(fgrid); hypre_BoxArray *grid_boxes; hypre_Box *grid_box, cbox; hypre_Box *intersect_box, boxman_entry_box; hypre_BoxManEntry **boxman_entries; HYPRE_Int nboxman_entries; hypre_BoxArrayArray *send_boxes; HYPRE_Int **send_processes; HYPRE_Int **send_remote_boxnums; hypre_Index ilower, iupper, index; HYPRE_Int myproc, proc; HYPRE_Int cnt; HYPRE_Int i, j; hypre_ClearIndex(index); hypre_MPI_Comm_rank(comm, &myproc); sendinfo_data= hypre_CTAlloc(hypre_SStructSendInfoData, 1); /*------------------------------------------------------------------------ * Create the structured sendbox patterns. * * send_boxes are obtained by intersecting this proc's fgrid boxes * with cgrid's box_man. Intersecting BoxManEntries not on this proc * will give boxes that we will need to send data to- i.e., we scan * through the boxes of grid and find the processors that own a chunk * of it. *------------------------------------------------------------------------*/ intersect_box = hypre_CTAlloc(hypre_Box, 1); grid_boxes = hypre_StructGridBoxes(fgrid); send_boxes= hypre_BoxArrayArrayCreate(hypre_BoxArraySize(grid_boxes)); send_processes= hypre_CTAlloc(HYPRE_Int *, hypre_BoxArraySize(grid_boxes)); send_remote_boxnums= hypre_CTAlloc(HYPRE_Int *, hypre_BoxArraySize(grid_boxes)); hypre_ForBoxI(i, grid_boxes) { grid_box= hypre_BoxArrayBox(grid_boxes, i); /*--------------------------------------------------------------------- * Find the boxarray that must be sent. BoxManIntersect returns * the full extents of the boxes that intersect with the given box. * We further need to intersect each box in the list with the given * box to determine the actual box that needs to be sent. *---------------------------------------------------------------------*/ hypre_SStructIndexScaleF_C(hypre_BoxIMin(grid_box), index, rfactor, hypre_BoxIMin(&cbox)); hypre_SStructIndexScaleF_C(hypre_BoxIMax(grid_box), index, rfactor, hypre_BoxIMax(&cbox)); hypre_BoxManIntersect(cboxman, hypre_BoxIMin(&cbox), hypre_BoxIMax(&cbox), &boxman_entries, &nboxman_entries); cnt= 0; for (j= 0; j< nboxman_entries; j++) { hypre_SStructBoxManEntryGetProcess(boxman_entries[j], &proc); if (proc != myproc) { cnt++; } } send_processes[i] = hypre_CTAlloc(HYPRE_Int, cnt); send_remote_boxnums[i]= hypre_CTAlloc(HYPRE_Int, cnt); cnt= 0; for (j= 0; j< nboxman_entries; j++) { hypre_SStructBoxManEntryGetProcess(boxman_entries[j], &proc); /* determine the chunk of the boxman_entries[j] box that is needed */ hypre_BoxManEntryGetExtents(boxman_entries[j], ilower, iupper); hypre_BoxSetExtents(&boxman_entry_box, ilower, iupper); hypre_IntersectBoxes(&boxman_entry_box, &cbox, &boxman_entry_box); if (proc != myproc) { send_processes[i][cnt] = proc; hypre_SStructBoxManEntryGetBoxnum(boxman_entries[j], &send_remote_boxnums[i][cnt]); hypre_AppendBox(&boxman_entry_box, hypre_BoxArrayArrayBoxArray(send_boxes, i)); cnt++; } } hypre_TFree(boxman_entries); } /* hypre_ForBoxI(i, grid_boxes) */
HYPRE_Int HYPRE_SStructSplitSetup( HYPRE_SStructSolver solver, HYPRE_SStructMatrix A, HYPRE_SStructVector b, HYPRE_SStructVector x ) { hypre_SStructVector *y; HYPRE_Int nparts; HYPRE_Int *nvars; void ****smatvec_data; HYPRE_Int (***ssolver_solve)(); HYPRE_Int (***ssolver_destroy)(); void ***ssolver_data; HYPRE_Int ssolver = (solver -> ssolver); MPI_Comm comm; hypre_SStructGrid *grid; hypre_SStructPMatrix *pA; hypre_SStructPVector *px; hypre_SStructPVector *py; hypre_StructMatrix *sA; hypre_StructVector *sx; hypre_StructVector *sy; HYPRE_StructMatrix sAH; HYPRE_StructVector sxH; HYPRE_StructVector syH; HYPRE_Int (*ssolve)(); HYPRE_Int (*sdestroy)(); void *sdata; HYPRE_Int part, vi, vj; comm = hypre_SStructVectorComm(b); grid = hypre_SStructVectorGrid(b); HYPRE_SStructVectorCreate(comm, grid, &y); HYPRE_SStructVectorInitialize(y); HYPRE_SStructVectorAssemble(y); nparts = hypre_SStructMatrixNParts(A); nvars = hypre_TAlloc(HYPRE_Int, nparts); smatvec_data = hypre_TAlloc(void ***, nparts); ssolver_solve = (HYPRE_Int (***)()) hypre_MAlloc((sizeof(HYPRE_Int (**)()) * nparts)); ssolver_destroy = (HYPRE_Int (***)()) hypre_MAlloc((sizeof(HYPRE_Int (**)()) * nparts)); ssolver_data = hypre_TAlloc(void **, nparts); for (part = 0; part < nparts; part++) { pA = hypre_SStructMatrixPMatrix(A, part); px = hypre_SStructVectorPVector(x, part); py = hypre_SStructVectorPVector(y, part); nvars[part] = hypre_SStructPMatrixNVars(pA); smatvec_data[part] = hypre_TAlloc(void **, nvars[part]); ssolver_solve[part] = (HYPRE_Int (**)()) hypre_MAlloc((sizeof(HYPRE_Int (*)()) * nvars[part])); ssolver_destroy[part] = (HYPRE_Int (**)()) hypre_MAlloc((sizeof(HYPRE_Int (*)()) * nvars[part])); ssolver_data[part] = hypre_TAlloc(void *, nvars[part]); for (vi = 0; vi < nvars[part]; vi++) { smatvec_data[part][vi] = hypre_TAlloc(void *, nvars[part]); for (vj = 0; vj < nvars[part]; vj++) { sA = hypre_SStructPMatrixSMatrix(pA, vi, vj); sx = hypre_SStructPVectorSVector(px, vj); smatvec_data[part][vi][vj] = NULL; if (sA != NULL) { smatvec_data[part][vi][vj] = hypre_StructMatvecCreate(); hypre_StructMatvecSetup(smatvec_data[part][vi][vj], sA, sx); } } sA = hypre_SStructPMatrixSMatrix(pA, vi, vi); sx = hypre_SStructPVectorSVector(px, vi); sy = hypre_SStructPVectorSVector(py, vi); sAH = (HYPRE_StructMatrix) sA; sxH = (HYPRE_StructVector) sx; syH = (HYPRE_StructVector) sy; switch(ssolver) { default: /* If no solver is matched, use Jacobi, but throw and error */ if (ssolver != HYPRE_Jacobi) { hypre_error(HYPRE_ERROR_GENERIC); } /* don't break */ case HYPRE_Jacobi: HYPRE_StructJacobiCreate(comm, (HYPRE_StructSolver *)&sdata); HYPRE_StructJacobiSetMaxIter(sdata, 1); HYPRE_StructJacobiSetTol(sdata, 0.0); if (solver -> zero_guess) { HYPRE_StructJacobiSetZeroGuess(sdata); } HYPRE_StructJacobiSetup(sdata, sAH, syH, sxH); ssolve = HYPRE_StructJacobiSolve; sdestroy = HYPRE_StructJacobiDestroy; break; case HYPRE_SMG: HYPRE_StructSMGCreate(comm, (HYPRE_StructSolver *)&sdata); HYPRE_StructSMGSetMemoryUse(sdata, 0); HYPRE_StructSMGSetMaxIter(sdata, 1); HYPRE_StructSMGSetTol(sdata, 0.0); if (solver -> zero_guess) { HYPRE_StructSMGSetZeroGuess(sdata); } HYPRE_StructSMGSetNumPreRelax(sdata, 1); HYPRE_StructSMGSetNumPostRelax(sdata, 1); HYPRE_StructSMGSetLogging(sdata, 0); HYPRE_StructSMGSetPrintLevel(sdata, 0); HYPRE_StructSMGSetup(sdata, sAH, syH, sxH); ssolve = HYPRE_StructSMGSolve; sdestroy = HYPRE_StructSMGDestroy; break; case HYPRE_PFMG: HYPRE_StructPFMGCreate(comm, (HYPRE_StructSolver *)&sdata); HYPRE_StructPFMGSetMaxIter(sdata, 1); HYPRE_StructPFMGSetTol(sdata, 0.0); if (solver -> zero_guess) { HYPRE_StructPFMGSetZeroGuess(sdata); } HYPRE_StructPFMGSetRelaxType(sdata, 1); HYPRE_StructPFMGSetNumPreRelax(sdata, 1); HYPRE_StructPFMGSetNumPostRelax(sdata, 1); HYPRE_StructPFMGSetLogging(sdata, 0); HYPRE_StructPFMGSetPrintLevel(sdata, 0); HYPRE_StructPFMGSetup(sdata, sAH, syH, sxH); ssolve = HYPRE_StructPFMGSolve; sdestroy = HYPRE_StructPFMGDestroy; break; } ssolver_solve[part][vi] = ssolve; ssolver_destroy[part][vi] = sdestroy; ssolver_data[part][vi] = sdata; } } (solver -> y) = y; (solver -> nparts) = nparts; (solver -> nvars) = nvars; (solver -> smatvec_data) = smatvec_data; (solver -> ssolver_solve) = ssolver_solve; (solver -> ssolver_destroy) = ssolver_destroy; (solver -> ssolver_data) = ssolver_data; if ((solver -> tol) > 0.0) { hypre_SStructMatvecCreate(&(solver -> matvec_data)); hypre_SStructMatvecSetup((solver -> matvec_data), A, x); } return hypre_error_flag; }