Skip to content

Commit

Permalink
Eliminate some dead _GPU_ code.
Browse files Browse the repository at this point in the history
Specifically, parts of GradNonLinearForm and GradFaceIntegrator had
_GPU_ specific code that was superseded by device-specific code in the
Gradients class and thus was never called.  This code has been
eliminated to improve hygiene in preparation for changes to address #198.
  • Loading branch information
trevilo committed May 9, 2023
1 parent 958d8d1 commit cec7f03
Show file tree
Hide file tree
Showing 3 changed files with 12 additions and 103 deletions.
31 changes: 7 additions & 24 deletions src/faceGradientIntegration.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -60,18 +60,17 @@ void GradFaceIntegrator::AssembleFaceVector(const FiniteElement &el1, const Fini
elvect.SetSize((dof1 + dof2) * num_equation * dim);
elvect = 0.0;

// #ifdef _GPU_
// DenseMatrix elfun1_mat(elfun.GetData(), dof1, num_equation);
// DenseMatrix elfun2_mat(elfun.GetData()+dof1*num_equation,
// dof2, num_equation);
#ifndef _GPU_
DenseMatrix elfun1_mat(elfun.GetData(), dof1, num_equation * dim);
DenseMatrix elfun2_mat(elfun.GetData() + dof1 * num_equation * dim, dof2, num_equation * dim);
// NB: Incoming Vector &elfun is in gradient space. The first
// component is the state, and the rest is zero. See
// GradNonLinearForm::Mult in gradNonLinearForm.cpp for details.
DenseMatrix elfun1_mat(elfun.GetData(), dof1, num_equation);
DenseMatrix elfun2_mat(elfun.GetData() + dof1 * num_equation * dim, dof2, num_equation);


DenseMatrix elvect1_mat(elvect.GetData(), dof1, num_equation * dim);
DenseMatrix elvect2_mat(elvect.GetData() + dof1 * num_equation * dim, dof2, num_equation * dim);
elvect1_mat = 0.;
elvect2_mat = 0.;
#endif

// Integration order calculation from DGTraceIntegrator
int intorder;
Expand Down Expand Up @@ -104,20 +103,12 @@ void GradFaceIntegrator::AssembleFaceVector(const FiniteElement &el1, const Fini
for (int eq = 0; eq < num_equation; eq++) {
double sum = 0.;
for (int k = 0; k < dof1; k++) {
#ifdef _GPU_
sum += shape1(k) * elfun(k + eq * dof1);
#else
sum += shape1[k] * elfun1_mat(k, eq);
#endif
}
iUp1(eq) = sum;
sum = 0.;
for (int k = 0; k < dof2; k++) {
#ifdef _GPU_
sum += shape2(k) * elfun(k + eq * dof2 + dof1 * num_equation);
#else
sum += shape2[k] * elfun2_mat(k, eq);
#endif
}
iUp2(eq) = sum;
mean(eq) = (iUp1(eq) + iUp2(eq)) / 2.;
Expand All @@ -133,18 +124,10 @@ void GradFaceIntegrator::AssembleFaceVector(const FiniteElement &el1, const Fini
const double du1n = (mean(eq) - iUp1(eq)) * nor(d);
const double du2n = (mean(eq) - iUp2(eq)) * nor(d);
for (int k = 0; k < dof1; k++) {
#ifdef _GPU_
elvect(k + eq + d * num_equation) += du1n * shape1(k);
#else
elvect1_mat(k, eq + d * num_equation) += du1n * shape1(k);
#endif
}
for (int k = 0; k < dof2; k++) {
#ifdef _GPU_
elvect(k + eq + d * num_equation + dof1 * num_equation * dim) -= du2n * shape2(k);
#else
elvect2_mat(k, eq + d * num_equation) -= du2n * shape2(k);
#endif
}
}
}
Expand Down
78 changes: 5 additions & 73 deletions src/gradNonLinearForm.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -39,82 +39,14 @@ GradNonLinearForm::GradNonLinearForm(ParFiniteElementSpace *_vfes, IntegrationRu

void GradNonLinearForm::Mult(const ParGridFunction *Up, Vector &y) {
Vector x;

#ifdef _GPU_
// TestNonLinearForm::setToZero_gpu(y,y.Size()); already done in gradient.cpd
x.UseDevice(true);
if (fnfi.Size()) {
x.NewMemoryAndSize(Up->GetMemory(), vfes->GetNDofs() * num_equation, false);
Mult_gpu(x, y);
}
#else
const double *data = Up->GetData();

// NB: Setting x here accounts for the fact that the space Up lives
// in is different from the space of the gradient.
x.SetSize(vfes->GetNDofs() * num_equation * dim);
x = 0.;
for (int i = 0; i < vfes->GetNDofs() * num_equation; i++) x(i) = data[i];
ParNonlinearForm::Mult(x, y);
#endif
}

#ifdef _GPU_
void GradNonLinearForm::Mult_gpu(const Vector &x, Vector &y) {
// INTEGRATION SHARED FACES
const Vector &px = Prolongate(x);

const int dofs = vfes->GetNDofs();

// Terms over shared interior faces in parallel.
ParFiniteElementSpace *pfes = ParFESpace();
ParMesh *pmesh = pfes->GetParMesh();
FaceElementTransformations *tr;
const FiniteElement *fe1, *fe2;
Array<int> vdofs1, vdofs2;
Vector el_x, el_y;
el_x.UseDevice(true);
el_y.UseDevice(true);

// if (bfnfi.Size()==0) y.HostReadWrite();
X.MakeRef(aux1, 0); // aux1 contains P.x
X.ExchangeFaceNbrData();
const int n_shared_faces = pmesh->GetNSharedFaces();
for (int i = 0; i < n_shared_faces; i++) {
tr = pmesh->GetSharedFaceTransformations(i, true);
int Elem2NbrNo = tr->Elem2No - pmesh->GetNE();

fe1 = pfes->GetFE(tr->Elem1No);
fe2 = pfes->GetFaceNbrFE(Elem2NbrNo);

pfes->GetElementVDofs(tr->Elem1No, vdofs1);
pfes->GetFaceNbrElementVDofs(Elem2NbrNo, vdofs2);

el_x.SetSize(vdofs1.Size() + vdofs2.Size());
X.GetSubVector(vdofs1, el_x.GetData());
X.FaceNbrData().GetSubVector(vdofs2, el_x.GetData() + vdofs1.Size());

const int elDof1 = fe1->GetDof();
const int elDof2 = fe2->GetDof();
for (int k = 0; k < fnfi.Size(); k++) {
fnfi[k]->AssembleFaceVector(*fe1, *fe2, *tr, el_x, el_y);
// y.AddElementVector(vdofs1, el_y.GetData());
GradNonLinearForm::IndexedAddToGlobalMemory(vdofs1, y, el_y, dofs, elDof1, num_equation, dim);
}
}
}

void GradNonLinearForm::IndexedAddToGlobalMemory(const Array<int> &vdofs, Vector &y, const Vector &el_y, const int &dof,
const int &elDof, const int &num_equation, const int &dim) {
double *dy = y.ReadWrite();
const double *d_ely = el_y.Read();
auto d_vdofs = vdofs.Read();

MFEM_FORALL(i, elDof, {
const int index = d_vdofs[i];
for (int eq = 0; eq < num_equation; eq++) {
for (int d = 0; d < dim; d++) {
dy[index + eq * dof + d * num_equation * dof] += d_ely[i + eq * elDof + d * num_equation * elDof];
}
}
});
// After setting x as above, base class Mult does the right thing
ParNonlinearForm::Mult(x, y);
}

#endif // _GPU_
6 changes: 0 additions & 6 deletions src/gradNonLinearForm.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -51,12 +51,6 @@ class GradNonLinearForm : public ParNonlinearForm {
GradNonLinearForm(ParFiniteElementSpace *f, IntegrationRules *intRules, const int dim, const int num_equation);

void Mult(const ParGridFunction *Up, Vector &y);

#ifdef _GPU_
void Mult_gpu(const Vector &x, Vector &y);
static void IndexedAddToGlobalMemory(const Array<int> &vdofs, Vector &y, const Vector &el_y, const int &dof,
const int &elDof, const int &num_equation, const int &dim);
#endif
};

#endif // GRADNONLINEARFORM_HPP_

0 comments on commit cec7f03

Please sign in to comment.