NGSolve
5.3
|
This is an overview of some representative functions of NGSolve.
For instruction purposes, the functions are simplified.
Element matrix assembling for BDB integrators:
Many bilinear-forms are of the structure
where is some differential operator, i.e., the gradient, and
is the material tensor. Then, the finite elemen matrix
is computed as
The elements of the matrix are the values of the differential operator applied to the shape functions in the integration point, i.e.,
.
The BDBIntegrator is a derived class from BilinearFormIntegrator. Here, the differential operator and the material tensor are defined by template arguments.
template <class DIFFOP, class DMATOP> void BDBIntegrator<DIFFOP,DMATOP> :: AssembleElementMatrix (const FiniteElement & fel, const ElementTransformation & eltrans, Matrix & elmat, LocalHeap & locheap) const { enum { DIM_SPACE = DIFFOP::DIM_SPACE }; enum { DIM_ELEMENT = DIFFOP::DIM_ELEMENT }; enum { DIM_DMAT = DIFFOP::DIM_DMAT }; enum { DIM = DIFFOP::DIM }; int ndof = fel.GetNDof(); elmat.AssignMemory (ndof*DIM, ndof*DIM, locheap); elmat = 0; MatrixFixHeight<DIM_DMAT> bmat (ndof * DIM, locheap); MatrixFixHeight<DIM_DMAT> dbmat (ndof * DIM, locheap); Mat<DIM_DMAT,DIM_DMAT> dmat; const IntegrationRule & ir = GetIntegrationRule (fel); // optimized memory management void * heapp = locheap.GetPointer(); for (int i = 0; i < ir.GetNIP(); i++) { // the mapped point, including Jacobian SpecificIntegrationPoint<DIM_ELEMENT,DIM_SPACE> sip(ir, i, eltrans, locheap); DIFFOP::GenerateMatrix (fel, sip, bmat, locheap); dmatop.GenerateMatrix (fel, sip, dmat, locheap); double fac = fabs (sip.GetJacobiDet()) * sip.IP().Weight(); dbmat = dmat * bmat; elmat += fac * (Trans (bmat) * dbmat); locheap.CleanUp (heapp); } }
Global matrix assembling The bilinear-form consists of a couple of integrators acting on the domain or on the boundary. It has a reference to the finite element space it is defined on. The sparse matrices are stored in the bilinear-form object.
template <typename SCAL = double> void BilinearForm :: Assemble () { Array<int> dnums; Matrix<SCAL> Matrix elmat; ElementTransformation eltrans; int ndof = fespace.GetNDof(); BaseMatrix & mat = GetMatrix(); mat = 0.0; LocalHeap locheap (1000000); // assembling of volume terms int ne = ma.GetNE(); for (int i = 0; i < ne; i++) { locheap.CleanUp(); if (!fespace.DefinedOn (ma.GetElIndex (i))) continue; ma.GetElementTransformation (i, eltrans); const FiniteElement & fel = fespace.GetFE (i); fespace.GetDofNrs (i, dnums); for (int j = 0; j < parts.Size(); j++) { const BilinearFormIntegrator & bfi = *parts.Get(j); if (bfi.BoundaryForm()) continue; bfi.AssembleElementMatrix (fel, eltrans, elmat, &locheap); fespace.TransformMatrix (i, elmat); mat->AddElementMatrix (dnums, elmat); } } // assembling of surface terms int nse = ma.GetNSE(); for (int i = 0; i < nse; i++) { locheap.CleanUp(); if (!fespace.DefinedOnBoundary (ma.GetSElIndex (i))) continue; ma.GetSurfaceElementTransformation (i, eltrans); const FiniteElement & fel = fespace.GetSFE (i); fespace.GetSDofNrs (i, dnums); for (int j = 1; j <= parts.Size(); j++) { const BilinearFormIntegrator & bfi = *parts.Get(j); if (!bfi.BoundaryForm()) continue; bfi.AssembleElementMatrix (fel, eltrans, elmat, &locheap); fespace.TransformSurfMatrix (i, elmat); mat->AddElementMatrix (dnums, elmat); } } }