From 065dd2d18ddeb15420bb81777c502a64b1922236 Mon Sep 17 00:00:00 2001 From: Hugh Carson Date: Fri, 21 Jul 2023 16:17:25 -0400 Subject: [PATCH 1/2] Fixing a range of compiler warnings for shadow/unused variables - Most of the OOP/functional clashes were resolved by static where possible, or a free function - SurfacePostOperator involved breaking out private struct definitions into the palace namespace these struct definitions had already bled into the global namespace through the use of accessor functions. - Source of the issue being methods of internal structs believing mat_op and local_to_shared were mirroring the containing classes variables, despite this being impossible. --- palace/drivers/eigensolver.cpp | 5 +- palace/fem/coefficient.hpp | 11 ++-- palace/linalg/feast.cpp | 30 ++++----- palace/linalg/feast.hpp | 4 +- palace/linalg/ksp.cpp | 4 +- palace/linalg/superlu.cpp | 2 +- palace/linalg/superlu.hpp | 1 - palace/models/romoperator.cpp | 2 + palace/models/romoperator.hpp | 4 +- palace/models/surfacepostoperator.cpp | 56 +++++++++-------- palace/models/surfacepostoperator.hpp | 88 +++++++++++++-------------- palace/models/waveportoperator.cpp | 53 ++++++++-------- palace/models/waveportoperator.hpp | 11 ++-- 13 files changed, 138 insertions(+), 133 deletions(-) diff --git a/palace/drivers/eigensolver.cpp b/palace/drivers/eigensolver.cpp index 11369c158..bb16100ae 100644 --- a/palace/drivers/eigensolver.cpp +++ b/palace/drivers/eigensolver.cpp @@ -164,7 +164,7 @@ void EigenSolver::Solve(std::vector> &mesh, std::unique_ptr ksp; std::unique_ptr pc; #if defined(PALACE_WITH_SLEPC) - auto *feast = dynamic_cast(eigen.get()); + auto * const feast = dynamic_cast(eigen.get()); if (feast) { // Configure the FEAST integration contour. The linear solvers are set up inside the @@ -339,8 +339,7 @@ void EigenSolver::Solve(std::vector> &mesh, #if defined(PALACE_WITH_SLEPC) if (!ksp) { - const auto &feast = dynamic_cast(*eigen); - SaveMetadata(feast.GetTotalKspMult(), feast.GetTotalKspIter()); + SaveMetadata(feast->GetTotalKspMult(), feast->GetTotalKspIter()); } else #endif diff --git a/palace/fem/coefficient.hpp b/palace/fem/coefficient.hpp index 18b53de55..8119f0759 100644 --- a/palace/fem/coefficient.hpp +++ b/palace/fem/coefficient.hpp @@ -217,8 +217,7 @@ class DielectricInterfaceCoefficient : public mfem::Coefficient, const mfem::Vector side; mutable mfem::Vector C1, V, nor; - int Initialize(mfem::ElementTransformation &T, const mfem::IntegrationPoint &ip, - mfem::Vector &V) + int Initialize(mfem::ElementTransformation &T, const mfem::IntegrationPoint &ip) { // Get neighboring elements. mfem::ElementTransformation *T1, *T2; @@ -277,7 +276,7 @@ inline double DielectricInterfaceCoefficient::Eval( mfem::ElementTransformation &T, const mfem::IntegrationPoint &ip) { // Get single-sided solution and neighboring element attribute. - Initialize(T, ip, V); + Initialize(T, ip); GetNormal(T, ip, nor); // Metal-air interface: 0.5 * t / ϵ_MA * |E_n|² . @@ -290,7 +289,7 @@ inline double DielectricInterfaceCoefficient::Eval( mfem::ElementTransformation &T, const mfem::IntegrationPoint &ip) { // Get single-sided solution and neighboring element attribute. - int attr = Initialize(T, ip, V); + int attr = Initialize(T, ip); GetNormal(T, ip, nor); // Metal-substrate interface: 0.5 * t * (ϵ_S)² / ϵ_MS * |E_n|² . @@ -304,7 +303,7 @@ inline double DielectricInterfaceCoefficient::Eval( mfem::ElementTransformation &T, const mfem::IntegrationPoint &ip) { // Get single-sided solution and neighboring element attribute. - Initialize(T, ip, V); + Initialize(T, ip); GetNormal(T, ip, nor); // Substrate-air interface: 0.5 * t * (ϵ_SA * |E_t|² + 1 / ϵ_MS * |E_n|²) . @@ -318,7 +317,7 @@ inline double DielectricInterfaceCoefficient:: mfem::ElementTransformation &T, const mfem::IntegrationPoint &ip) { // Get single-sided solution and neighboring element attribute. - Initialize(T, ip, V); + Initialize(T, ip); // No specific interface, use full field evaluation: 0.5 * t * ϵ * |E|² . return 0.5 * ts * epsilon * (V * V); diff --git a/palace/linalg/feast.cpp b/palace/linalg/feast.cpp index 24d1d7a6c..1008bfbde 100644 --- a/palace/linalg/feast.cpp +++ b/palace/linalg/feast.cpp @@ -140,7 +140,7 @@ class FeastLinearSolver Q.RestoreColumn(j, q); if (M > 1) { - petsc::PetscParVector q = Q.GetColumn(j + m0); + q = Q.GetColumn(j + m0); q.AXPY(zk / gamma, v); Q.RestoreColumn(j + m0, q); } @@ -182,8 +182,8 @@ FeastEigenSolver::FeastEigenSolver(MPI_Comm comm, const IoData &iodata, print = print_lvl; info = 0; nev = m0 = mQ = 0; - M = iodata.solver.eigenmode.feast_moments; - MFEM_VERIFY(M == 1 || M == 2, + N_moments = iodata.solver.eigenmode.feast_moments; + MFEM_VERIFY(N_moments == 1 || N_moments == 2, "FEAST eigensolver only supports up to 2 subspace moments!"); rtol = 0.0; max_it = 0; @@ -268,7 +268,7 @@ void FeastEigenSolver::SetNumModes(int numeig, int numvec) m0 = std::max(nev + 3, nev + (nev + 1) / 2); } } - mQ = 2 * M * m0; // Real-valued basis splitting leads to factor of 2 + mQ = 2 * N_moments * m0; // Real-valued basis splitting leads to factor of 2 } void FeastEigenSolver::SetTol(double tol) @@ -381,7 +381,7 @@ int FeastEigenSolver::SolveInternal(RG rg) for (PetscInt j = 1; j < mQ / 2; j++) { // Ensure homogeneous Dirichlet BC are satisfied by the starting subspace. - petsc::PetscParVector q = Q.GetColumn(j); + q = Q.GetColumn(j); q.PointwiseMult(*r0, false); Q.RestoreColumn(j, q); } @@ -458,27 +458,27 @@ int FeastEigenSolver::SolveInternal(RG rg) x.Normalize(); } GetResidual(sigma, x, r); - PetscReal res = r.Norml2() / (x.Norml2() * PetscAbsScalar(sigma)); + PetscReal local_res = r.Norml2() / (x.Norml2() * PetscAbsScalar(sigma)); // PetscReal res = r.Norml2()/x.Norml2(); X->RestoreColumn(j, x); R.RestoreColumn(j, r); - if (res < rtol) + if (local_res < rtol) { // Mark converged even for eigenvalues outside the contour. converged[j] = true; nconv++; - if (res > rmax) + if (local_res > rmax) { - rmax = res; + rmax = local_res; jmax = j; } } else { converged[j] = false; - if (res < rmin) + if (local_res < rmin) { - rmin = res; + rmin = local_res; jmin = j; } } @@ -493,7 +493,7 @@ int FeastEigenSolver::SolveInternal(RG rg) // Debug // Mpi::Print(comm, " res[{:d}] = {:e} (eig = {:+e}{:+e}i, inside = {:d})\n", - // j, res, PetscRealPart(sigma), + // j, local_res, PetscRealPart(sigma), // PetscImaginaryPart(sigma), inside[j]); } if (print > 0) @@ -1198,10 +1198,10 @@ void FeastPEPSolver::SolveProjectedProblem(const petsc::PetscDenseMatrix &Q_, for (PetscInt i = 0; i < m0; i++) { eig_[i] = alpha[sort[i]]; - const PetscScalar *pXQ = XQ->GetArrayRead(); + const PetscScalar* const local_pXQ = XQ->GetArrayRead(); PetscScalar *pXQ0 = XQ0->GetArray(); - PalacePetscCall(PetscArraycpy(pXQ0 + mQ * i, pXQ + 2 * mQ * sort[i], mQ)); - XQ->RestoreArrayRead(pXQ); + PalacePetscCall(PetscArraycpy(pXQ0 + mQ * i, local_pXQ + 2 * mQ * sort[i], mQ)); + XQ->RestoreArrayRead(local_pXQ); XQ0->RestoreArray(pXQ0); } diff --git a/palace/linalg/feast.hpp b/palace/linalg/feast.hpp index 63547c422..e4684e66d 100644 --- a/palace/linalg/feast.hpp +++ b/palace/linalg/feast.hpp @@ -55,7 +55,7 @@ class FeastEigenSolver : public EigenSolverBase PetscInt nev, m0, mQ; // Number of moments to consider for subspace construction. - PetscInt M; + PetscInt N_moments; // Relative eigenvalue error convergence tolerance for the solver. PetscReal rtol; @@ -192,7 +192,7 @@ class FeastEigenSolver : public EigenSolverBase { MFEM_ABORT("SetWhichEigenpairs not defined for FeastEigenSolver!"); } - void SetShiftInvert(double tr, double ti, bool precond = false) override + void SetShiftInvert(double tr_in, double ti_in, bool precond = false) override { MFEM_ABORT("SetShiftInvert not defined for FeastEigenSolver!"); } diff --git a/palace/linalg/ksp.cpp b/palace/linalg/ksp.cpp index 4a0946f72..65489f06c 100644 --- a/palace/linalg/ksp.cpp +++ b/palace/linalg/ksp.cpp @@ -104,8 +104,10 @@ void KspSolver::Configure(const IoData &iodata) } } -void KspSolver::ConfigureVerbose(int print, const std::string &prefix) +void KspSolver::ConfigureVerbose(int print_input, const std::string &prefix) { + print = print_input; + // Manage debugging output. if (print > 0) { diff --git a/palace/linalg/superlu.cpp b/palace/linalg/superlu.cpp index ba601c0bc..0355e44ab 100644 --- a/palace/linalg/superlu.cpp +++ b/palace/linalg/superlu.cpp @@ -37,7 +37,7 @@ int GetNpDep(int np, bool use_3d) SuperLUSolver::SuperLUSolver(MPI_Comm comm, config::LinearSolverData::SymFactType reorder, bool use_3d, int print) - : mfem::Solver(), comm(comm), A(nullptr), solver(comm, GetNpDep(Mpi::Size(comm), use_3d)) + : mfem::Solver(), A(nullptr), solver(comm, GetNpDep(Mpi::Size(comm), use_3d)) { // Configure the solver. if (print > 1) diff --git a/palace/linalg/superlu.hpp b/palace/linalg/superlu.hpp index 1daf86631..f5c52c1c4 100644 --- a/palace/linalg/superlu.hpp +++ b/palace/linalg/superlu.hpp @@ -20,7 +20,6 @@ namespace palace class SuperLUSolver : public mfem::Solver { private: - MPI_Comm comm; std::unique_ptr A; mfem::SuperLUSolver solver; diff --git a/palace/models/romoperator.cpp b/palace/models/romoperator.cpp index e369319f6..9c402dffd 100644 --- a/palace/models/romoperator.cpp +++ b/palace/models/romoperator.cpp @@ -467,6 +467,7 @@ double RomOperator::ComputeError(double omega) return num / den; } +// static void RomOperator::BVMatProjectInternal(petsc::PetscDenseMatrix &V, petsc::PetscParMatrix &A, petsc::PetscDenseMatrix &Ar, petsc::PetscParVector &r, int n0, int n) @@ -515,6 +516,7 @@ void RomOperator::BVMatProjectInternal(petsc::PetscDenseMatrix &V, petsc::PetscP } } +// static void RomOperator::BVDotVecInternal(petsc::PetscDenseMatrix &V, petsc::PetscParVector &b, petsc::PetscParVector &br, int n0, int n) { diff --git a/palace/models/romoperator.hpp b/palace/models/romoperator.hpp index 616d7ad4b..c72f986fa 100644 --- a/palace/models/romoperator.hpp +++ b/palace/models/romoperator.hpp @@ -65,10 +65,10 @@ class RomOperator double ComputeError(double omega); // Helper functions for reduced-order matrix or vector construction/update. - void BVMatProjectInternal(petsc::PetscDenseMatrix &V, petsc::PetscParMatrix &A, + static void BVMatProjectInternal(petsc::PetscDenseMatrix &V, petsc::PetscParMatrix &A, petsc::PetscDenseMatrix &Ar, petsc::PetscParVector &r, int n0, int n); - void BVDotVecInternal(petsc::PetscDenseMatrix &V, petsc::PetscParVector &b, + static void BVDotVecInternal(petsc::PetscDenseMatrix &V, petsc::PetscParVector &b, petsc::PetscParVector &br, int n0, int n); public: diff --git a/palace/models/surfacepostoperator.cpp b/palace/models/surfacepostoperator.cpp index 5ea49b447..c2b7b5a83 100644 --- a/palace/models/surfacepostoperator.cpp +++ b/palace/models/surfacepostoperator.cpp @@ -12,7 +12,7 @@ namespace palace { -SurfacePostOperator::InterfaceDielectricData::InterfaceDielectricData( +InterfaceDielectricData::InterfaceDielectricData( const config::InterfaceDielectricData &data, mfem::ParMesh &mesh) : ts(data.ts), tandelta(data.tandelta) { @@ -101,7 +101,7 @@ SurfacePostOperator::InterfaceDielectricData::InterfaceDielectricData( } std::unique_ptr -SurfacePostOperator::InterfaceDielectricData::GetCoefficient( +InterfaceDielectricData::GetCoefficient( int i, const mfem::ParGridFunction &U, const MaterialOperator &mat_op, const std::map &local_to_shared) const { @@ -124,21 +124,21 @@ SurfacePostOperator::InterfaceDielectricData::GetCoefficient( return {}; // For compiler warning } -SurfacePostOperator::SurfaceChargeData::SurfaceChargeData( +SurfaceChargeData::SurfaceChargeData( const config::CapacitanceData &data, mfem::ParMesh &mesh) { attr_markers.emplace_back(); mesh::AttrToMarker(mesh.bdr_attributes.Max(), data.attributes, attr_markers.back()); } -std::unique_ptr SurfacePostOperator::SurfaceChargeData::GetCoefficient( +std::unique_ptr SurfaceChargeData::GetCoefficient( int i, const mfem::ParGridFunction &U, const MaterialOperator &mat_op, const std::map &local_to_shared) const { return std::make_unique(U, mat_op, local_to_shared); } -SurfacePostOperator::SurfaceFluxData::SurfaceFluxData(const config::InductanceData &data, +SurfaceFluxData::SurfaceFluxData(const config::InductanceData &data, mfem::ParMesh &mesh) { // Check inputs. @@ -177,13 +177,36 @@ SurfacePostOperator::SurfaceFluxData::SurfaceFluxData(const config::InductanceDa mesh::AttrToMarker(mesh.bdr_attributes.Max(), data.attributes, attr_markers.back()); } -std::unique_ptr SurfacePostOperator::SurfaceFluxData::GetCoefficient( +std::unique_ptr SurfaceFluxData::GetCoefficient( int i, const mfem::ParGridFunction &U, const MaterialOperator &mat_op, const std::map &local_to_shared) const { return std::make_unique(U, direction, local_to_shared); } +namespace +{ +double GetSurfaceIntegral(const SurfaceData &data, + const mfem::ParGridFunction &U, + const mfem::ParGridFunction &ones, + const MaterialOperator &mat_op, + const std::map &local_to_shared) +{ + // Integrate the coefficient over the boundary attributes making up this surface index. + std::vector> fb; + mfem::ParLinearForm s(ones.ParFESpace()); + for (int i = 0; i < static_cast(data.attr_markers.size()); i++) + { + fb.emplace_back(data.GetCoefficient(i, U, mat_op, local_to_shared)); + s.AddBoundaryIntegrator(new BoundaryLFIntegrator(*fb.back()), data.attr_markers[i]); + } + s.UseFastAssembly(true); + s.Assemble(); + return s(ones); +} + +} // namespace + SurfacePostOperator::SurfacePostOperator(const IoData &iodata, const MaterialOperator &mat, const std::map &l2s, mfem::ParFiniteElementSpace &h1_fespace) @@ -219,7 +242,7 @@ SurfacePostOperator::GetInterfaceElectricFieldEnergy(int idx, auto it = eps_surfs.find(idx); MFEM_VERIFY(it != eps_surfs.end(), "Unknown dielectric loss postprocessing surface index requested!"); - return GetSurfaceIntegral(it->second, E); + return GetSurfaceIntegral(it->second, E, ones, mat_op, local_to_shared); } double SurfacePostOperator::GetInterfaceLossTangent(int idx) const @@ -236,7 +259,7 @@ double SurfacePostOperator::GetSurfaceElectricCharge(int idx, auto it = charge_surfs.find(idx); MFEM_VERIFY(it != charge_surfs.end(), "Unknown capacitance postprocessing surface index requested!"); - return GetSurfaceIntegral(it->second, E); + return GetSurfaceIntegral(it->second, E, ones, mat_op, local_to_shared); } double SurfacePostOperator::GetSurfaceMagneticFlux(int idx, @@ -245,23 +268,8 @@ double SurfacePostOperator::GetSurfaceMagneticFlux(int idx, auto it = flux_surfs.find(idx); MFEM_VERIFY(it != flux_surfs.end(), "Unknown inductance postprocessing surface index requested!"); - return GetSurfaceIntegral(it->second, B); + return GetSurfaceIntegral(it->second, B, ones, mat_op, local_to_shared); } -double SurfacePostOperator::GetSurfaceIntegral(const SurfaceData &data, - const mfem::ParGridFunction &U) const -{ - // Integrate the coefficient over the boundary attributes making up this surface index. - std::vector> fb; - mfem::ParLinearForm s(ones.ParFESpace()); - for (int i = 0; i < static_cast(data.attr_markers.size()); i++) - { - fb.emplace_back(data.GetCoefficient(i, U, mat_op, local_to_shared)); - s.AddBoundaryIntegrator(new BoundaryLFIntegrator(*fb.back()), data.attr_markers[i]); - } - s.UseFastAssembly(true); - s.Assemble(); - return s(ones); -} } // namespace palace diff --git a/palace/models/surfacepostoperator.hpp b/palace/models/surfacepostoperator.hpp index 697aca505..7c511c204 100644 --- a/palace/models/surfacepostoperator.hpp +++ b/palace/models/surfacepostoperator.hpp @@ -25,6 +25,50 @@ struct InductanceData; } // namespace config + +struct SurfaceData +{ + mutable std::vector> attr_markers; + + virtual ~SurfaceData() = default; + + virtual std::unique_ptr + GetCoefficient(int i, const mfem::ParGridFunction &U, const MaterialOperator &mat_op, + const std::map &local_to_shared) const = 0; +}; +struct InterfaceDielectricData : public SurfaceData +{ + DielectricInterfaceType type; + double epsilon, ts, tandelta; + std::vector sides; + + InterfaceDielectricData(const config::InterfaceDielectricData &data, + mfem::ParMesh &mesh); + + std::unique_ptr + GetCoefficient(int i, const mfem::ParGridFunction &U, const MaterialOperator &mat_op, + const std::map &local_to_shared) const override; +}; +struct SurfaceChargeData : public SurfaceData +{ + SurfaceChargeData(const config::CapacitanceData &data, mfem::ParMesh &mesh); + + std::unique_ptr + GetCoefficient(int i, const mfem::ParGridFunction &U, const MaterialOperator &mat_op, + const std::map &local_to_shared) const override; +}; +struct SurfaceFluxData : public SurfaceData +{ + mfem::Vector direction; + + SurfaceFluxData(const config::InductanceData &data, mfem::ParMesh &mesh); + + std::unique_ptr + GetCoefficient(int i, const mfem::ParGridFunction &U, const MaterialOperator &mat_op, + const std::map &local_to_shared) const override; +}; + + // // A class handling boundary surface postprocessing. // @@ -33,47 +77,6 @@ class SurfacePostOperator private: // Mapping from surface index to data structure containing surface postprocessing // information for surface loss, charge, or magnetic flux. - struct SurfaceData - { - mutable std::vector> attr_markers; - - virtual ~SurfaceData() = default; - - virtual std::unique_ptr - GetCoefficient(int i, const mfem::ParGridFunction &U, const MaterialOperator &mat_op, - const std::map &local_to_shared) const = 0; - }; - struct InterfaceDielectricData : public SurfaceData - { - DielectricInterfaceType type; - double epsilon, ts, tandelta; - std::vector sides; - - InterfaceDielectricData(const config::InterfaceDielectricData &data, - mfem::ParMesh &mesh); - - std::unique_ptr - GetCoefficient(int i, const mfem::ParGridFunction &U, const MaterialOperator &mat_op, - const std::map &local_to_shared) const override; - }; - struct SurfaceChargeData : public SurfaceData - { - SurfaceChargeData(const config::CapacitanceData &data, mfem::ParMesh &mesh); - - std::unique_ptr - GetCoefficient(int i, const mfem::ParGridFunction &U, const MaterialOperator &mat_op, - const std::map &local_to_shared) const override; - }; - struct SurfaceFluxData : public SurfaceData - { - mfem::Vector direction; - - SurfaceFluxData(const config::InductanceData &data, mfem::ParMesh &mesh); - - std::unique_ptr - GetCoefficient(int i, const mfem::ParGridFunction &U, const MaterialOperator &mat_op, - const std::map &local_to_shared) const override; - }; std::map eps_surfs; std::map charge_surfs; std::map flux_surfs; @@ -86,9 +89,6 @@ class SurfacePostOperator // Unit function used for computing surface integrals. mfem::ParGridFunction ones; - - double GetSurfaceIntegral(const SurfaceData &data, const mfem::ParGridFunction &U) const; - public: SurfacePostOperator(const IoData &iodata, const MaterialOperator &mat, const std::map &l2s, diff --git a/palace/models/waveportoperator.cpp b/palace/models/waveportoperator.cpp index ed5dc6963..ed26affca 100644 --- a/palace/models/waveportoperator.cpp +++ b/palace/models/waveportoperator.cpp @@ -539,6 +539,7 @@ void WavePortData::GetTrueDofs(const mfem::Array &dbc_marker, } } +// static void WavePortData::GetInitialSpace(int nt, int nn, petsc::PetscParVector &y0) { // Initial space chosen as such that B v₀ = y₀, with y₀ = [y₀ₜ, 0, ... 0]ᵀ ⟂ null(A) @@ -557,37 +558,35 @@ void WavePortData::GetInitialSpace(int nt, int nn, petsc::PetscParVector &y0) y0.RestoreArray(py0); } -std::complex WavePortData::Solve(petsc::PetscParVector &y0, - petsc::PetscParVector &e0, - petsc::PetscParVector &e, - petsc::PetscScatter &scatter) + +std::complex WavePortData::Solve() { double eig[2]; if (A) // Only on root { // The y0 and e0 vectors are still parallel vectors, but with all data on root. We want // true sequential vectors. - PetscScalar *pe0 = e0.GetArray(); - petsc::PetscParVector e0s(e0.GetSize(), pe0); + PetscScalar *pe0 = e0->GetArray(); + petsc::PetscParVector e0s(e0->GetSize(), pe0); // Set starting vector. { - PetscScalar *py0 = y0.GetArray(); - petsc::PetscParVector y0s(y0.GetSize(), py0); + PetscScalar *py0 = y0->GetArray(); + petsc::PetscParVector y0s(y0->GetSize(), py0); eigen->SetInitialSpace(y0s); - y0.RestoreArray(py0); + y0->RestoreArray(py0); } #if 0 // Alternatively, use B-orthogonal initial space. Probably want to call SetBMat for // the eigensolver in this case. { - PetscScalar *py0 = y0.GetArray(); - petsc::PetscParVector y0s(y0.GetSize(), py0); + PetscScalar *py0 = y0->GetArray(); + petsc::PetscParVector y0s(y0->GetSize(), py0); petsc::PetscParVector v0s(y0s); ksp->Mult(y0s, v0s); eigen->SetInitialSpace(v0s); - y0.RestoreArray(py0); + y0->RestoreArray(py0); } #endif @@ -598,12 +597,12 @@ std::complex WavePortData::Solve(petsc::PetscParVector &y0, MFEM_VERIFY(num_conv >= mode_idx, "Wave port eigensolver did not converge!"); eigen->GetEigenvalue(mode_idx - 1, eig[0], eig[1]); eigen->GetEigenvector(mode_idx - 1, e0s); - e0.RestoreArray(pe0); + e0->RestoreArray(pe0); } // Scatter the result to all processors. - scatter.Reverse(e0, e); - Mpi::Broadcast(2, eig, 0, e.GetComm()); + scatter->Reverse(*e0, *e); + Mpi::Broadcast(2, eig, 0, e->GetComm()); return {eig[0], eig[1]}; } @@ -631,7 +630,7 @@ void WavePortData::Initialize(double omega) } // Configure and solve the eigenvalue problem for the desired boundary mode. - std::complex lambda = Solve(*y0, *e0, *e, *scatter); + std::complex lambda = Solve(); // Extract the eigenmode solution and postprocess. The extracted eigenvalue is λ = // Θ²/(Θ²-kₙ²). @@ -747,8 +746,8 @@ std::complex WavePortData::GetSParameter(mfem::ParComplexGridFunction &E return {-(*sr)(E.real()) - (*si)(E.imag()), -(*sr)(E.imag()) + (*si)(E.real())}; } -std::complex WavePortData::GetPower(mfem::ParComplexGridFunction &E, - mfem::ParComplexGridFunction &B, +std::complex WavePortData::GetPower(mfem::ParComplexGridFunction &Ein, + mfem::ParComplexGridFunction &Bin, const MaterialOperator &mat_op, const std::map &local_to_shared) const { @@ -756,9 +755,9 @@ std::complex WavePortData::GetPower(mfem::ParComplexGridFunction &E, // using the computed E and H = μ⁻¹ B fields. The linear form is reconstructed from // scratch each time due to changing H. The BdrCurrentVectorCoefficient computes -n x H, // where n is an outward normal. - auto &nd_fespace = *E.ParFESpace(); - BdrCurrentVectorCoefficient nxHr_func(B.real(), mat_op, local_to_shared); - BdrCurrentVectorCoefficient nxHi_func(B.imag(), mat_op, local_to_shared); + auto &nd_fespace = *Ein.ParFESpace(); + BdrCurrentVectorCoefficient nxHr_func(Bin.real(), mat_op, local_to_shared); + BdrCurrentVectorCoefficient nxHi_func(Bin.imag(), mat_op, local_to_shared); mfem::ParLinearForm pr(&nd_fespace), pi(&nd_fespace); pr.AddBoundaryIntegrator(new VectorFEBoundaryLFIntegrator(nxHr_func), attr_marker); pi.AddBoundaryIntegrator(new VectorFEBoundaryLFIntegrator(nxHi_func), attr_marker); @@ -766,7 +765,7 @@ std::complex WavePortData::GetPower(mfem::ParComplexGridFunction &E, pi.UseFastAssembly(true); pr.Assemble(); pi.Assemble(); - return {pr(E.real()) + pi(E.imag()), pr(E.imag()) - pi(E.real())}; + return {pr(Ein.real()) + pi(Ein.imag()), pr(Ein.imag()) - pi(Ein.real())}; } WavePortOperator::WavePortOperator(const IoData &iod, const MaterialOperator &mat, @@ -777,13 +776,11 @@ WavePortOperator::WavePortOperator(const IoData &iod, const MaterialOperator &ma // Set up wave port boundary conditions. MFEM_VERIFY(nd_fespace.GetParMesh() == h1_fespace.GetParMesh(), "Mesh mismatch in WavePortOperator FE spaces!"); - SetUpBoundaryProperties(iodata, mat_op, nd_fespace, h1_fespace); - PrintBoundaryInfo(iodata, *nd_fespace.GetParMesh()); + SetUpBoundaryProperties(nd_fespace, h1_fespace); + PrintBoundaryInfo(*nd_fespace.GetParMesh()); } -void WavePortOperator::SetUpBoundaryProperties(const IoData &iodata, - const MaterialOperator &mat_op, - mfem::ParFiniteElementSpace &nd_fespace, +void WavePortOperator::SetUpBoundaryProperties(mfem::ParFiniteElementSpace &nd_fespace, mfem::ParFiniteElementSpace &h1_fespace) { // Check that wave port boundary attributes have been specified correctly. @@ -861,7 +858,7 @@ void WavePortOperator::SetUpBoundaryProperties(const IoData &iodata, } } -void WavePortOperator::PrintBoundaryInfo(const IoData &iodata, mfem::ParMesh &mesh) +void WavePortOperator::PrintBoundaryInfo(mfem::ParMesh &mesh) { // Print out BC info for all port attributes. if (ports.empty()) diff --git a/palace/models/waveportoperator.hpp b/palace/models/waveportoperator.hpp index 1769ffb74..a0ee90a66 100644 --- a/palace/models/waveportoperator.hpp +++ b/palace/models/waveportoperator.hpp @@ -72,9 +72,9 @@ class WavePortData mfem::Array &h1_tdof_list); // Configure and solve the linear eigenvalue problem for the boundary mode. - void GetInitialSpace(int nt, int nn, petsc::PetscParVector &y0); - std::complex Solve(petsc::PetscParVector &y0, petsc::PetscParVector &e0, - petsc::PetscParVector &e, petsc::PetscScatter &scatter); + static void GetInitialSpace(int nt, int nn, petsc::PetscParVector &y0); + + std::complex Solve(); public: WavePortData(const config::WavePortData &data, const MaterialOperator &mat_op, @@ -146,10 +146,9 @@ class WavePortOperator // Mapping from port index to data structure containing port information. std::map ports; mfem::Array port_marker; - void SetUpBoundaryProperties(const IoData &iodata, const MaterialOperator &mat_op, - mfem::ParFiniteElementSpace &nd_fespace, + void SetUpBoundaryProperties(mfem::ParFiniteElementSpace &nd_fespace, mfem::ParFiniteElementSpace &h1_fespace); - void PrintBoundaryInfo(const IoData &iodata, mfem::ParMesh &mesh); + void PrintBoundaryInfo(mfem::ParMesh &mesh); // Compute boundary modes for all wave port boundaries at the specified frequency. void Initialize(double omega); From d70060b53fc5ce2945065e6c85c09e0042503cc9 Mon Sep 17 00:00:00 2001 From: Hugh Carson Date: Tue, 25 Jul 2023 14:13:55 -0400 Subject: [PATCH 2/2] style --- palace/drivers/eigensolver.cpp | 2 +- palace/linalg/feast.cpp | 2 +- palace/models/romoperator.hpp | 6 ++--- palace/models/surfacepostoperator.cpp | 36 +++++++++++++-------------- palace/models/surfacepostoperator.hpp | 14 +++++------ palace/models/waveportoperator.cpp | 1 - 6 files changed, 28 insertions(+), 33 deletions(-) diff --git a/palace/drivers/eigensolver.cpp b/palace/drivers/eigensolver.cpp index bb16100ae..f9b9b9283 100644 --- a/palace/drivers/eigensolver.cpp +++ b/palace/drivers/eigensolver.cpp @@ -164,7 +164,7 @@ void EigenSolver::Solve(std::vector> &mesh, std::unique_ptr ksp; std::unique_ptr pc; #if defined(PALACE_WITH_SLEPC) - auto * const feast = dynamic_cast(eigen.get()); + auto *const feast = dynamic_cast(eigen.get()); if (feast) { // Configure the FEAST integration contour. The linear solvers are set up inside the diff --git a/palace/linalg/feast.cpp b/palace/linalg/feast.cpp index 1008bfbde..09b40ea44 100644 --- a/palace/linalg/feast.cpp +++ b/palace/linalg/feast.cpp @@ -1198,7 +1198,7 @@ void FeastPEPSolver::SolveProjectedProblem(const petsc::PetscDenseMatrix &Q_, for (PetscInt i = 0; i < m0; i++) { eig_[i] = alpha[sort[i]]; - const PetscScalar* const local_pXQ = XQ->GetArrayRead(); + const PetscScalar *const local_pXQ = XQ->GetArrayRead(); PetscScalar *pXQ0 = XQ0->GetArray(); PalacePetscCall(PetscArraycpy(pXQ0 + mQ * i, local_pXQ + 2 * mQ * sort[i], mQ)); XQ->RestoreArrayRead(local_pXQ); diff --git a/palace/models/romoperator.hpp b/palace/models/romoperator.hpp index c72f986fa..cbda23b65 100644 --- a/palace/models/romoperator.hpp +++ b/palace/models/romoperator.hpp @@ -66,10 +66,10 @@ class RomOperator // Helper functions for reduced-order matrix or vector construction/update. static void BVMatProjectInternal(petsc::PetscDenseMatrix &V, petsc::PetscParMatrix &A, - petsc::PetscDenseMatrix &Ar, petsc::PetscParVector &r, int n0, - int n); + petsc::PetscDenseMatrix &Ar, petsc::PetscParVector &r, + int n0, int n); static void BVDotVecInternal(petsc::PetscDenseMatrix &V, petsc::PetscParVector &b, - petsc::PetscParVector &br, int n0, int n); + petsc::PetscParVector &br, int n0, int n); public: RomOperator(const IoData &iodata, SpaceOperator &sp, int nmax); diff --git a/palace/models/surfacepostoperator.cpp b/palace/models/surfacepostoperator.cpp index c2b7b5a83..ea6353984 100644 --- a/palace/models/surfacepostoperator.cpp +++ b/palace/models/surfacepostoperator.cpp @@ -101,9 +101,9 @@ InterfaceDielectricData::InterfaceDielectricData( } std::unique_ptr -InterfaceDielectricData::GetCoefficient( - int i, const mfem::ParGridFunction &U, const MaterialOperator &mat_op, - const std::map &local_to_shared) const +InterfaceDielectricData::GetCoefficient(int i, const mfem::ParGridFunction &U, + const MaterialOperator &mat_op, + const std::map &local_to_shared) const { switch (type) { @@ -124,22 +124,22 @@ InterfaceDielectricData::GetCoefficient( return {}; // For compiler warning } -SurfaceChargeData::SurfaceChargeData( - const config::CapacitanceData &data, mfem::ParMesh &mesh) +SurfaceChargeData::SurfaceChargeData(const config::CapacitanceData &data, + mfem::ParMesh &mesh) { attr_markers.emplace_back(); mesh::AttrToMarker(mesh.bdr_attributes.Max(), data.attributes, attr_markers.back()); } -std::unique_ptr SurfaceChargeData::GetCoefficient( - int i, const mfem::ParGridFunction &U, const MaterialOperator &mat_op, - const std::map &local_to_shared) const +std::unique_ptr +SurfaceChargeData::GetCoefficient(int i, const mfem::ParGridFunction &U, + const MaterialOperator &mat_op, + const std::map &local_to_shared) const { return std::make_unique(U, mat_op, local_to_shared); } -SurfaceFluxData::SurfaceFluxData(const config::InductanceData &data, - mfem::ParMesh &mesh) +SurfaceFluxData::SurfaceFluxData(const config::InductanceData &data, mfem::ParMesh &mesh) { // Check inputs. MFEM_VERIFY(data.direction.length() == 2 && @@ -177,19 +177,18 @@ SurfaceFluxData::SurfaceFluxData(const config::InductanceData &data, mesh::AttrToMarker(mesh.bdr_attributes.Max(), data.attributes, attr_markers.back()); } -std::unique_ptr SurfaceFluxData::GetCoefficient( - int i, const mfem::ParGridFunction &U, const MaterialOperator &mat_op, - const std::map &local_to_shared) const +std::unique_ptr +SurfaceFluxData::GetCoefficient(int i, const mfem::ParGridFunction &U, + const MaterialOperator &mat_op, + const std::map &local_to_shared) const { return std::make_unique(U, direction, local_to_shared); } namespace { -double GetSurfaceIntegral(const SurfaceData &data, - const mfem::ParGridFunction &U, - const mfem::ParGridFunction &ones, - const MaterialOperator &mat_op, +double GetSurfaceIntegral(const SurfaceData &data, const mfem::ParGridFunction &U, + const mfem::ParGridFunction &ones, const MaterialOperator &mat_op, const std::map &local_to_shared) { // Integrate the coefficient over the boundary attributes making up this surface index. @@ -205,7 +204,7 @@ double GetSurfaceIntegral(const SurfaceData &data, return s(ones); } -} // namespace +} // namespace SurfacePostOperator::SurfacePostOperator(const IoData &iodata, const MaterialOperator &mat, const std::map &l2s, @@ -271,5 +270,4 @@ double SurfacePostOperator::GetSurfaceMagneticFlux(int idx, return GetSurfaceIntegral(it->second, B, ones, mat_op, local_to_shared); } - } // namespace palace diff --git a/palace/models/surfacepostoperator.hpp b/palace/models/surfacepostoperator.hpp index 7c511c204..691f79a47 100644 --- a/palace/models/surfacepostoperator.hpp +++ b/palace/models/surfacepostoperator.hpp @@ -25,7 +25,6 @@ struct InductanceData; } // namespace config - struct SurfaceData { mutable std::vector> attr_markers; @@ -34,7 +33,7 @@ struct SurfaceData virtual std::unique_ptr GetCoefficient(int i, const mfem::ParGridFunction &U, const MaterialOperator &mat_op, - const std::map &local_to_shared) const = 0; + const std::map &local_to_shared) const = 0; }; struct InterfaceDielectricData : public SurfaceData { @@ -42,12 +41,11 @@ struct InterfaceDielectricData : public SurfaceData double epsilon, ts, tandelta; std::vector sides; - InterfaceDielectricData(const config::InterfaceDielectricData &data, - mfem::ParMesh &mesh); + InterfaceDielectricData(const config::InterfaceDielectricData &data, mfem::ParMesh &mesh); std::unique_ptr GetCoefficient(int i, const mfem::ParGridFunction &U, const MaterialOperator &mat_op, - const std::map &local_to_shared) const override; + const std::map &local_to_shared) const override; }; struct SurfaceChargeData : public SurfaceData { @@ -55,7 +53,7 @@ struct SurfaceChargeData : public SurfaceData std::unique_ptr GetCoefficient(int i, const mfem::ParGridFunction &U, const MaterialOperator &mat_op, - const std::map &local_to_shared) const override; + const std::map &local_to_shared) const override; }; struct SurfaceFluxData : public SurfaceData { @@ -65,10 +63,9 @@ struct SurfaceFluxData : public SurfaceData std::unique_ptr GetCoefficient(int i, const mfem::ParGridFunction &U, const MaterialOperator &mat_op, - const std::map &local_to_shared) const override; + const std::map &local_to_shared) const override; }; - // // A class handling boundary surface postprocessing. // @@ -89,6 +86,7 @@ class SurfacePostOperator // Unit function used for computing surface integrals. mfem::ParGridFunction ones; + public: SurfacePostOperator(const IoData &iodata, const MaterialOperator &mat, const std::map &l2s, diff --git a/palace/models/waveportoperator.cpp b/palace/models/waveportoperator.cpp index ed26affca..3c32c4064 100644 --- a/palace/models/waveportoperator.cpp +++ b/palace/models/waveportoperator.cpp @@ -558,7 +558,6 @@ void WavePortData::GetInitialSpace(int nt, int nn, petsc::PetscParVector &y0) y0.RestoreArray(py0); } - std::complex WavePortData::Solve() { double eig[2];