Skip to content

Commit

Permalink
Fixes to better handle normalization
Browse files Browse the repository at this point in the history
  • Loading branch information
@hughcars
hughcars committed Oct 3, 2023
1 parent bf73494 commit a97349c
Show file tree
Hide file tree
Showing 10 changed files with 64 additions and 39 deletions.
11 changes: 6 additions & 5 deletions palace/drivers/basesolver.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -567,7 +567,8 @@ void BaseSolver::PostprocessFields(const PostOperator &postop, int step, double
}

void BaseSolver::PostprocessErrorIndicators(const std::string &name, int step, double time,
const ErrorIndicators &indicators) const
const ErrorIndicators &indicators,
bool normalized) const
{
if (post_dir.length() == 0)
{
Expand All @@ -584,10 +585,10 @@ void BaseSolver::PostprocessErrorIndicators(const std::string &name, int step, d
// clang-format off
output.print("{:>{}s},{:>{}s},{:>{}s},{:>{}s},{:>{}s},{:>{}s}\n",
name, table.w1,
"Global Rel.", table.w,
"Element Min.", table.w,
"Element Max.", table.w,
"Element Mean", table.w,
normalized ? "Sum Rel." : "Sum Abs.", table.w,
normalized ? "Min Elem Rel." : "Min Elem Abs.", table.w,
normalized ? "Max Elem Rel." : "Max Elem Abs.", table.w,
normalized ? "Mean Elem Rel." : "Mean Elem Abs.", table.w,
"Normalization", table.w);
// clang-format on
}
Expand Down
7 changes: 4 additions & 3 deletions palace/drivers/basesolver.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -67,10 +67,11 @@ class BaseSolver
double time) const;
void PostprocessFields(const PostOperator &postop, int step, double time) const;

// Postprocess granular error indicator to file.
// Postprocess granular error indicator to file. The argument normalized indicates if
// supplied indicators have already been normalized.
void PostprocessErrorIndicators(const std::string &name, int step, double time,
const ErrorIndicators &indicators) const;
// Append combined error indicator to file.
const ErrorIndicators &indicators, bool normalized) const;
// Write a string labeled error indicator. Used for writing statistics.
void PostprocessErrorIndicators(const std::string &name,
const ErrorIndicators &indicators) const;

Expand Down
12 changes: 8 additions & 4 deletions palace/drivers/drivensolver.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -151,13 +151,15 @@ ErrorIndicators DrivenSolver::SweepUniform(SpaceOperator &spaceop, PostOperator
&postop](const auto &E, int step, double f)
{
BlockTimer bt0(Timer::ESTSOLVE);
auto estimate = estimator(E);
constexpr bool normalized = true;
auto estimate = estimator(E, normalized);
BlockTimer bt1(Timer::POSTPRO);
// Write the indicator for this mode.
postop.SetIndicatorGridFunction(estimate.indicators);
PostprocessErrorIndicators(
"f (GHz)", step, f,
ErrorIndicators{estimate, indicators.GlobalTrueVSize(), indicators.GetComm()});
ErrorIndicators{estimate, indicators.GlobalTrueVSize(), indicators.GetComm()},
normalized);
ErrorReducer(indicators, std::move(estimate));
};

Expand Down Expand Up @@ -285,12 +287,14 @@ ErrorIndicators DrivenSolver::SweepAdaptive(SpaceOperator &spaceop, PostOperator
&ErrorReducer](const auto &E, int step, double frequency)
{
BlockTimer bt0(Timer::ESTSOLVE);
auto estimate = estimator(E);
constexpr bool normalized = true;
auto estimate = estimator(E, normalized);
BlockTimer bt1(Timer::POSTPRO);
// Write the indicator for this mode.
PostprocessErrorIndicators(
"f (GHz)", step, frequency,
ErrorIndicators{estimate, indicators.GlobalTrueVSize(), indicators.GetComm()});
ErrorIndicators{estimate, indicators.GlobalTrueVSize(), indicators.GetComm()},
normalized);
ErrorReducer(indicators, std::move(estimate));
};

Expand Down
6 changes: 4 additions & 2 deletions palace/drivers/eigensolver.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -272,13 +272,15 @@ EigenSolver::Solve(const std::vector<std::unique_ptr<mfem::ParMesh>> &mesh) cons
[this, &estimator, &indicators, &ErrorReducer, &postop](const auto &E, int i)
{
BlockTimer bt0(Timer::ESTSOLVE);
auto ind = estimator(E);
constexpr bool normalized = true;
auto ind = estimator(E, normalized);
BlockTimer bt1(Timer::POSTPRO);
postop.SetIndicatorGridFunction(ind.indicators);
// Write the indicator for this mode.
PostprocessErrorIndicators(
"m", i, i + 1,
ErrorIndicators{ind, indicators.GlobalTrueVSize(), indicators.GetComm()});
ErrorIndicators{ind, indicators.GlobalTrueVSize(), indicators.GetComm()},
normalized);
ErrorReducer(indicators, std::move(ind));
};

Expand Down
6 changes: 4 additions & 2 deletions palace/drivers/electrostaticsolver.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -108,14 +108,16 @@ ErrorIndicators ElectrostaticSolver::Postprocess(LaplaceOperator &laplaceop,
&postop](const auto &V, int i, double idx)
{
BlockTimer bt0(Timer::ESTSOLVE);
auto estimate = estimator(V);
constexpr bool normalized = true;
auto estimate = estimator(V, normalized);
ErrorReducer(indicators, estimate);
BlockTimer bt1(Timer::POSTPRO);
// Write the indicator for this mode.
postop.SetIndicatorGridFunction(estimate.indicators);
PostprocessErrorIndicators(
"i", i, idx,
ErrorIndicators{estimate, indicators.GlobalTrueVSize(), indicators.GetComm()});
ErrorIndicators{estimate, indicators.GlobalTrueVSize(), indicators.GetComm()},
normalized);
ErrorReducer(indicators, std::move(estimate));
};
if (iodata.solver.electrostatic.n_post > 0)
Expand Down
6 changes: 4 additions & 2 deletions palace/drivers/magnetostaticsolver.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -110,13 +110,15 @@ ErrorIndicators MagnetostaticSolver::Postprocess(CurlCurlOperator &curlcurlop,
&postop](const auto &A, int i, double idx)
{
BlockTimer bt0(Timer::ESTSOLVE);
auto estimate = estimator(A);
constexpr bool normalized = true;
auto estimate = estimator(A, normalized);
BlockTimer bt1(Timer::POSTPRO);
// Write the indicator for this mode.
postop.SetIndicatorGridFunction(estimate.indicators);
PostprocessErrorIndicators(
"i", i, idx,
ErrorIndicators{estimate, indicators.GlobalTrueVSize(), indicators.GetComm()});
ErrorIndicators{estimate, indicators.GlobalTrueVSize(), indicators.GetComm()},
normalized);
ErrorReducer(indicators, std::move(estimate));
};
if (iodata.solver.magnetostatic.n_post > 0)
Expand Down
16 changes: 8 additions & 8 deletions palace/drivers/transientsolver.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -90,12 +90,15 @@ TransientSolver::Solve(const std::vector<std::unique_ptr<mfem::ParMesh>> &mesh)
&postop](const auto &E, int step, double time)
{
BlockTimer bt0(Timer::ESTSOLVE);
auto estimate = estimator(E);
// Initial flux of zero would return nan.
bool constexpr normalized = false;
auto estimate = estimator(E, normalized);
BlockTimer bt1(Timer::POSTPRO);
postop.SetIndicatorGridFunction(estimate.indicators);
PostprocessErrorIndicators(
"t (ns)", step, time,
ErrorIndicators{estimate, indicators.GlobalTrueVSize(), indicators.GetComm()});
ErrorIndicators{estimate, indicators.GlobalTrueVSize(), indicators.GetComm()},
normalized);
ErrorReducer(indicators, std::move(estimate));
};

Expand Down Expand Up @@ -140,11 +143,8 @@ TransientSolver::Solve(const std::vector<std::unique_ptr<mfem::ParMesh>> &mesh)
E_elec + E_mag);
}

if (step > 0)
{
// Calculate and record the error indicators.
UpdateErrorIndicators(E, step, t);
}
// Calculate and record the error indicators.
UpdateErrorIndicators(E, step, t);

// Postprocess port voltages/currents and optionally write solution to disk.
Postprocess(postop, spaceop.GetLumpedPortOp(), spaceop.GetSurfaceCurrentOp(), step, t,
Expand All @@ -154,7 +154,7 @@ TransientSolver::Solve(const std::vector<std::unique_ptr<mfem::ParMesh>> &mesh)
step++;
}
SaveMetadata(timeop.GetLinearSolver());
PostprocessErrorIndicators("Mean", step, t, indicators);
PostprocessErrorIndicators("Mean", indicators);
return indicators;
}
std::function<double(double)> TransientSolver::GetTimeExcitation(bool dot) const
Expand Down
30 changes: 21 additions & 9 deletions palace/linalg/errorestimator.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -47,7 +47,8 @@ CurlFluxErrorEstimator::CurlFluxErrorEstimator(
}
}

IndicatorsAndNormalization CurlFluxErrorEstimator::operator()(const ComplexVector &v) const
IndicatorsAndNormalization CurlFluxErrorEstimator::operator()(const ComplexVector &v,
bool normalize) const
{
mfem::ParComplexGridFunction field(&fes);
field.real().SetFromTrueDofs(v.Real());
Expand Down Expand Up @@ -152,12 +153,16 @@ IndicatorsAndNormalization CurlFluxErrorEstimator::operator()(const ComplexVecto

Mpi::GlobalSum(1, &normalization, field.ParFESpace()->GetComm());
normalization = std::sqrt(normalization);
std::for_each(estimates.begin(), estimates.end(),
[&normalization](auto &x) { x /= normalization; });
if (normalize)
{
std::for_each(estimates.begin(), estimates.end(),
[&normalization](auto &x) { x /= normalization; });
}
return {estimates, normalization};
}

IndicatorsAndNormalization CurlFluxErrorEstimator::operator()(const Vector &v) const
IndicatorsAndNormalization CurlFluxErrorEstimator::operator()(const Vector &v,
bool normalize) const
{
mfem::ParGridFunction field(&fes);
field.SetFromTrueDofs(v);
Expand Down Expand Up @@ -236,8 +241,11 @@ IndicatorsAndNormalization CurlFluxErrorEstimator::operator()(const Vector &v) c

Mpi::GlobalSum(1, &normalization, field.ParFESpace()->GetComm());
normalization = std::sqrt(normalization);
std::for_each(estimates.begin(), estimates.end(),
[&normalization](auto &x) { x /= normalization; });
if (normalize)
{
std::for_each(estimates.begin(), estimates.end(),
[&normalization](auto &x) { x /= normalization; });
}
return {estimates, normalization};
}

Expand Down Expand Up @@ -276,7 +284,8 @@ GradFluxErrorEstimator::GradFluxErrorEstimator(
}
}

IndicatorsAndNormalization GradFluxErrorEstimator::operator()(const Vector &v) const
IndicatorsAndNormalization GradFluxErrorEstimator::operator()(const Vector &v,
bool normalize) const
{
mfem::ParGridFunction field(&fes);
field.SetFromTrueDofs(v);
Expand Down Expand Up @@ -388,8 +397,11 @@ IndicatorsAndNormalization GradFluxErrorEstimator::operator()(const Vector &v) c

paraview.Save();
}
std::for_each(estimates.begin(), estimates.end(),
[&normalization](auto &x) { x /= normalization; });
if (normalize)
{
std::for_each(estimates.begin(), estimates.end(),
[&normalization](auto &x) { x /= normalization; });
}
return {estimates, normalization};
}
} // namespace palace
7 changes: 4 additions & 3 deletions palace/linalg/errorestimator.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -46,10 +46,11 @@ class CurlFluxErrorEstimator
mfem::ParFiniteElementSpace &fes);

// Compute elemental error indicators given a complex vector of true DOF.
IndicatorsAndNormalization operator()(const ComplexVector &v) const;
IndicatorsAndNormalization operator()(const ComplexVector &v,
bool normalize = true) const;

// Compute elemental error indicators given a vector of true DOF.
IndicatorsAndNormalization operator()(const Vector &v) const;
IndicatorsAndNormalization operator()(const Vector &v, bool normalize = true) const;
};

// Class used for computing grad flux error estimate,
Expand Down Expand Up @@ -82,7 +83,7 @@ class GradFluxErrorEstimator
mfem::ParFiniteElementSpace &fes);

// Compute elemental error indicators given a vector of true DOF.
IndicatorsAndNormalization operator()(const Vector &v) const;
IndicatorsAndNormalization operator()(const Vector &v, bool normalize = true) const;
};

} // namespace palace
Expand Down
2 changes: 1 addition & 1 deletion palace/main.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -188,7 +188,7 @@ int main(int argc, char *argv[])
// Run the problem driver.
auto solver_output = solver->Solve(mesh);

Mpi::Print(world_comm, "Normalized Error Indicator: {:.3e}",
Mpi::Print(world_comm, "Error Estimate: {:.3e}\n",
solver_output.GetGlobalErrorIndicator());

// Print timing summary.
Expand Down

0 comments on commit a97349c

Please sign in to comment.