View source for Hp-FEM

{{Short description|Generalization of finite element method}}
{{lowercase title}}

'''hp-FEM''' is a generalization of the [[Finite element method|finite element method]] (FEM) for solving [[Partial differential equation|partial differential equation]]s [[Numerical analysis|numerically]] based on [[Piecewise|piecewise-polynomial]] approximations. hp-FEM originates from the discovery by Barna A. Szabó and [[Biography:Ivo Babuška|Ivo Babuška]] that the finite element method converges ''exponentially fast'' when the mesh is refined using a suitable combination of h-refinements (dividing elements into smaller ones) and [[P-FEM|p-refinements]] (increasing their polynomial degree)<ref>B. A. Szabó, A. K. Mehta: p-Convergent Finite Element Approximations in Fracture Mechanics, Int. J. Num. Meth. Engng, Volume 12, pp. 551-560, 1978.</ref><ref>[[Biography:Ivo Babuška|I. Babuška]], B. A. Szabó and I. N. Katz: The p-Version of the Finite Element Method, SIAM J. Numer. Anl., Volume 18, pp. 515-544, 1981.</ref><ref>[[Biography:Ivo Babuška|I. Babuška]], B. A. Szabó, On the Rates of Convergence of the Finite Element Method, Int. J. Numer. Meth.Engng., Volume 18, pp. 323-341, 1982.</ref><ref>[[Biography:Ivo Babuška|I. Babuška]]: The p- and hp-Versions of the Finite Element Method: the State of the Art, Finite Elements: Theory and Applications, edited by D. L. Dwoyer, M. Y. Hussaini and R. G. Voigt, New York, Springer-Verlag, 1988.</ref><ref>B. A. Szabó, [[Biography:Ivo Babuška|I. Babuška]]: Finite Element Analysis, John Wiley & Sons, {{ISBN|978-0-471-50273-9}}, 1991.</ref><ref>[[Biography:Ivo Babuška|I. Babuška]], B.Q. Guo: The h, p and h-p version of the finite element method: basis theory and applications, Advances in Engineering Software, Volume 15, Issue 3-4, 1992.</ref>.The exponential convergence of hp-FEM has been observed by numerous independent researchers.<ref>J.M. Melenk: hp-Finite Element Methods for Singular Perturbations, Springer, 2002</ref><ref>C. Schwab: p- and hp- Finite Element Methods: Theory and Applications in Solid and Fluid Mechanics, Oxford University Press, 1998</ref><ref>P. Solin: Partial Differential Equations and the Finite Element Method, J. Wiley & Sons, 2005</ref>

==Differences from standard FEM==
The hp-FEM differs from the standard (lowest-order) FEM in many aspects.<ref>P. Solin, K. Segeth, I. Dolezel: Higher-Order Finite Element Methods, Chapman & Hall/CRC Press, 2003</ref>
* '''Choice of higher-order shape functions'''{{Examples|date=September 2011}}: The higher-degree polynomials in elements can be generated using different sets of shape functions. The choice of such a set can influence dramatically the conditioning of the stiffness matrix, and in turn the entire solution process.<ref>I. Babuska, M. Griebel and J. Pitkaranta, The problem of selecting the shape functions for a p-type finite element, Internat. J. Numer. Methods Engrg. (1989), pp. 1891–1908</ref>
* '''Automatic hp-adaptivity''': In hp-FEM, an element can be hp-refined in many different ways, such as: Increasing its polynomial degree without subdividing it in space, or subdividing the element geometrically, where various polynomial degrees can be applied to the sub-elements. The number of element refinement candidates easily reaches 100 in two dimensions and 1000 in three dimensions. One number indicating the size of the error in an element is not enough to guide automatic hp-adaptivity (as opposed to adaptivity in standard FEM). Other techniques such as ''reference solutions'' or ''analyticity considerations'' must be employed to obtain more information about the ''shape of error'' in every element.<ref>L. Demkowicz, W. Rachowicz, and Ph. Devloo:  A Fully Automatic hp-Adaptivity, Journal of Scientific Computing, 17, Nos 1–3 (2002), 127–155</ref>
* '''The ratio of assembling and solution CPU times''': In standard FEM, the stiffness matrix usually is assembled quickly but it is quite large. Typically, the solution of the discrete problem consumes the largest part of the overall computing time. By contrast, the stiffness matrices in hp-FEM typically are much smaller, but (for the same matrix size) their assembly takes more time than in standard FEM. This is mainly due to the computational cost of numerical quadrature, which must have higher precision, and therefore be of higher order, compared to standard FEM to take advantage of the faster convergence rates.
* '''Programming challenges''': It is much harder to implement a hp-FEM solver than standard FEM code. The multiple issues that need to be overcome include (but are not limited to): higher-order quadrature formulas, higher-order shape functions, connectivity and orientation information relating shape functions on the reference domain with basis functions in the physical domain, etc.<ref>L. Demkowicz, J. Kurtz, D. Pardo, W. Rachowicz, M. Paszynski, A. Zdunek: Computing with hp-Adaptive Finite Elements, Chapman & Hall/CRC Press, 2007</ref>

==Solution to Fichera problem==
The Fichera problem (also called the Fichera corner problem) is a standard benchmark problem for adaptive FEM codes. One can use it to show the dramatic difference in the performance of standard FEM and hp-FEM. The problem geometry is a cube with missing corners The exact solution has a singular gradient (an analogy of infinite stress) at the centre. The knowledge of the exact solution makes it possible to calculate the approximation error exactly and thus compare various numerical methods. For illustration, the problem was solved using three different versions of adaptive FEM: linear elements, quadratic elements, and hp-FEM.

<gallery>
Image:grad fichera.png|Fichera problem: singular gradient.
Image:conv fichera.png|Fichera problem: convergence comparison.
</gallery>

The convergence graphs show the approximation error as a function of the number of [[Degrees of freedom|degrees of freedom]] (DOF). DOF refers to unknown parameters that are needed to define the approximation, and the number of DOF equals the size of the stiffness matrix. The reader can see in the graphs that the convergence of the hp-FEM is much faster than the convergence of both other methods. The performance gap is large enough that the linear FEM might not converge at all (in reasonable time) and the quadratic FEM would need hundreds of thousands or perhaps millions of DOF to reach the accuracy that hp-FEM attained with approximately 17,000 DOF. Obtaining very accurate results using relatively few degrees of freedom is the main strength of hp-FEM.

==Efficiency of hp-FEM==
Smooth functions can be approximated much more efficiently using large high-order elements than small piecewise-linear ones. This is illustrated in the figure below, where a one-dimensional Poisson equation with zero Dirichlet boundary conditions is solved on two different meshes. The exact solution is the sine function.

* Left: mesh consisting of two linear elements.
* Right: mesh consisting of one quadratic element.

[[Image:grad sin h.png|Piecewise-linear approximation.]]
[[Image:conv sin hp.png|Quadratic approximation.]]

While the number of unknowns is the same in both cases (1 DOF), the errors in the corresponding norm are 0.68 and 0.20, respectively. This means that the quadratic approximation was roughly 3.5 times more efficient than the piecewise linear one. When we proceed one step further and compare (a) four linear elements to (b) one quartic element (p=4), then both discrete problems will have three DOFs, but the quartic approximation will be approximately 40 times more efficient.

On the contrary, small low-order elements can capture small-scale features such as singularities much better than large high-order ones. hp-FEM is based on an optimal combination of these two approaches which leads to exponential convergence. Note that this exponential convergence is expressed in the axis of error vs. degrees of freedom. For real-life applications, we usually consider the computational time needed to reach the same level of accuracy. For this performance indicator h- and hp-refinement can provide similar results.<ref>{{Cite web|url=http://hpfem.org/wp-content/uploads/doc-web/doc-examples/src/hermes2d/examples/maxwell/microwave-oven.html|title=Microwave Oven — Hermes Examples Guide|access-date=2019-06-26|archive-date=2019-06-26|archive-url=https://web.archive.org/web/20190626013508/http://hpfem.org/wp-content/uploads/doc-web/doc-examples/src/hermes2d/examples/maxwell/microwave-oven.html|url-status=live|language=en}}</ref> As soon as it is harder to program and [[Parallel computing|parallelize]] hp-FEM compared to h-FEM, the convergence excellence of hp-refinement may become impractical.

==Hp-adaptivity==
Some FEM sites describe hp-adaptivity as a combination of h-adaptivity (splitting elements in space while keeping their polynomial degree fixed) and p-adaptivity (only increasing their polynomial degree).{{Citation needed|date=June 2022}} This is not entirely accurate, as hp-adaptivity is significantly different from both h- and p-adaptivity as the hp-refinement of an element can be done in many different ways. Besides a p-refinement, the element can be subdivided in space (as in h-adaptivity), but there are many combinations for the polynomial degrees on the sub-elements. This is illustrated in the figure on the right. For example, if a triangular or quadrilateral element is subdivided into four sub-elements where the polynomial degrees are allowed to vary by at most two, then this yields 3^4 = 81 refinement candidates (not considering polynomials anisotropic candidates). Analogously, splitting a hexahedron into eight sub-elements and varying their polynomial degrees by at most two yields 3^8 = 6,561 refinement candidates. Standard FEM error estimates providing one constant number per element is not enough to guide automatic hp-adaptivity.

==Higher-order shape functions==
In standard FEM one only works with shape functions associated with grid vertices (the so-called ''vertex functions''). In contrast, when using hp-FEM, one moreover regards ''edge functions'' (associated with element edges), ''face functions'' (corresponding to element faces – 3D only),
and ''bubble functions'' (higher-order polynomials which vanish nonelement boundaries). The following images show these functions restricted to a single element. All these functions are defined in the entire element interior.

<gallery>
Image:vertex new.jpg|Vertex function.
Image:edge new.jpg|Edge function.
Image:face new.jpg|Face function.
Image:bubble new.jpg|Bubble function.
</gallery>

==Open source hp-FEM codes==
* [[Software:Deal.II|Deal.II]]: Deal.II is a free, open-source library to solve partial differential equations using the finite element method.
* [http://www.concepts.math.ethz.ch/ Concepts]: C/C++ hp-FEM/DGFEM/BEM library for elliptic equations developed at SAM, [[Organization:ETH Zurich|ETH Zurich]] (Switzerland), and in the group of K. Schmidt at [[Organization:Technical University of Berlin|TU Berlin]] (Germany).{{Dead link|date=November 2023}}
* 2dhp90, 3dhp90: Fortran codes for elliptic problems and Maxwell's equations developed by L. Demkowicz at ICES, UT Austin.
* PHAML: The Parallel Hierarchical Adaptive Multi-level Project. Finite element software developed at the [[Organization:National Institute of Standards and Technology|National Institute for Standards and Technology]], USA, for the numerical solution of 2D elliptic partial differential equations on distributed memory parallel computers and multi-core computers using adaptive mesh refinement and multi-grid solution techniques.
* [[Software:Hermes Project|Hermes Project]]: C/C++/Python library for rapid prototyping of space- and space-time adaptive hp-FEM solvers for a large variety of PDEs and Multiphysics PDE systems, developed by the hp-FEM group at the University of Nevada, Reno (USA), Institute of Thermo-mechanics, Prague (Czech Republic), and the University of West Bohemia in Pilsen (Czech Republic) – with the [[Software:Agros2D|Agros2D]] engineering software built on top of the Hermes library.
* [http://lsec.cc.ac.cn/phg/index_en.htm PHG]: PHG is a toolbox for developing parallel adaptive finite element programs. It's suitable for h-, p- and hp-fem. PHG is currently under active development at the State Key Laboratory of Scientific and Engineering Computing, Institute of Computational Mathematics, and Scientific/Engineering Computing of the Chinese Academy of Sciences(LSEC, CAS, China). PHG deals with conforming tetrahedral meshes and uses bisection for adaptive local mesh refinement and MPI for message passing. PHG has an object-oriented design that hides parallelization details and provides common operations on meshes and finite element functions in an abstract way, allowing the users to concentrate on their numerical algorithms.
* [http://mofem.eng.gla.ac.uk MoFEM] is a finite element analysis code tailored for the solution of multi-physics problems with arbitrary levels of approximation and different levels of mesh refinement and optimized for high-performance computing. It is designed to be able to manage complexities related to a heterogeneous order of approximations for  L2, H1, H-div, and H-curl spaces.
* [http://www.sparselizard.org Sparselizard] is a multi-physics, hp-adaptive, user-friendly, open-source C++ finite element library currently developed at Tampere University, Finland. It combines 3D tetrahedral & 2D triangle/quadrangle conformal adaptive mesh refinement with arbitrary order hierarchical H1 & H-curl function spaces for general static and transient hp-FEM.

==Commercial hp-FEM software==
* [https://www.esrd.com/ StressCheck] is a finite element analysis tool with capabilities oriented toward detailed structural analysis.

==References==
{{reflist|30em}}
{{Numerical PDE}}

[[Category:Finite element method]]

{{Sourceattribution|Hp-FEM}}