hplab2D

*********************************************************************
********               hplab2D     (release 1)               ********
*********************************************************************

March 8, 2010
Tomas Vejchodsky
www.math.cas.cz/vejchod/hplab.html

An hp-FEM solver for Matlab. It solves 2D linear elliptic problems in
the form
    -div( a grad u) + c u = f    in Omega,
                        u = 0    on the boundary of Omega
by the hp-version of the finite element method (hp-FEM). This version
of the FEM is characterized by approximation of the solution u by a
piecewise polynomial function with various polynomial degrees on
various elements. Proper distribution of polynomial degrees in the
mesh enables to achieve exponential rate of convergence even for
problems with singularities.

The input are the coefficients a and c, the right-hand side f and the
hp-FEM mesh in the form of three lists: a list of vertices, a list of
triangles, and a list of polynomial degrees. The hp-FEM mesh defines
implicitly the polygonal domain Omega.

Quick Start
===========
- Unpack the archive  hplab2Drel1.tar.gz  with source files.
  Keep everything in the subdirectory hplab2D.
- Run Matlab.
- In Matlab, change current directory to hplab2D.
- Run the test computation.
  >> test1

User's Guide
============
The main function is "hpFEMsolver.m". It has three usual (visible)
input parameters:
  v, t, polydeg
and four hidden inputs:
  m-files: "a_coef.m", "c_coef.m", "f_coef.m", and
  a global variable "hpFEM_intlambda".

The parameters  v, t, polydeg  describe the hp-FEM mesh.
v ... A matrix Nvert-by-2 containing coordinates of the vertices in
      the triangulation. The value Nvert stands for the number of
      vertices.
t ... A matrix Ntri-by-3 defining the triangles in the mesh. Each row
      contains three indices to v.
polydeg ... A vector Ntri-by-1 with polynomial degrees of the
      elements. It is so-called interior polynomial degree as opposed
      to the edge polynomial degree.
      This can be a single number. Then all element have the same
      polynomial degree.
The m-files "a_coef.m", "c_coef.m", "f_coef.m" are in fact real
functions of two variables and they correspond to the coefficients a,
c, and f. They must be functions of two input values x and y
providing a single real number as an output. Internally, however, the
values in barycentres of elements only are evaluated and, hence, the
coefficients a, c, and f, are effectively treated as piecewise
constant.
Finally, the global variable "hpFEM_intlambda" can be ignored on the
input and the code will automatically compute the correct values
(certain integrals) in the run time. This, however, may lead to a
slowdown. It is recommended to initialize the global variable
"hpFEM_intlambda" at the beginning of computation by calling
"loadprecompint.m" (loading the values from a binary file)
or "precompint.m" (computing the values).

The output from "hpFEMsolver.m" is a piecewise polynomial function
uhp stored in an internal structure (see "hpFEMsolver.m" for more
detail). This structure can be used to obtain a value of uhp at any
point from the domain Omega by calling "evaluhploc.m" or to plot
a graph of uhp by calling "plotuhp.m".

For further information see the description directly in the
particular m-files.

The test example introduced in Quick Start describes well the ability
of the code. All data concerning the particular example are stored in
a specifice subdirectory - "ex2_circle" in this case. The script
"test1.m" just calls "ex2solve.m" from the subdirectory "ex2_circle".
The mesh is loaded from the data file "ex2mesh1A.dat" which was
created with the help of the Matlab PDE  toolbox. The PDE toolbox
uses arrays 'p', 'e', and 't' for the description of the
triangulation. The code hplab2D needs vertex array v and triangle
array t only. They can be obtained from the PDE toolbox arrays as:
>> v = p(1:2,:)';
>> t = t(1:3,:)';

Enjoy computation with hplab2D!