From 148d51e623d6b78bad372e7b997568348d549d45 Mon Sep 17 00:00:00 2001 From: Sebastian Grimberg Date: Mon, 17 Jul 2023 10:51:00 -0700 Subject: [PATCH] Fix full parallel assembly for AMS solver and wave ports --- palace/linalg/ams.cpp | 36 ++++----- palace/models/waveportoperator.cpp | 118 ++++++++++++++++------------- 2 files changed, 84 insertions(+), 70 deletions(-) diff --git a/palace/linalg/ams.cpp b/palace/linalg/ams.cpp index 72f72d209..875fcafb9 100644 --- a/palace/linalg/ams.cpp +++ b/palace/linalg/ams.cpp @@ -51,12 +51,13 @@ void HypreAmsSolver::ConstructAuxiliaryMatrices(mfem::ParFiniteElementSpace &nd_ // HypreAMS:Init. Start with the discrete gradient matrix. { // XX TODO: Partial assembly option? - auto grad = std::make_unique(&h1_fespace, &nd_fespace); - grad->AddDomainInterpolator(new mfem::GradientInterpolator); - grad->SetAssemblyLevel(mfem::AssemblyLevel::LEGACY); - grad->Assemble(); - grad->Finalize(); - ParOperator RAP_G(std::move(grad), h1_fespace, nd_fespace, true); + mfem::DiscreteLinearOperator grad(&h1_fespace, &nd_fespace); + grad.AddDomainInterpolator(new mfem::GradientInterpolator); + grad.SetAssemblyLevel(mfem::AssemblyLevel::LEGACY); + grad.Assemble(); + grad.Finalize(); + ParOperator RAP_G(std::unique_ptr(grad.LoseMat()), h1_fespace, + nd_fespace, true); G = RAP_G.StealParallelAssemble(); } @@ -112,18 +113,17 @@ void HypreAmsSolver::ConstructAuxiliaryMatrices(mfem::ParFiniteElementSpace &nd_ } else { - { - // XX TODO: Partial assembly option? - mfem::ParFiniteElementSpace h1d_fespace(&mesh, h1_fespace.FEColl(), space_dim, - mfem::Ordering::byVDIM); - auto pi = std::make_unique(&h1d_fespace, &nd_fespace); - pi->AddDomainInterpolator(new mfem::IdentityInterpolator); - pi->SetAssemblyLevel(mfem::AssemblyLevel::LEGACY); - pi->Assemble(); - pi->Finalize(); - ParOperator RAP_Pi(std::move(pi), h1d_fespace, nd_fespace, true); - Pi = RAP_Pi.StealParallelAssemble(); - } + // XX TODO: Partial assembly option? + mfem::ParFiniteElementSpace h1d_fespace(&mesh, h1_fespace.FEColl(), space_dim, + mfem::Ordering::byVDIM); + mfem::DiscreteLinearOperator pi(&h1d_fespace, &nd_fespace); + pi.AddDomainInterpolator(new mfem::IdentityInterpolator); + pi.SetAssemblyLevel(mfem::AssemblyLevel::LEGACY); + pi.Assemble(); + pi.Finalize(); + ParOperator RAP_Pi(std::unique_ptr(pi.LoseMat()), h1d_fespace, + nd_fespace, true); + Pi = RAP_Pi.StealParallelAssemble(); if (cycle_type >= 10) { // Get blocks of Pi corresponding to each component, and free Pi. diff --git a/palace/models/waveportoperator.cpp b/palace/models/waveportoperator.cpp index ff16a54cc..c537704e7 100644 --- a/palace/models/waveportoperator.cpp +++ b/palace/models/waveportoperator.cpp @@ -84,12 +84,13 @@ std::unique_ptr GetBtt(const MaterialOperator &mat_op, constexpr MaterialPropertyType MatType = MaterialPropertyType::INV_PERMEABILITY; constexpr MeshElementType ElemType = MeshElementType::BDR_SUBMESH; MaterialPropertyCoefficient muinv_func(mat_op); - auto btt = std::make_unique(&nd_fespace); - btt->AddDomainIntegrator(new mfem::VectorFEMassIntegrator(muinv_func)); - btt->SetAssemblyLevel(mfem::AssemblyLevel::LEGACY); - btt->Assemble(skip_zeros); - btt->Finalize(skip_zeros); - return std::make_unique(std::move(btt), nd_fespace); + mfem::SymmetricBilinearForm btt(&nd_fespace); + btt.AddDomainIntegrator(new mfem::VectorFEMassIntegrator(muinv_func)); + btt.SetAssemblyLevel(mfem::AssemblyLevel::LEGACY); + btt.Assemble(skip_zeros); + btt.Finalize(skip_zeros); + return std::make_unique(std::unique_ptr(btt.LoseMat()), + nd_fespace); } std::unique_ptr GetBtn(const MaterialOperator &mat_op, @@ -100,12 +101,13 @@ std::unique_ptr GetBtn(const MaterialOperator &mat_op, constexpr MaterialPropertyType MatType = MaterialPropertyType::INV_PERMEABILITY; constexpr MeshElementType ElemType = MeshElementType::BDR_SUBMESH; MaterialPropertyCoefficient muinv_func(mat_op); - auto btn = std::make_unique(&h1_fespace, &nd_fespace); - btn->AddDomainIntegrator(new mfem::MixedVectorGradientIntegrator(muinv_func)); - btn->SetAssemblyLevel(mfem::AssemblyLevel::LEGACY); - btn->Assemble(skip_zeros); - btn->Finalize(skip_zeros); - return std::make_unique(std::move(btn), h1_fespace, nd_fespace, false); + mfem::MixedBilinearForm btn(&h1_fespace, &nd_fespace); + btn.AddDomainIntegrator(new mfem::MixedVectorGradientIntegrator(muinv_func)); + btn.SetAssemblyLevel(mfem::AssemblyLevel::LEGACY); + btn.Assemble(skip_zeros); + btn.Finalize(skip_zeros); + return std::make_unique(std::unique_ptr(btn.LoseMat()), + h1_fespace, nd_fespace, false); } std::array, 3> GetBnn(const MaterialOperator &mat_op, @@ -115,38 +117,44 @@ std::array, 3> GetBnn(const MaterialOperator &mat_o constexpr MaterialPropertyType MatTypeMuInv = MaterialPropertyType::INV_PERMEABILITY; constexpr MeshElementType ElemType = MeshElementType::BDR_SUBMESH; MaterialPropertyCoefficient muinv_func(mat_op); - auto bnn1 = std::make_unique(&h1_fespace); - bnn1->AddDomainIntegrator(new mfem::DiffusionIntegrator(muinv_func)); - bnn1->SetAssemblyLevel(mfem::AssemblyLevel::LEGACY); - bnn1->Assemble(skip_zeros); - bnn1->Finalize(skip_zeros); + mfem::SymmetricBilinearForm bnn1(&h1_fespace); + bnn1.AddDomainIntegrator(new mfem::DiffusionIntegrator(muinv_func)); + bnn1.SetAssemblyLevel(mfem::AssemblyLevel::LEGACY); + bnn1.Assemble(skip_zeros); + bnn1.Finalize(skip_zeros); constexpr MaterialPropertyType MatTypeEpsReal = MaterialPropertyType::PERMITTIVITY_REAL; NormalProjectedCoefficient epsilon_func( std::make_unique>(mat_op)); - auto bnn2r = std::make_unique(&h1_fespace); - bnn2r->AddDomainIntegrator(new mfem::MixedScalarMassIntegrator(epsilon_func)); - bnn2r->SetAssemblyLevel(mfem::AssemblyLevel::LEGACY); - bnn2r->Assemble(skip_zeros); - bnn2r->Finalize(skip_zeros); + mfem::SymmetricBilinearForm bnn2r(&h1_fespace); + bnn2r.AddDomainIntegrator(new mfem::MixedScalarMassIntegrator(epsilon_func)); + bnn2r.SetAssemblyLevel(mfem::AssemblyLevel::LEGACY); + bnn2r.Assemble(skip_zeros); + bnn2r.Finalize(skip_zeros); // Contribution for loss tangent: ε -> ε * (1 - i tan(δ)). if (!mat_op.HasLossTangent()) { - return {std::make_unique(std::move(bnn1), h1_fespace), - std::make_unique(std::move(bnn2r), h1_fespace), nullptr}; + return {std::make_unique( + std::unique_ptr(bnn1.LoseMat()), h1_fespace), + std::make_unique( + std::unique_ptr(bnn2r.LoseMat()), h1_fespace), + nullptr}; } constexpr MaterialPropertyType MatTypeEpsImag = MaterialPropertyType::PERMITTIVITY_IMAG; NormalProjectedCoefficient negepstandelta_func( std::make_unique>(mat_op)); - auto bnn2i = std::make_unique(&h1_fespace); - bnn2i->AddDomainIntegrator(new mfem::MixedScalarMassIntegrator(negepstandelta_func)); - bnn2i->SetAssemblyLevel(mfem::AssemblyLevel::LEGACY); - bnn2i->Assemble(skip_zeros); - bnn2i->Finalize(skip_zeros); - return {std::make_unique(std::move(bnn1), h1_fespace), - std::make_unique(std::move(bnn2r), h1_fespace), - std::make_unique(std::move(bnn2i), h1_fespace)}; + mfem::SymmetricBilinearForm bnn2i(&h1_fespace); + bnn2i.AddDomainIntegrator(new mfem::MixedScalarMassIntegrator(negepstandelta_func)); + bnn2i.SetAssemblyLevel(mfem::AssemblyLevel::LEGACY); + bnn2i.Assemble(skip_zeros); + bnn2i.Finalize(skip_zeros); + return {std::make_unique(std::unique_ptr(bnn1.LoseMat()), + h1_fespace), + std::make_unique( + std::unique_ptr(bnn2r.LoseMat()), h1_fespace), + std::make_unique( + std::unique_ptr(bnn2i.LoseMat()), h1_fespace)}; } std::array, 3> GetAtt(const MaterialOperator &mat_op, @@ -157,36 +165,42 @@ std::array, 3> GetAtt(const MaterialOperator &mat_o constexpr MeshElementType ElemType = MeshElementType::BDR_SUBMESH; NormalProjectedCoefficient muinv_func( std::make_unique>(mat_op)); - auto att1 = std::make_unique(&nd_fespace); - att1->AddDomainIntegrator(new mfem::CurlCurlIntegrator(muinv_func)); - att1->SetAssemblyLevel(mfem::AssemblyLevel::LEGACY); - att1->Assemble(skip_zeros); - att1->Finalize(skip_zeros); + mfem::SymmetricBilinearForm att1(&nd_fespace); + att1.AddDomainIntegrator(new mfem::CurlCurlIntegrator(muinv_func)); + att1.SetAssemblyLevel(mfem::AssemblyLevel::LEGACY); + att1.Assemble(skip_zeros); + att1.Finalize(skip_zeros); constexpr MaterialPropertyType MatTypeEpsReal = MaterialPropertyType::PERMITTIVITY_REAL; MaterialPropertyCoefficient epsilon_func(mat_op); - auto att2r = std::make_unique(&nd_fespace); - att2r->AddDomainIntegrator(new mfem::VectorFEMassIntegrator(epsilon_func)); - att2r->SetAssemblyLevel(mfem::AssemblyLevel::LEGACY); - att2r->Assemble(skip_zeros); - att2r->Finalize(skip_zeros); + mfem::SymmetricBilinearForm att2r(&nd_fespace); + att2r.AddDomainIntegrator(new mfem::VectorFEMassIntegrator(epsilon_func)); + att2r.SetAssemblyLevel(mfem::AssemblyLevel::LEGACY); + att2r.Assemble(skip_zeros); + att2r.Finalize(skip_zeros); // Contribution for loss tangent: ε -> ε * (1 - i tan(δ)). if (!mat_op.HasLossTangent()) { - return {std::make_unique(std::move(att1), nd_fespace), - std::make_unique(std::move(att2r), nd_fespace), nullptr}; + return {std::make_unique( + std::unique_ptr(att1.LoseMat()), nd_fespace), + std::make_unique( + std::unique_ptr(att2r.LoseMat()), nd_fespace), + nullptr}; } constexpr MaterialPropertyType MatTypeEpsImag = MaterialPropertyType::PERMITTIVITY_IMAG; MaterialPropertyCoefficient negepstandelta_func(mat_op); - auto att2i = std::make_unique(&nd_fespace); - att2i->AddDomainIntegrator(new mfem::VectorFEMassIntegrator(negepstandelta_func)); - att2i->SetAssemblyLevel(mfem::AssemblyLevel::LEGACY); - att2i->Assemble(skip_zeros); - att2i->Finalize(skip_zeros); - return {std::make_unique(std::move(att1), nd_fespace), - std::make_unique(std::move(att2r), nd_fespace), - std::make_unique(std::move(att2i), nd_fespace)}; + mfem::SymmetricBilinearForm att2i(&nd_fespace); + att2i.AddDomainIntegrator(new mfem::VectorFEMassIntegrator(negepstandelta_func)); + att2i.SetAssemblyLevel(mfem::AssemblyLevel::LEGACY); + att2i.Assemble(skip_zeros); + att2i.Finalize(skip_zeros); + return {std::make_unique(std::unique_ptr(att1.LoseMat()), + nd_fespace), + std::make_unique( + std::unique_ptr(att2r.LoseMat()), nd_fespace), + std::make_unique( + std::unique_ptr(att2i.LoseMat()), nd_fespace)}; } std::array, 6>