From 7f79eaa6958bc4fb9bcc0a303b65c349f8e00e86 Mon Sep 17 00:00:00 2001 From: Sebastian Aland <sebastian.aland@yahoo.com> Date: Thu, 21 Feb 2013 13:35:12 +0000 Subject: [PATCH] PetscSolverNSCH is now fit for 3D --- AMDiS/src/parallel/PetscSolverNSCH.cc | 35 +++++++++++++-------------- AMDiS/src/parallel/PetscSolverNSCH.h | 1 + 2 files changed, 18 insertions(+), 18 deletions(-) diff --git a/AMDiS/src/parallel/PetscSolverNSCH.cc b/AMDiS/src/parallel/PetscSolverNSCH.cc index e97da33a..607b24e0 100644 --- a/AMDiS/src/parallel/PetscSolverNSCH.cc +++ b/AMDiS/src/parallel/PetscSolverNSCH.cc @@ -31,18 +31,16 @@ namespace AMDiS { CahnHilliardData2* data = static_cast<CahnHilliardData2*>(ctx); /// extract vectors - Vec b1, b2, b34, b5, x1, x2, x34, x5, b3,b4; - data->globalMatrixSolver->extractVectorComponent(b, 0+3, &b1); - data->globalMatrixSolver->extractVectorComponent(b, 1+3, &b2); - data->globalMatrixSolver->extractVectorComponent(b, 0, &b34, 2); - data->globalMatrixSolver->extractVectorComponent(b, 2, &b5); + Vec b1, b2, b34, b5, x1, x2, x34, x5; + data->globalMatrixSolver->extractVectorComponent(b, data->dim+1, &b1); + data->globalMatrixSolver->extractVectorComponent(b, data->dim+2, &b2); + data->globalMatrixSolver->extractVectorComponent(b, 0, &b34, data->dim); + data->globalMatrixSolver->extractVectorComponent(b, data->dim, &b5); - data->globalMatrixSolver->extractVectorComponent(x, 0+3, &x1); - data->globalMatrixSolver->extractVectorComponent(x, 1+3, &x2); - data->globalMatrixSolver->extractVectorComponent(x, 0, &x34, 2); - data->globalMatrixSolver->extractVectorComponent(x, 2, &x5); - data->globalMatrixSolver->extractVectorComponent(b, 0, &b3, 1); - data->globalMatrixSolver->extractVectorComponent(b, 1, &b4, 1); + data->globalMatrixSolver->extractVectorComponent(x, data->dim+1, &x1); + data->globalMatrixSolver->extractVectorComponent(x, data->dim+2, &x2); + data->globalMatrixSolver->extractVectorComponent(x, 0, &x34, data->dim); + data->globalMatrixSolver->extractVectorComponent(x, data->dim, &x5); /// solve Cahn-Hilliard Preconditioner @@ -190,7 +188,7 @@ namespace AMDiS { double wtime = MPI::Wtime(); int dim = componentSpaces[0]->getMesh()->getDim(); pressureComponent=dim; - const FiniteElemSpace *cahnHilliardFeSpace = componentSpaces[0+3]; + const FiniteElemSpace *cahnHilliardFeSpace = componentSpaces[dim+1]; const FiniteElemSpace *velocityFeSpace= componentSpaces[0]; const FiniteElemSpace *pressureFeSpace = componentSpaces[pressureComponent]; @@ -199,6 +197,7 @@ namespace AMDiS { PCShellSetContext(pc, &matShellContext); matShellContext.globalMatrixSolver = (this); matShellContext.mpiCommGlobal= &(meshDistributor->getMpiComm(0)); + matShellContext.dim = dim; /// Init Cahn-Hilliard part @@ -222,16 +221,16 @@ namespace AMDiS { laplaceOp2.addTerm(&sot); // -delta*K massMatrixCH.assembleOperator(massOpCH); - massMatrixSolverCH = createSubSolver(0+3, "mass_"); + massMatrixSolverCH = createSubSolver(dim+1, "mass_"); massMatrixSolverCH->fillPetscMatrix(&massMatrixCH); // === matrix (M + eps*sqrt(delta)*K) === laplaceMatrixCH.assembleOperator(laplaceOpCH); - laplaceMatrixSolverCH = createSubSolver(0+3, "laplace_"); + laplaceMatrixSolverCH = createSubSolver(dim+1, "laplace_"); laplaceMatrixSolverCH->fillPetscMatrix(&laplaceMatrixCH); // === matrix (-delta*K) === deltaKMatrix.assembleOperator(laplaceOp2); - deltaKMatrixSolver = createSubSolver(0+3, "laplace2_"); + deltaKMatrixSolver = createSubSolver(dim+1, "laplace2_"); deltaKMatrixSolver->fillPetscMatrix(&deltaKMatrix); } @@ -317,9 +316,9 @@ namespace AMDiS { conDifMatrixSolver->fillPetscMatrix(&conDifMatrix); matShellContext.matConDif = conDifMatrixSolver->getMatInterior(); - extractMatrixComponent(mat[0][0], 2, 1, 0, 2, &(matShellContext.matDiv)); - extractMatrixComponent(mat[0][0], 0, 2, 2, 1, &(matShellContext.matGrad)); - extractMatrixComponent(mat[0][0], 0, 2, 0, 2, &(matShellContext.velocityMat)); + extractMatrixComponent(mat[0][0], dim, 1, 0, dim, &(matShellContext.matDiv)); + extractMatrixComponent(mat[0][0], 0, dim, dim, 1, &(matShellContext.matGrad)); + extractMatrixComponent(mat[0][0], 0, dim, 0, dim, &(matShellContext.velocityMat)); ///erstelle kspVelocity KSPCreate((meshDistributor->getMpiComm(0)), &(matShellContext.kspVelocity)); diff --git a/AMDiS/src/parallel/PetscSolverNSCH.h b/AMDiS/src/parallel/PetscSolverNSCH.h index c5a0b926..31187869 100644 --- a/AMDiS/src/parallel/PetscSolverNSCH.h +++ b/AMDiS/src/parallel/PetscSolverNSCH.h @@ -31,6 +31,7 @@ namespace AMDiS { Mat matMassCH, matMinusDeltaK, matGrad, matDiv, matConDif, matSchur, velocityMat; PetscSolverGlobalMatrix *globalMatrixSolver; double *eps, *delta; + int dim; MPI::Intracomm *mpiCommGlobal; }; -- GitLab