From 0066a64dc4731d1b579838671a401bdf220c8b0a Mon Sep 17 00:00:00 2001 From: Joachim Schoeberl Date: Tue, 10 Jan 2023 18:30:53 +0100 Subject: [PATCH] parallel diagonal assembling --- comp/bilinearform.cpp | 74 +++++++++++++++++++++++-------------------- comp/bilinearform.hpp | 18 +++++------ 2 files changed, 49 insertions(+), 43 deletions(-) diff --git a/comp/bilinearform.cpp b/comp/bilinearform.cpp index bfbad27f8..2f304bdaa 100644 --- a/comp/bilinearform.cpp +++ b/comp/bilinearform.cpp @@ -1641,49 +1641,55 @@ namespace ngcomp if (vb==VOL && diagonal) { - double prevtime = WallTime(); - Array dnums; + static Timer tdiag("AssembleDiag"); RegionTimer reg(tdiag); + + // double prevtime = WallTime(); + // Array dnums; - for (int i = 0; i < ne; i++) + // for (int i = 0; i < ne; i++) + IterateElements + (*fespace, vb, clh, [&] (FESpace::Element el, LocalHeap & lh) { - HeapReset hr(clh); - ElementId id(vb, i); - gcnt++; - + // HeapReset hr(clh); + // ElementId id(vb, i); + // gcnt++; + /* if (clock()-prevtime > 0.1 * CLOCKS_PER_SEC) { cout << IM(3) << "\rassemble element " << i << "/" << ne << flush; ma->SetThreadPercentage ( 100.0*gcnt / (loopsteps) ); prevtime = clock(); } - + */ // clh.CleanUp(heapp); - if (!fespace->DefinedOn (vb,ma->GetElIndex (id))) continue; + // if (!fespace->DefinedOn (vb,ma->GetElIndex (id))) continue; - const FiniteElement & fel = fespace->GetFE (id, clh); - ElementTransformation & eltrans = ma->GetTrafo (id, clh); - auto ei = ElementId(vb,i); - fespace->GetDofNrs (ei, dnums); + const FiniteElement & fel = fespace->GetFE (el, lh); + ElementTransformation & eltrans = ma->GetTrafo (el, lh); + // auto ei = ElementId(vb,i); + // fespace->GetDofNrs (ei, dnums); + FlatArray dnums = el.GetDofs(); - FlatVector sum_diag(dnums.Size()*fespace->GetDimension(), clh); + FlatVector sum_diag(dnums.Size()*fespace->GetDimension(), lh); sum_diag = 0; for (auto & bfip : VB_parts[vb]) { const BilinearFormIntegrator & bfi = *bfip; - if (!bfi.DefinedOn (ma->GetElIndex (id))) continue; - if (!bfi.DefinedOnElement (i)) continue; + if (!bfi.DefinedOn (el.GetIndex())) continue; + if (!bfi.DefinedOnElement (el.Nr())) continue; - FlatVector diag(sum_diag.Size(), clh); + FlatVector diag(sum_diag.Size(), lh); try { - bfi.CalcElementMatrixDiag (fel, eltrans, diag, clh); + bfi.CalcElementMatrixDiag (fel, eltrans, diag, lh); if (printelmat) { + lock_guard guard(printelmat_mutex); testout->precision(8); - (*testout) << "elnum= " << ElementId(vb,i) << endl; + (*testout) << "elnum= " << el << endl; (*testout) << "integrator " << bfi.Name() << endl; (*testout) << "dnums = " << endl << dnums << endl; (*testout) << "diag-elmat = " << endl << diag << endl; @@ -1705,14 +1711,14 @@ namespace ngcomp sum_diag += diag; } - AddDiagElementMatrix (dnums, sum_diag, 1, i, clh); + AddDiagElementMatrix (dnums, sum_diag, 1, el.Nr(), lh); if (check_unused) for (int j = 0; j < dnums.Size(); j++) if (IsRegularDof(dnums[j])) useddof[dnums[j]] = true; - } - cout << IM(3) << "\rassemble element " << ne << "/" << ne << endl; + }); + // cout << IM(3) << "\rassemble element " << ne << "/" << ne << endl; } else // not diagonal { @@ -6257,8 +6263,8 @@ namespace ngcomp template void S_BilinearForm :: - AddDiagElementMatrix (const Array & dnums1, - const FlatVector & diag, + AddDiagElementMatrix (FlatArray dnums1, + FlatVector diag, bool inner_element, int elnr, LocalHeap & lh) { @@ -6877,8 +6883,8 @@ namespace ngcomp template void T_BilinearFormDiagonal :: - AddDiagElementMatrix (const Array & dnums1, - const FlatVector & diag, + AddDiagElementMatrix (FlatArray dnums1, + FlatVector diag, bool inner_element, int elnr, LocalHeap & lh) { @@ -6904,10 +6910,10 @@ namespace ngcomp /// template <> void T_BilinearFormDiagonal:: - AddDiagElementMatrix (const Array & dnums1, - const FlatVector & diag, - bool inner_element, int elnr, - LocalHeap & lh) + AddDiagElementMatrix (FlatArray dnums1, + FlatVector diag, + bool inner_element, int elnr, + LocalHeap & lh) { for (int i = 0; i < dnums1.Size(); i++) if (IsRegularDof(dnums1[i])) @@ -6916,10 +6922,10 @@ namespace ngcomp /// template <> void T_BilinearFormDiagonal:: - AddDiagElementMatrix (const Array & dnums1, - const FlatVector & diag, - bool inner_element, int elnr, - LocalHeap & lh) + AddDiagElementMatrix (FlatArray dnums1, + FlatVector diag, + bool inner_element, int elnr, + LocalHeap & lh) { for (int i = 0; i < dnums1.Size(); i++) if (IsRegularDof(dnums1[i])) diff --git a/comp/bilinearform.hpp b/comp/bilinearform.hpp index d53dd0405..cef00dd56 100644 --- a/comp/bilinearform.hpp +++ b/comp/bilinearform.hpp @@ -578,8 +578,8 @@ namespace ngcomp */ // { cerr << "ApplyElementMatrix called for baseclass" << endl;} - virtual void AddDiagElementMatrix (const Array & dnums1, - const FlatVector & diag, + virtual void AddDiagElementMatrix (FlatArray dnums1, + FlatVector diag, bool inner_element, int elnr, LocalHeap & lh); @@ -835,20 +835,20 @@ namespace ngcomp const Flags & flags); virtual ~T_BilinearFormDiagonal (); - virtual void AllocateMatrix (); - virtual unique_ptr CreateRowVector() const; - virtual unique_ptr CreateColVector() const; + virtual void AllocateMatrix () override; + virtual unique_ptr CreateRowVector() const override; + virtual unique_ptr CreateColVector() const override; virtual void AddElementMatrix (FlatArray dnums1, FlatArray dnums2, BareSliceMatrix elmat, ElementId id, bool addatomic, - LocalHeap & lh); + LocalHeap & lh) override; - virtual void AddDiagElementMatrix (const Array & dnums1, - const FlatVector & diag, + virtual void AddDiagElementMatrix (FlatArray dnums1, + FlatVector diag, bool inner_element, int elnr, - LocalHeap & lh); + LocalHeap & lh) override; /* virtual void ApplyElementMatrix(const BaseVector & x, BaseVector & y,