152int main(
int argc,
char *argv[])
155 MPI_Init(&argc,&argv);
161 Galeri::core::Workspace::setNumDimensions(2);
163 Galeri::grid::Loadable domain, boundary;
165 int numGlobalElementsX = 4 * comm.NumProc();
166 int numGlobalElementsY = 4;
168 int mx = comm.NumProc();
171 Galeri::grid::Generator::
172 getSquareWithTriangles(comm, numGlobalElementsX, numGlobalElementsY,
173 mx, my, domain, boundary);
186 Galeri::problem::VectorLaplacian<MyVectorLaplacian> problem(numPDEs,
"Triangle");
188 Epetra_BlockMap matrixMap(domain.getNumGlobalVertices(), numPDEs, 0, comm);
190 problem.createGraph(domain, matrixGraph);
196 problem.integrate(domain, A,
RHS);
200 problem.imposeDirichletBoundaryConditions(boundary, A,
RHS, LHS);
208 AztecOO solver(linearProblem);
209 solver.SetAztecOption(AZ_solver, AZ_cg);
210 solver.SetAztecOption(AZ_precond, AZ_Jacobi);
211 solver.SetAztecOption(AZ_subdomain_solve, AZ_icc);
212 solver.SetAztecOption(AZ_output, 16);
214 solver.Iterate(1550, 1e-9);
216 Epetra_MultiVector* LHSComponent = Galeri::core::Workspace::createMultiVectorComponent(LHS);
218 for (
int ieq = 0; ieq < numPDEs; ++ieq)
220 Galeri::core::Workspace::extractMultiVectorComponent(LHS, ieq, *LHSComponent);
223 sprintf(fileName,
"sol%d", ieq);
224 Galeri::viz::MEDIT::write(domain, fileName, *LHSComponent);
226 problem.setEquation(ieq);
227 problem.computeNorms(domain, *LHSComponent);
static double getElementLHS(const double &x, const double &y, const double &z, const int &ieq, const int &jeq, const double &phi_i, const double &phi_i_derx, const double &phi_i_dery, const double &phi_i_derz, const double &phi_j, const double &phi_j_derx, const double &phi_j_dery, const double &phi_j_derz)