Provided by: librheolef-dev_5.93-2_i386 bug

NAME

       mass -- L2 scalar product

SYNOPSIS

               form(const space& V, const space& V, "mass");
               form(const space& M, const space& V, "mass");
               form (const space& V, const space& V, "mass", const domain& gamma);
               form_diag(const space& V, "mass");

DESCRIPTION

       Assembly  the  matrix associated to the L2 scalar product of the finite
       element space V.

                        /
                        |
               m(u,v) = | u v dx
                        |
                        / Omega

       The V space may be either a P0, P1, P2,  bubble,  P1d  and  P1d  finite
       element spaces for building a form see form(3).

       The  use  of  quadrature  formulae  is  sometime  usefull  for building
       diagonal matrix.  These approximate matrix are  eay  to  invert.   This
       procedure is available for P0 and P1 approximations.

       Notes  that  when dealing with discontinuous finite element space, i.e.
       P0 and P1d, the corresponding mass matrix is block  diagonal,  and  the
       inv_mass form may be usefull.

       When  two  different  space  M  and V are supplied, assembly the matrix
       associated to the projection operator from one finite element  space  M
       to space V.

                        /
                        |
               m(q,v) = | q v dx
                        |
                        / Omega

       for all q in M and v in V.

       This  form  is  usefull  for instance to convert discontinuous gradient
       components to a continuous approximation.  The transpose  operator  may
       also be usefull to performs the opposite operation.

       The   following  $V$  and  $M$  space  approximation  combinations  are
       supported for the mass form: P0-P1, P0-P1d, P1d-P2, P1-P1d and P1-P2.

EXAMPLE

       The following piece of code build the mass matrix associated to the  P1
       approximation:

               geo g("square");
               space V(g, "P1");
               form m(V, V, "mass");

       The use of lumped mass form write also:

               form_diag md(V, "mass");

       The following piece of code build the projection form:

               geo g("square");
               space V(g, "P1");
               space M(g, "P0");
               form m(M, V, "mass");

SCALAR PRODUCT ON THE BOUNDARY

       Assembly  the  matrix  associated to the L2 scalar product related to a
       boundary domain of a mesh and  a  specified  polynomial  approximation.
       These  forms are usefull when defining non-homogeneous Neumann or Robin
       boundary conditions.

       Let W be a space of  functions  defined  on  Gamma,  a  subset  of  the
       boundary of the whole domain Omega.

                        /
                        |
               m(u,v) = | u v dx
                        |
                        / Gamma

       for  all  u,  v  in W.  Let V a space of functions defined on Omega and
       gamma the trace operator from V into W.  For all u in W and v in V:

                         /
                         |
               mb(u,v) = | u gamma(v) dx
                         |
                         / Gamma

       For all u and v in V:

                         /
                         |
               ab(u,v) = | gamma(u) gamma(v) dx
                         |
                         / Gamma

EXAMPLE

       The following piece of code  build  forms  for  the  P1  approximation,
       assuming that the mesh contains a domain named boundary:

               geo omega ("square");
               domain gamma = omega.boundary();
               space V  (omega, "P1");
               space W  (omega, gamma, "P1");
               form  m  (W, W, "mass");
               form  mb (W, V, "mass");
               form  ab (V, V, "mass", gamma);

SEE ALSO

       form(3)