Solving anisotropic Poisson problems via discrete exterior calculus

<!-- metadata commented in wiki content
==Anisotropy in 2D Discrete Exterior Calculus==

''''''
-->
==Abstract==

We present a local formulation for 2D Discrete Exterior Calculus (DEC) similar to that of the Finite Element Method (FEM), which allows  a natural treatment of material heterogeneity (assigning material properties element by element). It also allows us to deduce, in a principled  manner, anisotropic fluxes and the DEC discretization of the pullback of 1-forms by the anisotropy tensor, i.e. we deduce the discrete action of the anisotropy tensor on primal 1-forms. Due to the local formulation,  the computational cost of DEC is similar to that of the Finite Element Method with Linear interpolating functions (FEML). The numerical  DEC solutions to the anisotropic Poisson equation show numerical convergence, are very close to those of FEML on fine meshes  and are slightly better than those of FEML on coarse meshes.
<!--
<span id="fn-1"></span>
<span style="text-align: center; font-size: 75%;">([[#fnc-1|<sup>1</sup>]]) Centro de Investigación en Matemáticas, Calle Jalisco s/n,  Guanajuato, GTO 36023, México. Email: esqueda,rherrera,botello@cimat.mx </span>

<span id="fn-2"></span>
<span style="text-align: center; font-size: 75%;">([[#fnc-2|<sup>2</sup>]]) Departamento de  Matemáticas, Universidad de Guanajuato,  Guanajuato, GTO 36023, México. Email: valerocar@gmail.com</span>
-->

==1. Introduction==

The theory of Discrete Exterior Calculus (DEC) is a relatively recent discretization  <span id='citeF-8'></span>[[#cite-8|[8]]] of the classical theory of  Exterior Differential Calculus developed by  E. Cartan <span id='citeF-3'></span>[[#cite-3|[3]]], which is a fundamental tool in Differential Geometry and Topology.  The aim of DEC is to solve partial differential equations preserving their geometrical and physical features as much as possible. There are only a few papers about implementions of DEC to solve certain PDEs, such as  the Darcy flow and Poisson's equation <span id='citeF-9'></span>[[#cite-9|[9]]],  the Navier-Stokes equations <span id='citeF-10'></span>[[#cite-10|[10]]], the simulation of  elasticity, plasticity and failure of isotropic materials <span id='citeF-5'></span>[[#cite-5|[5]]], some comparisons with  the finite differences and finite volume methods on regular flat meshes <span id='citeF-7'></span>[[#cite-7|[7]]], as well as applications in digital geometry processing <span id='citeF-4'></span>[[#cite-4|[4]]].

In this paper, we describe a local formulation of DEC which is reminiscent of that of the Finite Element Method  (FEM). Indeed,  once the local systems of equations have been established, they can be assembled into a global linear system. This local formulation is also efficient and helpful in understanding various features of DEC that can otherwise remain  unclear if one is dealing dealing with an entire mesh. Besides, we believe the local description to DEC will be accesible to a wide readership. We will, therefore, take a local approach when recalling all the objects required by 2D DEC <span id='citeF-6'></span>[[#cite-6|[6]]]. Our main results are the following:

* We present a local formulation of DEC analogous to that of FEM, which allows a natural treatment of heterogeneous material properties assigned to subdomains (element by element) and  eliminates the need of dealing with it through ad hoc modifications of the global discrete Hodge star operator.

* Guided by the local formulation, we deduce a natural way to approximate the flux/gradient-vector of a discretized function, as well as the anisotropic flux vector. We carry out a comparison of the formulas defining the flux in both DEC and Finite Element Method with linear interpolation functions (FEML).

* In order to understand how to discretize the anisotropic Poisson equation, we develop the discretization of the pull-back operator of 1-forms under an arbitrary linear trasformation or tensor. The pull-back operator is one of the basic ingredients in Exterior Differential Calculus in the smooth setting and is essential in all of its applications  in Topology [bott-Tu].  We discretize the ''pullback operator on primal 1-forms induced by an arbitrary tensor'' using Whitney interpolation (the Whitney map). The Whitney interpolation forms were introduced by Hassler Whitney in 1957 <span id='citeF-12'></span>[[#cite-12|[12]]] and  Bossavit <span id='citeF-1'></span>[[#cite-1|[1]]] explained their relevance in “mixed methods” of finite elements.  To our knowledge, this is the first time that the discrete version of this operator is presented in the DEC literature.   This has allowed us to discretize, in a principled manner, the anisotropic heat equation.

*  We carry out an analytic comparison of the DEC and FEML local formulations of the anisotropic Poisson equation.

* We present three numerical examples of the approximate solutions to the stationary anisotropic Poisson equation on different domains using DEC and FEML. The numerical DEC-solutions exhibit numerical convergence (see the error measurement tables) and a competitive performance, as well as a computational cost similar to that of FEML. In fact, the numerical solutions with both methods on fine meshes are identical, and DEC shows a slightly better performance than FEML on coarse meshes.


The paper is organized as follows. In Section [[#2 Preliminaries on DEC from a local viewpoint|2]], we describe the local versions of the discrete derivative operator, the dual mesh, the discrete Hodge star operator and the meaning of a continuous 1-form on the plane.  In Section [[#3 Pullback operator on primal 1-forms|3]], we deduce the discretization of the pullback operator on primal 1-forms. In Section [[#4 Flux and anisotropy|4]], we deduce the natural way of computing flux vectors in DEC  (which turns out to be equivalent to the FEML result), as well as the anisotropic flux vectors.  In Section [[#5 Anisotropic Poisson equation in 2D|5]],  we present the local DEC formulation of the 2D anisotropic Poisson equation and compare it with the local system of FEML, proving that the diffusion terms  are identical in both schemes, while the source terms are differentl due to a different area-weight assignment for the nodes.  In Section [[#6 Some remarks about DEC quantities|6]], we re-examine the geometry of some of the local DEC quantities. In Section [[#7 Numerical Examples|7]], we present and compare numerical examples of DEC and FEML approximate solutions to the 2D anisotropic  Poisson equation on different domains with meshes of various resolutions.  In Section [[#8 Conclusions|8]], we summarize the contributions of this paper.
<!--
''Acknowledgements''. The second named author was partially supported by a CONACyT grant, and would like to thank  the International Centre for Numerical Methods in Engineering (CIMNE) and the University of Swansea for their hospitality. We gratefully acknowledge the support of NVIDIA Corporation with the donation of the Titan X Pascal GPU used for this research.
-->

==2. Preliminaries on DEC from a local viewpoint==

In this section, we will recall the basic operators of Dicrete Exterior Calculus  restricting ourselves to a mesh made up of one simplex/triangle. The local results we derive in the  paper can be assembled, just as in the Finite elemnt Method,  due to the additivity of both differentiation and integration.

Let us consider a primal mesh made up of a single (positively oriented) triangle with vertices <math display="inline">v_1,v_2,v_3</math> (Figure 1).

<div id='img-1'></div>
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: auto;"
|-
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-Triangle01.png|228px|Triangle ]]
|- style="text-align: center; font-size: 75%;"
| colspan="1" style="padding-bottom:10px;"| '''Figure 1'''. Triangle 
|}

Such a mesh has

* one oriented 2-dimensional face

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>[v_1,v_2,v_3];</math>
|}
|}
* three oriented 1-dimensional edges

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>[v_1,v_2],\quad [v_1,v_3],\quad [v_2,v_3];</math>
|}
|}
* and three 0-dimensional vertices

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>[v_1],\quad [v_2],\quad [v_3].</math>
|}
|}

===2.1 Boundary operator===

There is a well known boundary operator

<span id="eq-1"></span>
{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>\partial _{2,1} [v_1,v_2,v_3]=[v_2,v_3]-[v_1,v_3]+[v_1,v_2], </math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (1)
|}

which describes the boundary of the triangle as an alternated sum of its ordered oriented edges <math display="inline">[v_1,v_2]</math>, <math display="inline">[v_1,v_3]</math> and <math display="inline">[v_2,v_3]</math>.

Similarly, one can compute the boundary of each edge

<span id="eq-2"></span>
{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>\partial _{1,0} [v_1,v_2] = [v_2]-[v_1], </math>
|-
| style="text-align: center;" | <math> \partial _{1,0} [v_1,v_3] = [v_3]-[v_1], </math>
| style="width: 5px;text-align: right;white-space: nowrap;" | (2)
|-
| style="text-align: center;" | <math> \partial _{1,0} [v_2,v_3] = [v_3]-[v_2].  </math>
|}
|}

If we consider

* the symbol <math display="inline">[v_1,v_2,v_3]</math> as a basis vector of a 1-dimensional vector space,
* the symbols <math display="inline">[v_1,v_2]</math>, <math display="inline">[v_1,v_3]</math>, <math display="inline">[v_2,v_3]</math> as an ordered basis of a 3-dimensional vector space,
* the symbols <math display="inline">[v_1]</math>, <math display="inline">[v_2]</math>, <math display="inline">[v_3]</math> as an ordered basis of a 3-dimensional vector space,


then the map ([[#eq-1|1]]), which sends the oriented triangle to a sum of its oriented edges, is represented by the matrix

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>[\partial _{2,1}]=\left(\begin{array}{r} 1 \\  -1 \\  1 \end{array}\right),</math>
|}
|}

while the map ([[#eq-2|2]]), which sends the oriented edges to sums of their oriented vertices, is represented by the matrix

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>[\partial _{1,0}]=\left(\begin{array}{rrr} -1 & -1 & 0 \\  1 & 0 & -1 \\  0 & 1 & 1 \end{array} \right).</math>
|}
|}

===2.2 Discrete derivative===

It has been argued that the DEC discretization of the differential of a function is given by the transpose of the matrix of the  boundary operator on edges  <span id='citeF-8'></span><span id='citeF-6'></span>[[#cite-8|[6,8]]]. More precisely, suppose we have a function discretized by its values at the vertices

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>f\sim \left(\begin{array}{l} f_1\\ f_2\\ f_3               \end{array} \right).</math>
|}
|}

Its discrete derivative, according to DEC, is

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>\left(\begin{array}{rrr}-1 & -1 & 0 \\  1 & 0 & -1 \\  0 & 1 & 1 \end{array} \right)^T \left(\begin{array}{l}f_1\\ f_2\\ f_3               \end{array} \right) =\left(\begin{array}{rrr}-1 & 1 & 0 \\  -1 & 0 & 1 \\  0 & -1 & 1 \end{array} \right) \left(\begin{array}{l}f_1\\ f_2\\ f_3               \end{array} \right)</math>
|-
| style="text-align: center;" | <math> =\left(\begin{array}{l}f_2-f_1\\ f_3-f_1\\ f_3-f_2               \end{array} \right). </math>
|}
|}

Indeed, such differences are rough approximations of the directional derivatives of <math display="inline">f</math> along the oriented edges.  For instance, <math display="inline">f_2-f_1</math> is a rough approximation of the directional derivative of <math display="inline">f</math> at <math display="inline">v_1</math> in the direction of the vector <math display="inline">v_2-v_1</math>, i.e.

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>f_2-f_1 \approx df_{v_1}(v_2-v_1) .</math>
|}
|}

It is precisely in this sense that, according to DEC,

* the value <math display="inline">f_2-f_1</math> is assigned to the edge <math display="inline">[v_1,v_2]</math>,
* the value <math display="inline">f_3-f_1</math> is assigned to the edge <math display="inline">[v_1,v_3]</math>,
* and the value <math display="inline">f_3-f_2</math> is assigned to the edge <math display="inline">[v_2,v_3]</math>.

Let

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>D_0:=\left(\begin{array}{rrr} -1 & 1 & 0 \\  -1 & 0 & 1 \\  0 & -1 & 1 \end{array} \right). </math>
|}
|}

===2.3 Dual mesh===

The dual mesh of the primal mesh consisting of a single triangle is constructured as follows:

* To the 2-dimensional triangular face <math display="inline">[v_1,v_2,v_3]</math> will correspond the 0-dimensional point given by the circumcenter <math display="inline">c</math> of the triangle (Figure 2).

<div id='img-2'></div>
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: auto;"
|-
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-Triangle02.png|240px|Circumcenter ]]
|- style="text-align: center; font-size: 75%;"
| colspan="1" style="padding-bottom:10px;"| '''Figure 2'''. Circumcenter 
|}

* To the 1-dimensional edge <math display="inline">[v_1,v_2]</math> will correspond the 1-dimensional straight line segment <math display="inline">[p_1,c]</math> joining the midpoint <math display="inline">p_1</math> of the edge <math display="inline">[v_1,v_2]</math> to the circumcenter <math display="inline">c</math> (Figure 3).   Similarly for the other edges.

<div id='img-3'></div>
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: auto;"
|-
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-Triangle03.png|240px|Dual segment ]]
|- style="text-align: center; font-size: 75%;"
| colspan="1" style="padding-bottom:10px;"| '''Figure 3'''. Dual segment 
|}

* To the 0-dimensional vertex/node <math display="inline">[v_1]</math> will correspond the oriented <math display="inline">2</math>-dimensional quadrilateral <math display="inline">[v_1,p_1,c,p_3]</math> (Figure 4).

<div id='img-4'></div>
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: bottom;"
|-
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-Triangle04.png|240px|Dual quadrilateral  ]]
|- style="text-align: center; font-size: 75%;"
| colspan="1" style="padding-bottom:10px;"| '''Figure 4'''. Dual quadrilateral  
|}

===2.4 Discrete Hodge star===

For the Poisson equation in 2D, we need two matrices: one relating original edges to dual edges, and another relating vertices to dual cells.

* The discrete Hodge star map <math display="inline">M_1</math> applied to the discrete differential of a discretized function <math display="inline">f\sim (f_1,f_2,f_3)</math>  is given as follows:

* the value <math display="inline">f_2-f_1</math> assigned to the edge <math display="inline">[v_1,v_2]</math> is changed to the new value

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>{{length}[p_1,c]\over {length}[v_1,v_2]}(f_2-f_1)</math>
|}
|}

assigned to the segment <math display="inline">[p_1,c]</math>;
* the value <math display="inline">f_3-f_1</math> assigned to the edge <math display="inline">[v_1,v_3]</math> is changed to the new value

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>{{length}[p_3,c]\over {length}[v_1,v_3]}(f_3-f_1)</math>
|}
|}

assigned to the edge <math display="inline">[p_3,c]</math>.
* the value <math display="inline">f_3-f_2</math> assigned to the edge <math display="inline">[v_2,v_3]</math> is changed to the new value

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>{{length}[p_2,c]\over {length}[v_2,v_3]}(f_3-f_2)</math>
|}
|}

assigned to the segment <math display="inline">[p_2,c]</math>;

In other words,

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>M_1=\left(\begin{array}{ccc} \displaystyle{{length}[p_1,c]\over {length}[v_1,v_2]} & 0 & 0 \\  0 & \displaystyle{{length}[p_3,c]\over {length}[v_1,v_3]} & 0 \\  0 & 0 & \displaystyle{{length}[p_2,c]\over {length}[v_2,v_3]} \end{array} \right).</math>
|}
|}

* Similarly, the discrete Hodge star map <math display="inline">M_0</math> on values on vertices is given as follows:

* the value <math display="inline">f_1</math> assigned to the vertex <math display="inline">[v_1]</math> is changed to the new value

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>{Area}[v_1,p_1,c,p_3]f_1</math>
|}
|}

assigned to the quadrilateral <math display="inline">[v_1,p_1,c,p_3]</math>;
* the value <math display="inline">f_2</math> assigned to the vertex <math display="inline">[v_2]</math> is changed to the new value

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>{Area}[v_2,p_2,c,p_1]f_2</math>
|}
|}

assigned to the quadrilateral <math display="inline">[v_2,p_2,c,p_1]</math>;
* the value <math display="inline">f_3</math> assigned to the vertex <math display="inline">[v_3]</math> is changed to the new value

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>{Area}[v_3,p_3,c,p_2]f_2</math>
|}
|}

assigned to the quadrilateral <math display="inline">[v_3,p_3,c,p_2]</math>.

In other words,

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>M_0=\left(\begin{array}{ccc} {Area}[v_1,p_1,c,p_3] & 0 & 0 \\  0 & {Area}[v_2,p_2,c,p_1] & 0 \\  0 & 0 & {Area}[v_3,p_3,c,p_2] \end{array} \right).</math>
|}
|}

===2.5 Continuous 1-forms===

A differential 1-form <math display="inline">\alpha </math> on <math display="inline">\mathbb{R}^2</math> can be thought of as a function whose arguments are  a point in <math display="inline">\mathbb{R} ^2</math> and a vector in <math display="inline">\mathbb{R} ^2</math>, i.e.

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>\alpha : \mathbb{R}^2 \times \mathbb{R}^2 \longrightarrow \mathbb{R},</math>
|}
|}

and which is linear on the vector argument.

The dot product of <math display="inline">\mathbb{R}^2</math> allows us to see continuous 1-forms in a more familiar form as vector fields in the following way: If <math display="inline">X</math> is a vector field on <math display="inline">\mathbb{R}^2</math>, at each point <math display="inline">p\in \mathbb{R}^2</math> we can consider the  1-form

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>(p,Y)\mapsto   \alpha (p,Y):=X(p)\cdot Y, </math>
|}
|}

where <math display="inline">Y\in \mathbb{R}^2</math>. On the other hand, from a continuos 1-form we can build a vector field as follows

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>X(p) := \alpha (p,e_1) e_1 + \alpha (p,e_2) e_2 = (\alpha (p,e_1),\alpha (p,e_2)),</math>
|}
|}

where

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math> e_1=\left(\begin{array}{c} 1 \\  0 \end{array} \right) \quad  \mbox{and} \quad  e_2=\left(\begin{array}{c} 0 \\  1 \end{array} \right). </math>
|}
|}

====2.5.1 Pullback of a conitnuous 1-form====

Given a tensor <math display="inline">S</math> (a linear transformation) and the vector field <math display="inline">X</math> associated to a 1-form <math display="inline">\alpha </math> as above, consider the vector field <math display="inline">S(X)</math>. For any vector <math display="inline">Y</math>, the value of the corresponding continuous 1-form is given by the product

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>S(X(p))\cdot Y     =X(p)\cdot S^T Y </math>
|-
| style="text-align: center;" | <math>    =\alpha (p, S^T Y) </math>
|-
| style="text-align: center;" | <math>    =(S^T)^*\alpha (p, Y),  </math>
|}
|}

which is the value of the 1-form called ''the pullback of <math>\alpha </math> by the tensor <math>S^T</math>''.

====2.5.2 Discretization of continuous 1-forms====

In order to produce a primal 1-form from such a continuous 1-form we need to inegrate it along the edges of the mesh (the de Rham map). In our local description, this produces a collection of 3 numbers associated to the 3 oriented edges of the triangle, i.e. for the oriented edge joining <math display="inline">v_i</math> to <math display="inline">v_j</math> we have

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>\alpha _{ij}:=\int _0^1 \alpha (v_i+t(v_j-v_i), v_j-v_i) dt.</math>
|}
|}

This is equivalent to calculating the line integral (circulation) of the corresponding vector field <math display="inline">X</math> along the path <math display="inline">\gamma _{ij}:=v_i+t(v_j-v_i)</math>, <math display="inline">0\leq t\leq 1</math>, parametrizing the edge <math display="inline">[v_i,v_j]</math>

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>\alpha _{ij}:=\int _{\gamma _{ij}} X. </math>
|}
|}

==3. Pullback operator on primal 1-forms==

In Section [[#2 Preliminaries on DEC from a local viewpoint|2]], we dealt with the primal discretization of a continuous 1-form. Now we wish to understand the pullback of a primal 1-form. In order to do this, we will interpolate a primal 1-form using the Whitney interpolation forms to produce a  continuous 1-form (Whitney map) to which we can apply the pullback operator by a tensor (linear transformation), and then intergrate along the edges (de Rham map).

The Whitney interpolation forms were introduced by  Whitney in 1957 <span id='citeF-12'></span>[[#cite-12|[12]]], and Bossavit <span id='citeF-1'></span>[[#cite-1|[1]]] explained their relevance to “mixed methods” of finite elements.  For the sake of simplicity, we will only define the Whitney 0 and 1-forms.

First consider the Finite Element linear basis functions (or Whitney interpolation functions)

<span id="eq-3"></span>
{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>N_1 = {1\over 2A}[(y_2-y_3)x+(x_3-x_2)y+x_2y_3-x_3y_2] ,</math>
|-
| style="text-align: center;" | <math> N_2 = {1\over 2A}[(y_3-y_1)x+(x_1-x_3)y+x_3y_1-x_1y_3] ,</math>
| style="width: 5px;text-align: right;white-space: nowrap;" | (3)
|-
| style="text-align: center;" | <math> N_3 = {1\over 2A}[(y_1-y_2)x+(x_2-x_1)y+x_1y_2-x_2y_1] , </math>
|}
|}

and their differentials

<span id="eq-4"></span>
{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>dN_1 = {1\over 2A}[(y_2-y_3)dx+(x_3-x_2)dy] ,</math>
|-
| style="text-align: center;" | <math> dN_2 = {1\over 2A}[(y_3-y_1)dx+(x_1-x_3)dy] ,</math>
| style="width: 5px;text-align: right;white-space: nowrap;" | (4)
|-
| style="text-align: center;" | <math> dN_3 = {1\over 2A}[(y_1-y_2)dx+(x_2-x_1)dy] . </math>
|}
|}

The Whitney interpolation 1-forms are

<span id="eq-5"></span>
{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>\phi _{12}=N_1dN_2-N_2dN_1,</math>
|-
| style="text-align: center;" | <math> \phi _{13}=N_1dN_3-N_3dN_1 ,</math>
| style="width: 5px;text-align: right;white-space: nowrap;" | (5)
|-
| style="text-align: center;" | <math> \phi _{23}=N_2dN_3-N_3dN_2. </math>
|}
|}

These differential 1-forms are such that, along the edge <math display="inline">[v_k,v_l]</math>

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>\int _{v_k}^{v_l}\phi _{ij} = \delta _{ik}\delta _{jl}-\delta _{il}\delta _{jk},</math>
|}
|}

i.e. <math display="inline">0</math> or <math display="inline">\pm 1</math>.

Now, suppose we have a primal 1-form given on our triangle by the collection of numbers

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>\beta _{12},\quad \beta _{13},\quad \beta _{23}.</math>
|}
|}

We form a continuous 1-form as follows

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>\beta :=\beta _{12}\phi _{12} + \beta _{13}\phi _{13} + \beta _{23}\phi _{23} </math>
|}
|}

which is such that, for <math display="inline">1\leq i <j \leq 3</math>

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>\int _{v_i}^{v_j} \beta = \beta _{ij}. </math>
|}
|}

Now we will calculate

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>\int _{v_i}^{v_j} S^*\beta </math>
|}
|}

for an arbitrary linear transformation

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>S:= \left( \begin{array}{cc} a & b \\  c & d \end{array}  \right)</math>
|}
|}

Let

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>w_1 = v_2-v_1,</math>
|-
| style="text-align: center;" | <math> w_2 = v_3-v_1,</math>
|-
| style="text-align: center;" | <math> w_3 = v_3-v_2, </math>
|}
|}

and

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>\eta _1  =           (y_2 - y_3)dx +  (x_3 - x_2)dy \quad =\quad  2A\,\, dN_1  \quad      =\quad \left|\begin{array}{cc}x_3-x_2 & y_3-y_2 \\      dx & dy     \end{array} \right|, </math>
|-
| style="text-align: center;" | <math>  \eta _2 =           (y_3- y_1 )dx +  (x_1 -  x_3)dy \quad =\quad  2A\,\, dN_2  \quad      =\quad \left|\begin{array}{cc}x_1-x_3 & y_1-y_3 \\      dx & dy     \end{array} \right|, </math>
|-
| style="text-align: center;" | <math>  \eta _3 =           (y_1 -  y_2)dx +  (x_2-x_1 )dy \quad =\quad  2A\,\, dN_3  \quad      =\quad \left|\begin{array}{cc}x_2-x_1 & y_2-y_1 \\      dx & dy     \end{array} \right|. </math>
|}
|}

As it turns out,

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>(S^*\beta )_{12}     =\int _{v_1}^{v_2} S^*\beta </math>
|-
| style="text-align: center;" | <math>    = {1\over 4A} \left[(\eta _2(S(w_1))-\eta _1(S(w_1)))\beta _{12} + \eta _3(S(w_1))\beta _{13} + \eta _3(S(w_1))\beta _{23}\right],</math>
|-
| style="text-align: center;" | <math> (S^*\beta )_{13}     =\int _{v_1}^{v_3} S^*\beta </math>
|-
| style="text-align: center;" | <math>    = {1\over 4A} \left[\eta _2(S(w_2))\beta _{12} + (\eta _3(S(w_2))-\eta _1(S(w_2)))\beta _{13}  - \eta _2(S(w_2))\beta _{23}\right],</math>
|-
| style="text-align: center;" | <math> (S^*\beta )_{23}     =\int _{v_2}^{v_3} S^*\beta </math>
|-
| style="text-align: center;" | <math>    = {1\over 4A} \left[- \eta _1(S(w_3))\beta _{12} - \eta _1(S(w_3))\beta _{13} (\eta _3(S(w_3))-\eta _2(S(w_3)))\beta _{23} \right], </math>
|}
|}

where

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>\eta _1(S(w_1))    = \left|\begin{array}{cc}x_3-x_2 & y_3-y_2 \\      a(x_2-x_1)+b(y_2-y_1) & c(x_2-x_1)+d(y_2-y_1)     \end{array} \right|,</math>
|-
| style="text-align: center;" | <math> \eta _1(S(w_2))    = \left|\begin{array}{cc}x_3-x_2 & y_3-y_2 \\      a(x_3-x_1)+b(y_3-y_1) & c(x_3-x_1)+d(y_3-y_1)     \end{array} \right|,</math>
|-
| style="text-align: center;" | <math> \eta _1(S(w_3))    = \left|\begin{array}{cc}x_3-x_2 & y_3-y_2 \\      a(x_3-x_2)+b(y_3-y_2) & c(x_3-x_2)+d(y_3-y_2)     \end{array} \right|,</math>
|-
| style="text-align: center;" | <math> \eta _2(S(w_1))    = \left|\begin{array}{cc}x_1-x_3 & y_1-y_3 \\      a(x_2-x_1)+b(y_2-y_1) & c(x_2-x_1)+d(y_2-y_1)     \end{array} \right|,</math>
|-
| style="text-align: center;" | <math> \eta _2(S(w_2))    = \left|\begin{array}{cc}x_1-x_3 & y_1-y_3 \\      a(x_3-x_1)+b(y_3-y_1) & c(x_3-x_1)+d(y_3-y_1)     \end{array} \right|,</math>
|-
| style="text-align: center;" | <math> \eta _2(S(w_3))    = \left|\begin{array}{cc}x_1-x_3 & y_1-y_3 \\      a(x_3-x_2)+b(y_3-y_2) & c(x_3-x_2)+d(y_3-y_2)     \end{array} \right|,</math>
|-
| style="text-align: center;" | <math> \eta _3(S(w_1))    = \left|\begin{array}{cc}x_2-x_1 & y_2-y_1 \\      a(x_2-x_1)+b(y_2-y_1) & c(x_2-x_1)+d(y_2-y_1)     \end{array} \right|,</math>
|-
| style="text-align: center;" | <math> \eta _3(S(w_2))    = \left|\begin{array}{cc}x_2-x_1 & y_2-y_1 \\      a(x_3-x_1)+b(y_3-y_1) & c(x_3-x_1)+d(y_3-y_1)     \end{array} \right|,</math>
|-
| style="text-align: center;" | <math> \eta _3(S(w_3))    = \left|\begin{array}{cc}x_2-x_1 & y_2-y_1 \\      a(x_3-x_2)+b(y_3-y_2) & c(x_3-x_2)+d(y_3-y_2)     \end{array} \right|. </math>
|}
|}

Anyhow, these quantities are areas of the parallelograms formed by one fixed edge of the original triangle and one edge transformed by <math display="inline">S</math>.

Thus, we have that the induced transformation on primal 1-forms is

<span id="eq-6"></span>
{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math> S^*\beta ={1\over 4A}\left(\begin{array}{ccc} \eta _2(S(w_1))-\eta _1(S(w_1)) & \eta _3(S(w_1)) & \eta _3(S(w_1)) \\  \eta _2(S(w_2)) & \eta _3(S(w_2))-\eta _1(S(w_2)) & -\eta _2(S(w_2)) \\  - \eta _1(S(w_3)) & - \eta _1(S(w_3)) & \eta _3(S(w_3))-\eta _2(S(w_3)) \end{array} \right) \left(\begin{array}{c} \beta _{12} \\  \beta _{13} \\  \beta _{23} \end{array} \right)   </math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (6)
|}

When the material is isotropic, i.e. <math display="inline">S</math> is a constant multiple of the <math display="inline">2\times 2</math> identity matrix,  we get a constant multiple of the <math display="inline">3\times 3</math> identity matrix above as discrete pull-back operator.

==4. Flux and anisotropy==

In this section, we deduce the DEC formulae for the local flux, the local anisotropic flux and the local anisotropy operator  for primal 1-forms.

===4.1 The flux in local DEC===

We wish to find a natural construction for the discrete flux (discrete gradient vector) of a discrete function. Recall from Vector Calculus that the directional derivative of a differentiable function <math display="inline">f:\mathbb{R}^2\longrightarrow \mathbb{R}</math> at a point <math display="inline">p\in \mathbb{R}^2</math> in the direction of <math display="inline">w\in \mathbb{R}^2</math> is defined by

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>df_p(w):=\lim _{t\rightarrow 0}{f(p+tw)-f(p)\over t}=\nabla f(p)\cdot w.</math>
|-
<!--| style="text-align: center;" | <math> =\nabla f(p)\cdot w. </math>-->
|}
|}

Thus, we have three Vector Calculus identities

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>df_{v_1}(v_2-v_1) = \nabla f(v_1)\cdot (v_2-v_1)  ,</math>
|-
| style="text-align: center;" | <math>  df_{v_2}(v_3-v_2) = \nabla f(v_2)\cdot (v_3-v_2),</math>
|-
| style="text-align: center;" | <math>  df_{v_3}(v_1-v_3) = \nabla f(v_3)\cdot (v_1-v_3) . </math>
|}
|}

As in Subsection [[#2.2 Discrete derivative|2.2]], the rough approximations to directional derivatives of a function <math display="inline">f</math> in the directions of the (oriented) edges are given as follows

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>df_{v_1}(v_2-v_1) \approx  f_2-f_1  ,</math>
|-
| style="text-align: center;" | <math>  df_{v_2}(v_3-v_2) \approx  f_3-f_2  ,</math>
|-
| style="text-align: center;" | <math>  df_{v_3}(v_1-v_3) \approx  f_1-f_3  . </math>
|}
|}

Thus, if we want to find a discrete gradient vector <math display="inline">W_1</math> of <math display="inline">f</math> at the point <math display="inline">v_1</math>,  we need to solve the equations

<span id="eq-7"></span>
<span id="eq-8"></span>
{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>W_1\cdot (v_2-v_1) = f_2-f_1 </math>
| style="width: 5px;text-align: right;white-space: nowrap;" | (7)
|-
| style="text-align: center;" | <math> W_1\cdot (v_3-v_1) = f_3-f_1. </math>
| style="width: 5px;text-align: right;white-space: nowrap;" | (8)
|}
|}

If

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>v_1=(x_1,y_1),</math>
|-
| style="text-align: center;" | <math> v_2=(x_2,y_2),</math>
|-
| style="text-align: center;" | <math> v_3=(x_3,y_3), </math>
|}
|}

then

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>W_1 =\left({f_1y_2-f_1y_3-f_2y_1+f_2y_3+f_3y_1-f_3y_2\over x_1y_2-x_1y_3-x_2y_1+x_2y_3+x_3y_1-x_3y_2} ,  -{f_1x_2-f_1x_3-f_2x_1+f_2x_3+f_3x_1-f_3x_2\over x_1y_2-x_1y_3-x_2y_1+x_2y_3+x_3y_1-x_3y_2}\right)^T </math>
|}
|}

Now, if we were to find  a discrete gradient vector <math display="inline">W_2</math> of <math display="inline">f</math> at the point <math display="inline">v_2</math>,  we need to solve the equations

<span id="eq-9"></span>
{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>W_2\cdot (v_1-v_2) = f_1-f_2</math>
| style="width: 5px;text-align: right;white-space: nowrap;" | (9)
|-
| style="text-align: center;" | <math> W_2\cdot (v_3-v_2) = f_3-f_2.</math>
|}
|}

The vectors <math display="inline">W_2</math> solving these equations is actually equal to <math display="inline">W_1</math>.  Indeed, consider

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>f_3-f_1 = W_1\cdot (v_3-v_1) </math>
|-
| style="text-align: center;" | <math> = W_1\cdot (v_3-v_2+v_2-v_1)</math>
|-
| style="text-align: center;" | <math> = W_1\cdot (v_3-v_2)+W_1\cdot (v_2-v_1)</math>
|-
| style="text-align: center;" | <math> = W_1\cdot (v_3-v_2)+f_2-f_1, </math>
|}
|}

so that

<span id="eq-10"></span>
{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>W_1\cdot (v_3-v_2) = f_3-f_2. </math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (10)
|}

Thus, adding up ([[#eq-7|7]]) and ([[#eq-9|9]]) we get

<span id="eq-11"></span>
{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>(W_1-W_2)\cdot (v_2-v_1)=0. </math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (11)
|}

Subtracting  ([[#eq-8|8]]) from ([[#eq-10|10]]) we get

<span id="eq-12"></span>
{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>(W_1-W_2)\cdot (v_3-v_2)=0. </math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (12)
|}

Since <math display="inline">v_2-v_1</math> and <math display="inline">v_3-v_2</math> are linearly independent  and the two inner products in ([[#eq-11|11]]) and ([[#eq-12|12]]) vanish,

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>W_1-W_2=0.</math>
|}
|}

Analogously, the corresponding gradient vector <math display="inline">W_3</math> of <math display="inline">f</math> at the vertex <math display="inline">v_3</math> is equal to <math display="inline">W_1</math>. This means that the three approximate gradient vectors at the three vertices coincide, i.e. we have a unique discrete flux vector <math display="inline">W=W_1=W_2=W_3</math> on the given element.  Note that the discrete flux <math display="inline">W</math> satisfies

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>W\cdot (v_2-v_1) = f_2-f_1,</math>
|-
| style="text-align: center;" | <math> W\cdot (v_3-v_1) = f_3-f_1,</math>
|-
| style="text-align: center;" | <math> W\cdot (v_3-v_2) = f_3-f_2  . </math>
|}
|}

This means that the primal 1-form discretizing <math display="inline">df</math> can be obtained by performing the dot products of the discrete flux vector <math display="inline">W</math>  with the vectors given by the triangle's oriented edges.

====4.1.1 Comparison of DEC and FEML local fluxes====

The local flux (gradient) of <math display="inline">f</math> in FEML is given by

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math> \left( \begin{array}{ccc} \displaystyle{\partial N_1\over \partial x} & \displaystyle\displaystyle{\partial N_2\over \partial x} & \displaystyle{\partial N_3\over \partial x} \\  \displaystyle{\partial N_1\over \partial y} & \displaystyle{\partial N_2\over \partial y} & \displaystyle{\partial N_3\over \partial y} \end{array}  \right) \left( \begin{array}{c} f_1 \\  f_2 \\  f_3 \end{array}  \right), </math>
|}
|}

where <math display="inline">N_1,N_2</math> and <math display="inline">N_3</math> are the basis functions given in ([[#eq-3|3]]). Explicitly

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>\left( \begin{array}{ccc}\displaystyle{\partial N_1\over \partial x} & \displaystyle{\partial N_2\over \partial x} & \displaystyle{\partial N_3\over \partial x} \\  \displaystyle{\partial N_1\over \partial y} & \displaystyle{\partial N_2\over \partial y} & \displaystyle{\partial N_3\over \partial y} \end{array}  \right) ={1\over 2A} \left( \begin{array}{ccc}y_2-y_3 & y_3-y_1 & y_1-y_2 \\  x_3-x_2 & x_1-x_3 & x_2-x_1 \end{array}  \right)  </math>
|}
|}

where

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>A = {1\over 2} [(x_2y_3-x_3y_2) -(x_1y_3-x_3y_1)+(x_1y_2-x_2y_1)] </math>
|}
|}

is the area of the triangle, so that the FEML flux is given by

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math> \left( {[(y_2-y_3)f_1+(y_3-y_1)f_2+(y_1-y_2)f_3]\over 2A} , {[(x_3-x_2)f_1+(x_1-x_3)f_2+(x_2-x_1)f_3]\over 2A}  \right)^T, </math>
|}
|}

and we can see that it coincides with the DEC flux.

===4.2 The anisotropic flux vector in local DEC===

We will now discuss how to discretize anisotropy in 2D DEC. Let <math display="inline">K</math> denote the symmetric anisotropy tensor

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>K=\left(\begin{array}{cc} k_{11} & k_{12} \\  k_{12} & k_{22} \end{array} \right)</math>
|}
|}

and recall the anisotropic Poisson equation

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>-\nabla \cdot (K\, \nabla f) = q.</math>
|}
|}

First recall that, in Exterior Differential Calculus, for any <math display="inline">w\in \mathbb{R}^2</math>,

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>(K\nabla f(p))\cdot w =\nabla f(p)\cdot (K^Tw) =\nabla f(p)\cdot (Kw) </math>
|-
| style="text-align: center;" | <math>    =df_p(Kw)=(df_p\circ K)(w) =:(K^*df_p)(w), </math>
|}
|}

where <math display="inline">K^*df_p</math> is the ''pullback'' of <math display="inline">df_p</math> by <math display="inline">K</math>.

Since we already have a discrete candidate <math display="inline">W</math> for <math display="inline">\nabla f</math> on the given element, the product of the matrix <math display="inline">K</math> and <math display="inline">W</math> gives us a candidate <math display="inline">W':=K W</math> for the anisptropic flux vector on such an element

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math> W' := KW =\left(\begin{array}{cc}k_{11} & k_{12} \\  k_{12} & k_{22} \end{array} \right) \left( \begin{array}{c}{[(y_2-y_3)f_1+(y_3-y_1)f_2+(y_1-y_2)f_3]\over 2A} \\  {[(x_3-x_2)f_1+(x_1-x_3)f_2+(x_2-x_1)f_3]\over 2A} \end{array}  \right)</math>
|-
| style="text-align: center;" | <math>    =\left( \begin{array}{c}{k_{11}[f_1(y_2-y_3)+f_2(y_3-y_1)+f_3(y_1-y_2)] +k_{12}[f_1(x_3-x_2)+f_2(x_1-x_3)+f_3(x_2-x_1)] \over x_1y_2-x_1y_3-x_2y_1+x_2y_3+x_3y_1-x_3y_2} \\ {k_{12}[f_1(y_2-y_3)+f_2(y_3-y_1)+f_3(y_1-y_2)] +k_{22}[f_1(x_3-x_2)+f_2(x_1-x_3)+f_3(x_2-x_1)] \over x_1y_2-x_1y_3-x_2y_1+x_2y_3+x_3y_1-x_3y_2} \end{array} \right). </math>
|}
|}

In order to produce a primal 1-form from <math display="inline">W'</math>, all we have to do is assume that <math display="inline">K</math> is constant on the given triangle and take  dot products with the edges of the triangle, i.e.

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>\left(\begin{array}{c}W'\cdot w_1 \\  W'\cdot w_2 \\  W'\cdot w_3 \end{array} \right). </math>
|}
|}

Nevertheless, in order to better understand the structure of the local system (as in FEML),  it is desirable to have a factorization of the result as a matrix product

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>\left(\begin{array}{c}W'\cdot w_1 \\  W'\cdot w_2 \\  W'\cdot w_3 \end{array} \right)= K^{DEC} D_0 \left(\begin{array}{l}f_1\\ f_2\\ f_3               \end{array} \right) </math>
|}
|}

where the matrix <math display="inline">K^{DEC}</math> is a real <math display="inline">3\times 3</math> matrix depending on the geometry of the triangle and the matrix <math display="inline">K</math>, and such that it is a multiple of the identity matrix when the domain is isotropic. In order to do this, we will use the discretization of the pullback opertator of Exterior Differential Calculus for arbitrary 1-forms from Section [[#3 Pullback operator on primal 1-forms|3]]. Thus, using ([[#eq-6|6]]), <math display="inline">K^{DEC}</math> is the following <math display="inline">3\times 3</math> matrix where we have set <math display="inline">S=K</math>, i.e.

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>K^{DEC}={1\over 4A}\left(\begin{array}{ccc} \eta ^K_{21}-\eta ^K_{11} & \eta ^K_{31} & \eta ^K_{31} \\  \eta ^K_{22} & \eta ^K_{32}-\eta ^K_{12} & -\eta ^K_{22} \\  - \eta ^K_{13} & - \eta ^K_{13} & \eta ^K_{33}-\eta ^K_{23} \end{array} \right) </math>
|}
|}

where

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>\eta ^K_{kl}=\eta _k(K(w_l)).</math>
|}
|}

====4.2.1 <span id='lb-4.2.1'></span>Geometric interpretation of the entries of <math>K^{DEC}</math>====

Consider the Figure [[#img-5|5]], where <math display="inline">J</math> denotes the <math display="inline">90^\circ </math> anti-clockwise rotation

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>J    =\left(\begin{array}{cc}0 & -1 \\    1 & 0   \end{array} \right). </math>
|}
|}

<div id='img-5'></div>
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: auto;"
|-
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-DiscretizedAnisotropyOperator03-v2.png|360px|Geometric interpretation of the entries of the anisotropy tensor discretization ]]
|- style="text-align: center; font-size: 75%;"
| colspan="1" style="padding:10px;"| '''Figure 5'''. Geometric interpretation of the entries of the anisotropy tensor discretization <math>K^{DEC}</math>
|}

We have

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>\eta ^K_{21}=-{J(w_2)\cdot K(w_1) \over 2A}= -{1\over 2A} |J(w_2)| |K(w_1)| \cos (\beta ) </math>
|-
| style="text-align: center;" | <math>    = -{1\over 2A} |w_2| |K(w_1)| \cos (\alpha{+\pi}/2) = {1\over 2A} |w_2| |K(w_1)| \sin (\alpha ) ={A'\over A}, </math>
|}
|}

where <math display="inline">A'</math> is the area of the red triangle.  Thus, the numbers <math display="inline">\eta ^K_{kl}</math> are quotients of areas of transformed triangles divided by the area of the original triangle.

==5. Anisotropic Poisson equation in 2D==

In this section, we describe the local DEC discretization of the 2D anisotropic Poisson equation and compare it to that of FEML.

===5.1 Local DEC discretization of the 2D anisotropic Poisson equation===

The anisotropic Poisson equation reads as follows

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>-\nabla \cdot (K\, \nabla f) = q,</math>
|}
|}

where <math display="inline">f</math> and <math display="inline">q</math> are two functions on a certain domain in <math display="inline">\mathbb{R}^2</math> and <math display="inline">K</math> is the anisotropy tensor. In terms of the exterior derivative <math display="inline">d</math> and the Hodge star operator <math display="inline">\star </math> it reads as follows

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>-\star d \star (K^*df) = q</math>
|}
|}

where <math display="inline">K^*df:=df\circ K</math> and <math display="inline">K=K^T</math>. Following the discretization of the discretized divergence operator <span id='citeF-6'></span>[[#cite-6|[6]]], the corresponding local DEC discretization of the anisotropic Poisson equation is

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>-M_0^{-1} \, \left(-D_0^T\right) \, M_1\, K^{DEC} \, D_0 \, [f] = [q], </math>
|}
|}

or equivalently

<span id="eq-13"></span>
{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>D_0^T  \, M_1\, K^{DEC}\, D_0 \, [f] = M_0 \, [q].  </math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (13)
|}

In order to simplify the notation, consider the lengths and areas defined in the Figure [[#img-6|6]].

<div id='img-6'></div>
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: auto;"
|-
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-WellCenteredTriangle02-eps-converted-to.png|210px|Triangle ]]
|- style="text-align: center; font-size: 75%;"
| colspan="1" style="padding-bottom:10px;"| '''Figure 6'''. Triangle 
|}

Now, the discretized equation ([[#eq-13|13]]) looks as follows:

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math> {1\over 4A} \left(\begin{array}{ccc} 1 & -1 & 0\\ -1 & 0 & -1\\ 0 & 1 & 1       \end{array} \right) \left(\begin{array}{ccc} {l_1\over L_1} & 0 & 0\\ 0 & {l_3\over L_3}  & 0\\ 0 & 0 & {l_2\over L_2}        \end{array} \right) \left(\begin{array}{ccc} \eta ^K_{21}-\eta ^K_{11} & \eta ^K_{31} & \eta ^K_{31} \\  \eta ^K_{22} & \eta ^K_{32}-\eta ^K_{12} & -\eta ^K_{22} \\  - \eta ^K_{13} & - \eta ^K_{13} & \eta ^K_{33}-\eta ^K_{23} \end{array} \right)\left(\begin{array}{rrr} -1 & 1 & 0 \\  -1 & 0 & 1 \\  0 & -1 & 1 \end{array} \right) \left(\begin{array}{c} f_1 \\  f_2 \\  f_3 \end{array} \right)    = \left(\begin{array}{c} A_1q_1 \\  A_2q_2 \\  A_3q_3 \end{array} \right). </math>
|}
|}

The diffusive term matrix  is actually symmetric (see Subsection [[#5.3.1 The diffusive term|5.3.1]]).

===5.2 Local FEML-Discretized 2D anisotropic Poisson equation===

The diffusive elemental matrix in FEM (frequently called “stiffness matrix”) on an element <math display="inline">e</math> is given by

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>K_e=\int B^tDBdA,</math>
|}
|}

where <math display="inline">D=K</math> is the matrix representing the anisotropic diffusion tensor,  and the matrix <math display="inline">B</math> is given explicitly by

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math> B=\left( \begin{array}{ccc} \displaystyle{\partial N_1\over \partial x} & \displaystyle{\partial N_2\over \partial x} & \displaystyle{\partial N_3\over \partial x} \\  \displaystyle{\partial N_1\over \partial y} & \displaystyle{\partial N_2\over \partial y} & \displaystyle{\partial N_3\over \partial y} \end{array}  \right) = {1\over 2A} \left( \begin{array}{ccc} y_2-y_3 & y_3-y_1 & y_1-y_2 \\  x_3-x_2 & x_1-x_3 & x_2-x_1 \end{array}  \right).  </math>
|}
|}

Since the matrix <math display="inline">B</math> is constant on an element of the mesh, the integral is easy to compute. Thus, the difussive matrix <math display="inline">K_e</math> for a linear triangular element (FEML) is given by

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>K_e    =\int B^T DB dA =B^T D B A_e </math>
|-
| style="text-align: center;" | <math>    ={1\over 4A_e} \left( \begin{array}{ccc}y_2-y_3 & x_3-x_2  \\  y_3-y_1 & x_1-x_3  \\ y_1-y_2 & x_2-x_1 \end{array}  \right)   \left(\begin{array}{cc}k_{11} & k_{12} \\    k_{12} & k_{22}   \end{array} \right) \left( \begin{array}{ccc}y_2-y_3 & y_3-y_1 & y_1-y_2 \\  x_3-x_2 & x_1-x_3 & x_2-x_1 \end{array}  \right) </math>
|}
|}

Now, let us consider the first diagonal entry of the local FEML anisotropic Poisson diffusive matrix <math display="inline">K_e</math>,

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>(K_e)_{11}      ={1\over 4A}(k_{11}(y_2-y_3)^2 +(k_{12}+k_{12})(y_2-y_3)(x_3-x_2) + k_{22}(x_3-x_2)^2)</math>
|-
| style="text-align: center;" | <math>    ={1\over 4A}  (\begin{array}{cc}-(y_3-y_2), & x_3-x_2   \end{array} )   \left(\begin{array}{cc}k_{11} & k_{12} \\    k_{12} & k_{22}   \end{array} \right)   \left(\begin{array}{c}-(y_3-y_2) \\    x_3-x_2   \end{array} \right)</math>
|-
| style="text-align: center;" | <math>    ={1\over 4A}  (\begin{array}{cc}x_3-x_2, & y_3-y_2   \end{array} )   \left(\begin{array}{cc}0 & 1 \\    -1 & 0   \end{array} \right)   \left(\begin{array}{cc}k_{11} & k_{12} \\    k_{12} & k_{22}   \end{array} \right)   \left(\begin{array}{cc}0 & -1 \\    1 & 0   \end{array} \right)   \left(\begin{array}{c}x_3-x_2 \\    y_3-y_2   \end{array} \right)</math>
|-
| style="text-align: center;" | <math>    ={1\over 4A}  (J(v_3-v_2))^T K J(v_3-v_2)={1\over 4A}  J(v_3-v_2)\cdot K (J(v_3-v_2)), </math>
|}
|}

where <math display="inline">J</math> is the <math display="inline">90^\circ </math> anticlockwise rotation,

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>J    =\left(\begin{array}{cc}0 & -1 \\    1 & 0   \end{array} \right). </math>
|}
|}

In this notation, the diffusive term in local FEML is given as follows <span style="text-align: center; font-size: 75%;">

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math> {1\over 4A} \left(\begin{array}{ccc} J(v_3-v_2)\cdot K(J(v_3-v_2)) & J(v_3-v_2)\cdot K(J(v_1-v_3)) & J(v_3-v_2)\cdot K(J(v_2-v_1)) \\  J(v_1-v_3)\cdot K(J(v_3-v_2)) & J(v_1-v_3)\cdot K(J(v_1-v_3)) & J(v_1-v_3)\cdot K(J(v_2-v_1)) \\  J(v_2-v_1)\cdot K(J(v_3-v_2)) & J(v_2-v_1)\cdot K(J(v_1-v_3)) & J(v_2-v_1)\cdot K(J(v_2-v_1)) \end{array}\right). </math>
|}
|}

</span>

===5.3 Comparison between local DEC and FEML discretizations===

For the sake of brevity,  we are only going to compare the entries of the first row and first column of each formulation. Consider the various lengths, areas and angles labeled in Figures [[#img-6|6]] and [[#img-7|7]].

<div id='img-7'></div>
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: auto;"
|-
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-WellCenteredTriangle03-eps-converted-to.png|198px|Circumscribed triangle.]]
|- style="text-align: center; font-size: 75%;"
| colspan="1" style="padding-bottom:10px;"| '''Figure 7'''. Circumscribed triangle.
|}

We have the following identities:

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>\pi = 2(\alpha _1+\alpha _2+\alpha _3),</math>
|-
| style="text-align: center;" | <math>  {2l_i\over L_i} = \tan (\alpha _i),</math>
|-
| style="text-align: center;" | <math>  {l_i\over R} = \sin (\alpha _i),</math>
|-
| style="text-align: center;" | <math>  {L_i\over 2R} = \cos (\alpha _i),</math>
|-
| style="text-align: center;" | <math>  A_1= {L_1l_1\over 4} + {L_3l_3\over 4},</math>
|-
| style="text-align: center;" | <math>  A_2= {L_1l_1\over 4} + {L_2l_2\over 4},</math>
|-
| style="text-align: center;" | <math>  A_3= {L_2l_2\over 4} + {L_3l_3\over 4}. </math>
|}
|}

====5.3.1 The diffusive term====

For instance, we claim that

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math> {J(v_3-v_2)\cdot K(J(v_3-v_2))\over 4A} =  {1\over 4A}\left({l_1(\eta ^K_{3,1}+\eta ^K_{2, 1}-\eta ^K_{1, 1})\over L_1}+{l_3(\eta ^K_{3, 2}+\eta ^K_{2,2}-\eta ^K_{1, 2})\over L_3}\right) </math>
|}
|}

Indeed,

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>\eta ^K_{3,1}+\eta ^K_{2, 1}-\eta ^K_{1, 1}    =\left|\begin{array}{cc}x_2-x_1 & y_2-y_1 \\    dx(K(w_1)) & dy(K(w_1))   \end{array} \right|   +   \left|\begin{array}{cc}x_1-x_3 & y_1-y_3 \\    dx(K(w_1)) & dy(K(w_1))   \end{array} \right|   -   \left|\begin{array}{cc}x_3-x_2 & y_3-y_2 \\    dx(K(w_1)) & dy(K(w_1))   \end{array} \right|    </math>
|-
| style="text-align: center;" | <math>    =-2  \left|\begin{array}{cc}x_3-x_2 & y_3-y_2 \\    dx(K(w_1)) & dy(K(w_1))   \end{array} \right|   =-2 J(v_3-v_2)\cdot K(v_2-v_1)    </math>
|-
| style="text-align: center;padding-top:20px;" | <math> \eta ^K_{3, 2}+\eta ^K_{2,2}-\eta ^K_{1, 2}    =\left|\begin{array}{cc}x_2-x_1 & y_2-y_1 \\    dx(K(w_2)) & dy(K(w_2))   \end{array} \right|   +   \left|\begin{array}{cc}x_1-x_3 & y_1-y_3 \\    dx(K(w_2)) & dy(K(w_2))   \end{array} \right|   -   \left|\begin{array}{cc}x_3-x_2 & y_3-y_2 \\    dx(K(w_2)) & dy(K(w_2))   \end{array} \right|    </math>
|-
| style="text-align: center;" | <math>    =-2  \left|\begin{array}{cc}x_3-x_2 & y_3-y_2 \\    dx(K(w_2)) & dy(K(w_2))   \end{array} \right|   =-2 J(v_3-v_2)\cdot K(v_3-v_1) </math>
|}
|}

Thus,

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>\begin{array}{l} & \displaystyle{1\over 4A}\left(\displaystyle{l_1(\eta ^K_{3,1}+\eta ^K_{2, 1}-\eta ^K_{1, 1})\over L_1}+\displaystyle{l_3(\eta ^K_{3, 2}+\eta ^K_{2,2}-\eta ^K_{1, 2})\over L_3}\right)    =\\
&\qquad = -\displaystyle{1 \over 2A }\left(\displaystyle{l_1\over L_1}J(v_3-v_2)\cdot K(v_2-v_1) +  \displaystyle{l_1\over L_1}J(v_3-v_2)\cdot K(v_3-v_1) \right)    \\     
&\qquad =     \displaystyle{J(v_3-v_2)\cdot K(v_1-v_3)\over 2A}{\tan (\alpha _3)\over 2}    -\displaystyle{J(v_3-v_2)\cdot K( v_2-v_1)\over 2A}\displaystyle{\tan (\alpha _1)\over 2}\\    
&\qquad =     \displaystyle{1\over 4A}J(v_3-v_2)\cdot K((v_1-v_3)\tan (\alpha _3)- (v_2-v_1)\tan (\alpha _1)). \end{array}</math>
|}
|}

All we have to do now is show that

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>(v_1-v_3)\tan (\alpha _3)- (v_2-v_1)\tan (\alpha _1) = J(v_3-v_2).</math>
|}
|}

Note that, since <math display="inline">J(v_3-v_2)</math> is orthogonal to <math display="inline">v_3-v_2</math>, <math display="inline">J(v_3-v_2)</math> must be parallel to <math display="inline">c -  {v_2+v_3\over 2}</math>. Thus,

<span id="eq-14"></span>
{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>J(v_3-v_2) = {L_2\over l_2}\left(c -  {v_2+v_3\over 2}\right).  </math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (14)
|}

Now we are going to express <math display="inline">c</math> in terms of <math display="inline">v_1,v_2,v_3</math>. Let us consider

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>c-v_1 = a (v_2-v_1) + b(v_3-v_1)</math>
|}
|}

where <math display="inline">a,b</math> are coefficients to be determined. Taking inner products with <math display="inline">(v_2-v_1)</math> and <math display="inline">(v_3-v_1)</math> we get the two equations

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>R\cos (\alpha _1)    =aL_1+bL_3\cos (\alpha _1+\alpha _3),</math>
|-
| style="text-align: center;" | <math> R\cos (\alpha _3)    =aL_1\cos (\alpha _1+\alpha _3)+bL_3.  </math>
|}
|}

Solving for <math display="inline">a</math> and <math display="inline">b</math>

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>a    ={\sin (\alpha _3)\over 2\cos (\alpha _1)\sin (\alpha _1+\alpha _3)},</math>
|-
| style="text-align: center;" | <math> b    ={\sin (\alpha _1)\over 2\cos (\alpha _3)\sin (\alpha _1+\alpha _3)}. </math>
|}
|}

Substituting all the relevant quantities in ([[#eq-14|14]]) we have, for instance, that the coefficient of <math display="inline">(v_2-v_1)</math> is

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>2{\cos (\alpha _2)\over \sin (\alpha _2)}\left({\sin (\alpha _3)\over 2\cos (\alpha _1)\sin (\alpha _1+\alpha _3)}-{1\over 2}\right)    ={\cos (\alpha _2)\over \sin (\alpha _2)}\left({\sin (\alpha _3)-\cos (\alpha _1)\sin (\alpha _1+\alpha _3)\over \cos (\alpha _1)\sin (\alpha _1+\alpha _3)}\right)</math>
|-
| style="text-align: center;" | <math>     ={\cos (\alpha _2)\over \sin (\alpha _2)}\left({\sin (\alpha _3)-\cos (\alpha _1)(\sin (\alpha _1)\cos (\alpha _3)+\sin (\alpha _3)\cos (\alpha _1))\over \cos (\alpha _1)\sin (\pi /2-\alpha _2)}\right)</math>
|-
| style="text-align: center;" | <math>     ={\sin (\alpha _1)\over \sin (\alpha _2)}\left({\sin (\alpha _3)\sin (\alpha _1)-\cos (\alpha _1)\cos (\alpha _3)\over \cos (\alpha _1)}\right)</math>
|-
| style="text-align: center;" | <math>     =\tan (\alpha _1){-\cos (\alpha _1+\alpha _3)\over \sin (\alpha _2)}=\tan (\alpha _1){-\cos (\pi /2-\alpha _2)\over \sin (\alpha _2)}</math>
|-
| style="text-align: center;" | <math>    =\tan (\alpha _1){-\sin (\alpha _2)\over \sin (\alpha _2)}=-\tan (\alpha _1),  </math>
|}
|}

and similarly for the coefficient of <math display="inline">(v_1-v_3)</math>. The calculations for the remaining entries are similar.

Thus, the local DEC and FEML diffusive terms of the 2D anisotropic Poisson equation coincide.

====5.3.2 The source term====

As already observed in <span id='citeF-6'></span>[[#cite-6|[6]]], the right hand sides of the local DEC and FEML systems are different

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>\left(\begin{array}{c} A_1q_1\\ A_2q_2\\ A_3q_3       \end{array} \right) \not = {A\over 3}\left(\begin{array}{l} q_1\\ q_2\\ q_3         \end{array} \right).</math>
|}
|}

While FEML uses a barycentric subdivision to calculate the areas associated to each node/vertex, DEC uses a circumcentric subdivision. Eventually, this leads the DEC discretization to a better approximation of the solution (on coarse meshes).

==6. Some remarks about DEC quantities==

===6.1 The discrete Hodge star quantities revisited===

The numbers appearing in the local DEC matrices can be expressed both in terms of determinants and in terms of trigonometric functions. More precisely,

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>A_1 = {1\over 4}\left[\det \left(\begin{array}{ccc}x_1 & y_1 & 1\\ x_c & y_c & 1\\ x_2 & y_2 & 1                             \end{array} \right) +\det \left(\begin{array}{ccc}x_3 & y_3 & 1\\ x_c & y_c & 1\\ x_1 & y_1 & 1            \end{array} \right) \right]  \quad =\quad {R^2\over 4}(\sin (2\alpha _1)+\sin (2\alpha _3)),</math>
|-
| style="text-align: center;" | <math>  A_2  = {1\over 4}\left[\det \left(\begin{array}{ccc}x_1 & y_1 & 1\\ x_c & y_c & 1\\ x_2 & y_2 & 1                             \end{array} \right) +\det \left(\begin{array}{ccc}x_2 & y_2 & 1\\ x_c & y_c & 1\\ x_3 & y_3 & 1            \end{array} \right) \right] \quad = \quad {R^2\over 4}(\sin (2\alpha _1)+\sin (2\alpha _2)),</math>
|-
| style="text-align: center;" | <math>  A_3  = {1\over 4}\left[\det \left(\begin{array}{ccc}x_2 & y_2 & 1\\ x_c & y_c & 1\\ x_3 & y_3 & 1                             \end{array} \right) +\det \left(\begin{array}{ccc}x_3 & y_3 & 1\\ x_c & y_c & 1\\ x_1 & y_1 & 1            \end{array} \right) \right] \quad = \quad {R^2\over 4}(\sin (2\alpha _2)+\sin (2\alpha _3)),</math>
|-
| style="text-align: center;" | <math> {l_1\over L_1}  ={1\over L_1^2} \det \left(\begin{array}{ccc}x_1 & y_1 & 1\\ x_c & y_c & 1\\ x_2 & y_2 & 1                             \end{array} \right) \quad = \quad   {\tan (\alpha _1)\over 2},</math>
|-
| style="text-align: center;" | <math> {l_2\over L_2}  ={1\over L_2^2} \det \left(\begin{array}{ccc}x_2 & y_2 & 1\\ x_c & y_c & 1\\ x_3 & y_3 & 1                             \end{array} \right)  \quad = \quad  {\tan (\alpha _2)\over 2},</math>
|-
| style="text-align: center;" | <math> {l_3\over L_3}  ={1\over L_3^2} \det \left(\begin{array}{ccc}x_3 & y_3 & 1\\ x_c & y_c & 1\\ x_1 & y_1 & 1                             \end{array} \right)   \quad = \quad {\tan (\alpha _3)\over 2}. </math>
|}
|}

These expressions are valid regardless of the location of the circumcenter and can, indeed, take negative values. The angles that are measured in the scheme can be negative as in the obtuse triangle of Figure [[#img-8|8]].

<div id='img-8'></div>
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: auto;"
|-
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-NonWellCenteredTriangle01a-eps-converted-to.png|180px|]]
|-
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-NegativeAngles01-eps-converted-to.png|180px|Negative (exterior) angles measured in an obtuse triangle.]]
|- style="text-align: center; font-size: 75%;"
| colspan="1" style="padding:10px;"| '''Figure 8'''. Negative (exterior) angles measured in an obtuse triangle.
|}

and some quantities can even be zero. For instance, if

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>\alpha _2={\pi \over 2}-2\alpha _1,</math>
|}
|}

then

{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math>A_1=0.</math>
|}
|}

===6.2 Area weights assigned to vertices===

In order to understand how local DEC assigns area weights to vertices differently from FEML, let us consider the obtuse triangle shown in Figure [[#img-8|8]]. Let <math display="inline">p_1,p_2,p_3</math> be the middle points of the segments <math display="inline">[v_1,v_2],[v_1,v_3],[v_2,v_3]</math> respectively.     As shown in Figure [[#img-9|9]], the triangle <math display="inline">[v_1,p_3,c]</math> lies completely outside of the triangle <math display="inline">[v_1,v_2,v_3]</math>. Geometrically,  this implies that its area must be assigned a negative sign, which is confirmed  by the determinant formulas of Subsection [[#6.1 The discrete Hodge star quantities revisited|6.1]]. On the other hand, the triangle <math display="inline">[v_1,p_1,c]</math> will have positive area.  Thus, their sum gives us the area <math display="inline">A_1</math> in Figure [[#img-9|9]].

<div id='img-9'></div>
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: auto;"
|-
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-NonWellCenteredTriangle03a-eps-converted-to.png|200px|]]
|-
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-NonWellCenteredTriangle04a-eps-converted-to.png|200px|]]
|-
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-NonWellCenteredTriangle05a-eps-converted-to.png|200px|Area weight assigned to ]]
|- style="text-align: center; font-size: 75%;"
| colspan="1" style="padding-bottom:10px;"| '''Figure 9'''. Area weight assigned to <math display="inline">v_1</math>
|}


The area <math display="inline">A_3</math> is computed similarly, where the triangle <math display="inline">[p_3,v_3,c]</math> is assigned negative area (Figure [[#img-10|10]]).

<div id='img-10'></div>
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: auto;"
|-
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-NonWellCenteredTriangle06a-eps-converted-to.png|200px|]]
|-
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-NonWellCenteredTriangle07a-eps-converted-to.png|200px|]]
|-
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-NonWellCenteredTriangle08a-eps-converted-to.png|200px|]]
|- style="text-align: center; font-size: 75%;"
| colspan="1" style="padding-bottom:10px;"| '''Figure 10'''. Area weight assigned to <math display="inline">v_3</math>
|}

Note that for <math display="inline">A_2</math>, the two triangles <math display="inline">[p_1,v_2,c]</math> and <math display="inline">[v_2,p_2,c]</math> both have positive areas (Figure [[#img-11|11]]).

<div id='img-11'></div>
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: auto;"
|-
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-NonWellCenteredTriangle09a-eps-converted-to.png|200px|]]
|-
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-NonWellCenteredTriangle10a-eps-converted-to.png|200px|]]
|-
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-NonWellCenteredTriangle11a-eps-converted-to.png|200px|Area weight assigned to ]]
|- style="text-align: center; font-size: 75%;"
| colspan="1" style="padding-bottom:10px;"| '''Figure 11'''. Area weight assigned to <math display="inline">v_2</math>
|}

==7. Numerical examples==

In this section, we present three numerical examples in order to illustrate the performance of DEC resulting from the local formulation and its implementation.  In all cases, we solve the anisotropic Poisson equation. The FEML methodology that we have used in the comparison can be consulted <span id='citeF-11'></span><span id='citeF-13'></span><span id='citeF-2'></span>[[#cite-11|[2,11,13]]].

===7.1 First example: Heterogeneity===

This example is intended to highlight how Local DEC deals effectively with heterogeneous materials.  Consider the region in the plane given in Figure [[#img-12|12]]. 

<div id='img-12'></div>
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: auto;"
|-
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-Fig_12_SquareWithConditions.png|220px|Square and inner circle with different conditions.]]
|- style="text-align: center; font-size: 75%;"
| colspan="1" style="padding:10px;"| '''Figure 12'''. Square and inner circle with different conditions
|}

* The difussion constant for the region labelled mat1 is <math display="inline">k=12</math> and its source term is <math display="inline">q=20</math>.
* The difussion constant for the region labelled mat2 is <math display="inline">k=6</math> and its source term is <math display="inline">q=5</math>.


The meshes used in this example are shown in Figure [[#img-13|13]] and vary from coarse to very fine. 

<div id='img-13a'></div>
<div id='img-13b'></div>
<div id='img-13c'></div>
<div id='img-13d'></div>
<div id='img-13e'></div>
<div id='img-13f'></div>
<div id='img-13'></div>
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: auto;"
|-
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-Fig_13_Square_m1.png|150px|]]
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-Fig_13_Square_m2.png|144px|]]
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-Fig_13_Square_m3.png|144px|]]
|- style="text-align: center; font-size: 75%;"
|  (a) 
| (b) 
| (c) 
|-
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-Fig_13_Square_m4.png|150px|]]
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-Fig_13_Square_m5.png|150px|]]
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-Fig_13_Square_m6.png|150px|]]
|- style="text-align: center; font-size: 75%;"
| (d) 
| (e) 
| (f) 
|- style="text-align: center; font-size: 75%;"
| colspan="3" style="padding:10px;"| '''Figure 13'''. Six of the meshes used in the first example
|}

The numerical results for the maximum temperature value are exemplified in Table [[#table-1|1]].

<div id='table-1'></div>
<div class="center" style="font-size: 75%;">'''Table 1'''. Numerical simulation results of the first example</div>

{|  class="floating_tableSCP wikitable" style="text-align: center; margin: 1em auto;min-width:50%;font-size:85%;"
|- style="border-top: 2px solid;border-bottom: 2px solid;"
| rowspan='2' style="border-left: 2px solid;border-right: 2px solid;" | Mesh
| rowspan='2' style="border-left: 2px solid;border-right: 2px solid;" |  nodes
| rowspan='2' style="border-left: 2px solid;border-right: 2px solid;" | elements
| colspan='2' style="border-right: 2px solid;border-left: 2px solid;border-right: 2px solid;" | Max. Temp. Value
| colspan='2' style="border-right: 2px solid;border-left: 2px solid;border-right: 2px solid;" | Max. Flux Magnitude
|-style="border-top: 2px solid;border-bottom: 2px solid;"
| style="border-left: 2px solid;border-right: 1px solid;" | DEC  
| style="border-left: 1px solid;border-right: 1px solid;" | FEML
| style="border-left: 1px solid;border-right: 1px solid;" | DEC  
| style="border-left: 1px solid;border-right: 2px solid;" | FEML  
|- style="border-top: 2px solid;"
| style="border-left: 2px solid;border-right: 1px solid;" |    Figure [[#img-13|13]](a) 
| style="border-left: 1px solid;border-right: 1px solid;" | 49 
| style="border-left: 1px solid;border-right: 1px solid;" | 80 
| style="border-left: 1px solid;border-right: 1px solid;" | 5.51836 
| style="border-left: 1px solid;border-right: 1px solid;" | 5.53345 
| style="border-left: 1px solid;border-right: 1px solid;" | 13.837 
| style="border-left: 1px solid;border-right: 2px solid;" | 13.453 
|-
| style="border-left: 2px solid;border-right: 1px solid;" | Figure [[#img-13|13]](b) 
| style="border-left: 1px solid;border-right: 1px solid;" | 98 
| style="border-left: 1px solid;border-right: 1px solid;" | 162 
| style="border-left: 1px solid;border-right: 1px solid;" | 5.65826 
| style="border-left: 1px solid;border-right: 1px solid;" | 5.66648 
| style="border-left: 1px solid;border-right: 1px solid;" | 14.137 
| style="border-left: 1px solid;border-right: 2px solid;" | 14.024 
|-
| style="border-left: 2px solid;border-right: 1px solid;" | Figure [[#img-13|13]](c) 
| style="border-left: 1px solid;border-right: 1px solid;" | 258 
| style="border-left: 1px solid;border-right: 1px solid;" | 466 
| style="border-left: 1px solid;border-right: 1px solid;" | 5.70585 
| style="border-left: 1px solid;border-right: 1px solid;" | 5.71709 
| style="border-left: 1px solid;border-right: 1px solid;" | 14.858 
| style="border-left: 1px solid;border-right: 2px solid;" | 14.770 
|-
| style="border-left: 2px solid;border-right: 1px solid;" | Figure [[#img-13|13]](d) 
| style="border-left: 1px solid;border-right: 1px solid;" | 1,010 
| style="border-left: 1px solid;border-right: 1px solid;" | 1,914 
| style="border-left: 1px solid;border-right: 1px solid;" | 5.72103 
| style="border-left: 1px solid;border-right: 1px solid;" | 5.72280 
| style="border-left: 1px solid;border-right: 1px solid;" | 15.008 
| style="border-left: 1px solid;border-right: 2px solid;" | 15.006 
|-
| style="border-left: 2px solid;border-right: 1px solid;" | Figure [[#img-13|13]](e) 
| style="border-left: 1px solid;border-right: 1px solid;" | 3,813 
| style="border-left: 1px solid;border-right: 1px solid;" | 7,424 
| style="border-left: 1px solid;border-right: 1px solid;" | 5.72725 
| style="border-left: 1px solid;border-right: 1px solid;" | 5.72725 
| style="border-left: 1px solid;border-right: 1px solid;" | 15.229 
| style="border-left: 1px solid;border-right: 2px solid;" | 15.228 
|-
| style="border-left: 2px solid;border-right: 1px solid;" | Figure [[#img-13|13]](f) 
| style="border-left: 1px solid;border-right: 1px solid;" | 13,911 
| style="border-left: 1px solid;border-right: 1px solid;" | 27,420 
| style="border-left: 1px solid;border-right: 1px solid;" | 5.72821 
| style="border-left: 1px solid;border-right: 1px solid;" | 5.72826 
| style="border-left: 1px solid;border-right: 1px solid;" | 15.342 
| style="border-left: 1px solid;border-right: 2px solid;" | 15.337 
|-
| style="border-left: 2px solid;border-right: 1px solid;" | 
| style="border-left: 1px solid;border-right: 1px solid;" | 50,950 
| style="border-left: 1px solid;border-right: 1px solid;" | 101,098 
| style="border-left: 1px solid;border-right: 1px solid;" | 5.72841 
| style="border-left: 1px solid;border-right: 1px solid;" | 5.72842 
| style="border-left: 1px solid;border-right: 1px solid;" | 15.395 
| style="border-left: 1px solid;border-right: 2px solid;" | 15.396 
|-
| style="border-left: 2px solid;border-right: 1px solid;" | 
| style="border-left: 1px solid;border-right: 1px solid;" | 135,519 
| style="border-left: 1px solid;border-right: 1px solid;" | 269,700 
| style="border-left: 1px solid;border-right: 1px solid;" | 5.72845 
| style="border-left: 1px solid;border-right: 1px solid;" | 5.72845 
| style="border-left: 1px solid;border-right: 1px solid;" | 15.420 
| style="border-left: 1px solid;border-right: 2px solid;" | 15.417 
|-
| style="border-left: 2px solid;border-right: 1px solid;" | 
| style="border-left: 1px solid;border-right: 1px solid;" | 298,299 
| style="border-left: 1px solid;border-right: 1px solid;" | 594,596 
| style="border-left: 1px solid;border-right: 1px solid;" | 5.72848 
| style="border-left: 1px solid;border-right: 1px solid;" | 5.72848 
| style="border-left: 1px solid;border-right: 1px solid;" | 15.430 
| style="border-left: 1px solid;border-right: 2px solid;" | 15.429 
|-
| style="border-left: 2px solid;border-right: 1px solid;" | 
| style="border-left: 1px solid;border-right: 1px solid;" | 600,594 
| style="border-left: 1px solid;border-right: 1px solid;" | 1,198,330 
| style="border-left: 1px solid;border-right: 1px solid;" | 5.72848 
| style="border-left: 1px solid;border-right: 1px solid;" | 5.72848 
| style="border-left: 1px solid;border-right: 1px solid;" | 15.433 
| style="border-left: 1px solid;border-right: 2px solid;" | 15.433 
|- style="border-bottom: 2px solid;"
| style="border-left: 2px solid;border-right: 1px solid;" | 
| style="border-left: 1px solid;border-right: 1px solid;" | 1,175,238 
| style="border-left: 1px solid;border-right: 1px solid;" | 2,346,474 
| style="border-left: 1px solid;border-right: 1px solid;" | 5.72849 
| style="border-left: 1px solid;border-right: 1px solid;" | 5.72849 
| style="border-left: 1px solid;border-right: 1px solid;" | 15.43724 
| style="border-left: 1px solid;border-right: 2px solid;" | 15.43724 
|}


The temperature and flux-magnitude distribution fields are shown in Figure [[#img-14|14]].     

<div id='img-14a'></div>
<div id='img-14b'></div>
<div id='img-14'></div>
{| style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: auto;max-width: auto;"
|-
|  style="text-align: center;padding:10px;"| [[Image:Esqueda_et_al_2020a-Fig_14_Square_contour_temp.png|240px|Contour Fill of Temperatures]]
|  style="text-align: center;padding:10px;"| [[Image:Esqueda_et_al_2020a-Fig_14_Square_contour_flux.png|240px|Contour Fill of Flux vectors on Elems]] 
|-
|  style="text-align: center;font-size: 75%;"|(a) Contour fill of temperatures
|  style="text-align: center;font-size: 75%;"|(b) Contour fill of flux vectors on Elems
|- style="text-align: center; font-size: 75%;"
| colspan="2" style="padding:10px;"| '''Figure 14'''. Temperature and flux-magnitude distribution fields of the first example
|}


Figure [[#img-15|15]] shows the graphs of the temperature and the flux-magnitude along a horizontal line crossing the inner circle for the first two meshes. 

<div id='img-15a'></div>
<div id='img-15b'></div>
<div id='img-15c'></div>
<div id='img-15d'></div>
<div id='img-15e'></div>
<div id='img-15f'></div>
<div id='img-15'></div>
{| style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: auto;max-width: auto;"
|-
|  style="text-align: center;padding:10px;"| [[Image:Esqueda_et_al_2020a-Fig_15_Square_m1_diametral_temp.png|300px|]]
|  style="text-align: center;padding:10px;"| [[Image:Esqueda_et_al_2020a-Fig_15_Square_m1_diametral_flux.png|300px|]]
|-
|  style="text-align: center;font-size: 75%;"|(a) 
|  style="text-align: center;font-size: 75%;"|(b) 
|-
|style="text-align: center;padding:10px;"|[[Image:Esqueda_et_al_2020a-Fig_15_Square_m2_diametral_temp.png|300px|]]
|style="text-align: center;padding:10px;"|[[Image:Esqueda_et_al_2020a-Fig_15_Square_m2_diametral_flux.png|300px|]]
|- 
| style="text-align: center;font-size: 75%;"|(c) 
| style="text-align: center;font-size: 75%;"|(d) 
|-
|style="text-align: center;padding:10px;"|[[Image:Esqueda_et_al_2020a-Fig_15_Square_m3_diametral_temp.png|300px|]]
|style="text-align: center;padding:10px;"|[[Image:Esqueda_et_al_2020a-Fig_15_Square_m3_diametral_flux.png|300px|]]
|- 
| style="text-align: center;font-size: 75%;"|(e) 
| style="text-align: center;font-size: 75%;"|(f) 
|- style="text-align: center; font-size: 75%;"
| colspan="2" style="padding:10px;"| '''Figure 15'''. Temperature and Flux magnitude graphs of the first example along a cross-section of the domain for different meshes
|}


Table [[#table-2|2]] shows some global error metrics for different meshes. Figure [[#img-16|16]] shows the error evolution in <math display="inline">L^2</math> norm for this example.

<div id='table-2'></div>
<div class="center" style="font-size: 75%;">'''Table 2'''. DEC <math>L^2</math> errors in the first example</div>

{|  class="floating_tableSCP wikitable" style="text-align: center; margin: 1em auto;min-width:50%;font-size:85%;"
|- style="border-top: 2px solid;"
| style="border-left: 2px solid;border-right: 2px solid;" |    '''Mesh''' 
| style="border-left: 2px solid;border-right: 2px solid;" | '''Nodes''' 
| style="border-left: 2px solid;border-right: 2px solid;" | '''<math>    \sum (u-u_i)^2 \over nodes    </math>''' 
| style="border-left: 2px solid;border-right: 2px solid;" | '''<math>L^2</math> norm''' 
|- style="border-top: 2px solid;"
| style="border-left: 2px solid;border-right: 1px solid;" |    1 
| style="border-left: 1px solid;border-right: 1px solid;" | 49 
| style="border-left: 1px solid;border-right: 1px solid;" | 1.0537e-02 
| style="border-left: 1px solid;border-right: 2px solid;" | 1.4438e-01 
|- 
| style="border-left: 2px solid;border-right: 1px solid;" |    2 
| style="border-left: 1px solid;border-right: 1px solid;" | 98 
| style="border-left: 1px solid;border-right: 1px solid;" | 4.2447e-03 
| style="border-left: 1px solid;border-right: 2px solid;" | 4.8558e-02 
|- 
| style="border-left: 2px solid;border-right: 1px solid;" |    3 
| style="border-left: 1px solid;border-right: 1px solid;" | 258 
| style="border-left: 1px solid;border-right: 1px solid;" | 6.9781e-04 
| style="border-left: 1px solid;border-right: 2px solid;" | 3.0390e-03 
|- 
| style="border-left: 2px solid;border-right: 1px solid;" |    4 
| style="border-left: 1px solid;border-right: 1px solid;" | 1,010 
| style="border-left: 1px solid;border-right: 1px solid;" | 8.8386e-05 
| style="border-left: 1px solid;border-right: 2px solid;" | 1.4877e-04 
|- 
| style="border-left: 2px solid;border-right: 1px solid;" |    5 
| style="border-left: 1px solid;border-right: 1px solid;" | 3,813 
| style="border-left: 1px solid;border-right: 1px solid;" | 1.0736e-05 
| style="border-left: 1px solid;border-right: 2px solid;" | 7.7369e-06 
|- 
| style="border-left: 2px solid;border-right: 1px solid;" |    6 
| style="border-left: 1px solid;border-right: 1px solid;" | 13,911 
| style="border-left: 1px solid;border-right: 1px solid;" | 1.4422e-06 
| style="border-left: 1px solid;border-right: 2px solid;" | 4.9791e-07 
|- 
| style="border-left: 2px solid;border-right: 1px solid;" |    7 
| style="border-left: 1px solid;border-right: 1px solid;" | 50,950 
| style="border-left: 1px solid;border-right: 1px solid;" | 1.7582e-07 
| style="border-left: 1px solid;border-right: 2px solid;" | 2.9608e-08 
|- 
| style="border-left: 2px solid;border-right: 1px solid;" |    8 
| style="border-left: 1px solid;border-right: 1px solid;" | 135,518 
| style="border-left: 1px solid;border-right: 1px solid;" | 3.2621e-08 
| style="border-left: 1px solid;border-right: 2px solid;" | 2.9233e-09 
|- 
| style="border-left: 2px solid;border-right: 1px solid;" |    9 
| style="border-left: 1px solid;border-right: 1px solid;" | 298,299 
| style="border-left: 1px solid;border-right: 1px solid;" | 7.3566e-09 
| style="border-left: 1px solid;border-right: 2px solid;" | 3.3610e-10 
|- style="border-top: 1px solid;border-bottom: 2px solid;"
| style="border-left: 2px solid;border-right: 1px solid;" |    10 
| style="border-left: 1px solid;border-right: 1px solid;" | 603,440 
| style="border-left: 1px solid;border-right: 1px solid;" | 1.8577e-09 
| style="border-left: 1px solid;border-right: 2px solid;" | 4.9496e-11 
|}


<div id='img-16'></div>
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: auto;"
|-
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-DEC_Error_Graph_Square_Circle.png|430px|DEC ]]
|- style="text-align: center; font-size: 75%;"
| colspan="1" style="padding-bottom:10px;"| '''Figure 16'''. DEC <math>L^2</math> error in the first example
|}

===7.2 Second example: Anisotropy===

Let us solve the Poisson equation in a circle of radius one centered at the origin <math display="inline">(0,0)</math> under the following conditions (see Figure [[#img-17|17]]):

* heat anisotropic diffusion constants <math display="inline">K_x = 1.5, K_y=1.0</math>;
* material angle <math display="inline">30^\circ </math>;
* source term <math display="inline">q= 1</math>;
* Dirichlet boundary condition <math display="inline">u=10</math>.

<div id='img-18'></div>
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: auto;"
|-
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-Fig_16_CircleWithConditions.png|180px|Disk of radius one.]]
|- style="text-align: center; font-size: 75%;"
| colspan="1" style="padding-bottom:10px;"| '''Figure 17'''. Disk of radius one.
|}


The meshes used in this example are shown in Figure [[#img-18|18]] and vary from very coarse to very fine. 

<div id='img-18a'></div>
<div id='img-18b'></div>
<div id='img-18c'></div>
<div id='img-18d'></div>
<div id='img-18e'></div>
<div id='img-18f'></div>
<div id='img-18'></div>
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: auto;"
|-
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-Fig_17_CircleMesh1.png|144px|]]
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-Fig_17_CircleMesh2.png|144px|]]
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-Fig_17_CircleMesh3.png|144px|]]
|- 
| style="text-align: center;font-size: 75%;"|(a) 
| style="text-align: center;font-size: 75%;"|(b) 
| style="text-align: center;font-size: 75%;"|(c) 
|-
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-Fig_17_CircleMesh4.png|144px|]]
|colspan="1" style="padding:10px;"|[[Image:Esqueda_et_al_2020a-Fig_17_CircleMesh5.png|144px|]]
|colspan="1" style="padding:10px;"|[[Image:Esqueda_et_al_2020a-Fig_17_CircleMesh6.png|144px|]]
|- 
| style="text-align: center;font-size: 75%;"|(d) 
|colspan="1"  style="text-align: center;font-size: 75%;"|(e) 
| colspan="1" style="text-align: center;font-size: 75%;"|(f) 
|- style="text-align: center; font-size: 75%;"
| colspan="3" style="padding:10px;" | '''Figure 18'''. First six meshes used for unit disk.
|}


The numerical results for the maximum temperature value (<math display="inline">u(0,0)=10.2</math>) are exemplified in Table [[#table-3|3]] where  a comparison with the Finite Element Method with linear interpolation functions (FEML) is also shown.

<div id='table-3'></div>
<div style="text-align: center; font-size: 75%; ">'''Table 3'''. Temperature value at the point (0,0) and Flux magnitude value at the point (-1, 0) of the numerical simulations for the second example</div>  

{|  class="floating_tableSCP wikitable" style="text-align: center; margin: 1em auto;min-width:50%;font-size:85%;"
|- style="border-top: 2px solid;border-bottom: 2px solid;"
| rowspan='2' style="border-left: 2px solid;border-right: 2px solid;" | Mesh
| rowspan='2' style="border-left: 2px solid;border-right: 2px solid;" | Nodes
| rowspan='2' style="border-left: 2px solid;border-right: 2px solid;" | Elements
| colspan='2' style="border-right: 2px solid;border-left: 2px solid;border-right: 2px solid;" | Temp. Value at <math>(0,0)</math>
| colspan='2' style="border-right: 2px solid;border-left: 2px solid;border-right: 2px solid;" | Flux Magnitude at <math>(-1,0)</math>
|-style="border-top: 2px solid;border-bottom: 2px solid;"
| style="border-left: 2px solid;border-right: 1px solid;" | DEC  
| style="border-left: 1px solid;border-right: 1px solid;" | FEML
| style="border-left: 1px solid;border-right: 1px solid;" | DEC  
| style="border-left: 1px solid;border-right: 2px solid;" | FEML  
|- 
| style="border-left: 2px solid;border-right: 1px solid;" |    Figure [[#img-18|18]](a)
| style="border-left: 1px solid;border-right: 1px solid;" | 17 
| style="border-left: 1px solid;border-right: 1px solid;" | 20  
| style="border-left: 1px solid;border-right: 1px solid;" | 10.20014
| style="border-left: 1px solid;border-right: 1px solid;" | 10.19002  
| style="border-left: 1px solid;border-right: 1px solid;" | 0.42133 
| style="border-left: 1px solid;border-right: 2px solid;" | 0.43865
|-
| style="border-left: 2px solid;border-right: 1px solid;" | Figure [[#img-18|18]](b)
| style="border-left: 1px solid;border-right: 1px solid;" | 41 
| style="border-left: 1px solid;border-right: 1px solid;" | 56  
| style="border-left: 1px solid;border-right: 1px solid;" | 10.20007
| style="border-left: 1px solid;border-right: 1px solid;" | 10.19678  
| style="border-left: 1px solid;border-right: 1px solid;" | 0.48544 
| style="border-left: 1px solid;border-right: 2px solid;" | 0.49387
|-
| style="border-left: 2px solid;border-right: 1px solid;" | Figure [[#img-18|18]](c)
| style="border-left: 1px solid;border-right: 1px solid;" | 201 
| style="border-left: 1px solid;border-right: 1px solid;" | 344  
| style="border-left: 1px solid;border-right: 1px solid;" | 10.20012
| style="border-left: 1px solid;border-right: 1px solid;" | 10.20158  
| style="border-left: 1px solid;border-right: 1px solid;" | 0.52470 
| style="border-left: 1px solid;border-right: 2px solid;" | 0.52428
|-
| style="border-left: 2px solid;border-right: 1px solid;" | Figure [[#img-18|18]](d)
| style="border-left: 1px solid;border-right: 1px solid;" | 713 
| style="border-left: 1px solid;border-right: 1px solid;" | 1304 
| style="border-left: 1px solid;border-right: 1px solid;" | 10.20000
| style="border-left: 1px solid;border-right: 1px solid;" | 10.19969  
| style="border-left: 1px solid;border-right: 1px solid;" | 0.54143 
| style="border-left: 1px solid;border-right: 2px solid;" | 0.54224
|-
| style="border-left: 2px solid;border-right: 1px solid;" | Figure [[#img-18|18]](e)
| style="border-left: 1px solid;border-right: 1px solid;" | 2455
| style="border-left: 1px solid;border-right: 1px solid;" | 4660 
| style="border-left: 1px solid;border-right: 1px solid;" | 10.20000
| style="border-left: 1px solid;border-right: 1px solid;" | 10.19990  
| style="border-left: 1px solid;border-right: 1px solid;" | 0.54971 
| style="border-left: 1px solid;border-right: 2px solid;" | 0.55138
|-
| style="border-left: 2px solid;border-right: 1px solid;" | Figure [[#img-18|18]](f)
| style="border-left: 1px solid;border-right: 1px solid;" | 8180
| style="border-left: 1px solid;border-right: 1px solid;" | 15862 
| style="border-left: 1px solid;border-right: 1px solid;" | 10.20000
| style="border-left: 1px solid;border-right: 1px solid;" | 10.20002  
| style="border-left: 1px solid;border-right: 1px solid;" | 0.55326 
| style="border-left: 1px solid;border-right: 2px solid;" | 0.55409
|-
| style="border-left: 2px solid;border-right: 1px solid;" | 
| style="border-left: 1px solid;border-right: 1px solid;" | 20016 
| style="border-left: 1px solid;border-right: 1px solid;" | 39198 
| style="border-left: 1px solid;border-right: 1px solid;" | 10.20000
| style="border-left: 1px solid;border-right: 1px solid;" | 10.19999  
| style="border-left: 1px solid;border-right: 1px solid;" | 0.55470 
| style="border-left: 1px solid;border-right: 2px solid;" | 0.55520
|- style="border-bottom: 2px solid;"
| style="border-left: 2px solid;border-right: 1px solid;" | 
| style="border-left: 1px solid;border-right: 1px solid;" | 42306 
| style="border-left: 1px solid;border-right: 1px solid;" | 83362 
| style="border-left: 1px solid;border-right: 1px solid;" | 10.20000
| style="border-left: 1px solid;border-right: 1px solid;" | 10.20000  
| style="border-left: 1px solid;border-right: 1px solid;" | 0.55540 
| style="border-left: 1px solid;border-right: 2px solid;" | 0.55572
|}


The temperature distribution and Flux magnitude fields for the finest mesh are shown in Figure [[#img-19|19]]. 

<div id='img-19a'></div>
<div id='img-19b'></div>
<div id='img-19'></div>
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: auto;"
|-
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-Fig_18_CircleContourTemp.png|234px|Contour Fill of Temperatures]]
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-Fig_18_CircleContourFlux.png|234px|Contour Fill of Flux vectors on Elems]]
|- style="text-align: center; font-size: 75%;"
| (a) Contour Fill of Temperatures
| (b) Contour Fill of Flux vectors on Elems
|- style="text-align: center; font-size: 75%;"
| colspan="2" style="padding:10px;"| '''Figure 19'''. Temperature distribution and Flux magnitude fields for the finest mesh of the second example.
|}


Figures [[#img-20|20]](a), (b) and (c)  show the graphs of the temperature and flux magnitude values along a diameter of the circle for the different meshes of Figures [[#img-19|18]](a), (b) and (c), respectively. 

<div id='img-20a'></div>
<div id='img-20b'></div>
<div id='img-20c'></div>
<div id='img-20d'></div>
<div id='img-20e'></div>
<div id='img-20f'></div>
<div id='img-20'></div>
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: auto;"
|-
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-Fig_19_CircleTempCrossSection01.png|300px|]]
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-Fig_19_CircleFluxCrossSection01.png|300px|]]
|- style="text-align: center; font-size: 75%;"
| (a) 
| (b) 
|-
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-Fig_19_CircleTempCrossSection02.png|300px|]]
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-Fig_19_CircleFluxCrossSection02.png|300px|]]
|- style="text-align: center; font-size: 75%;"
| (c) 
| (d) 
|-
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-Fig_19_CircleTempCrossSection03.png|300px|]]
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-Fig_19_CircleFluxCrossSection03.png|300px|]]
|- style="text-align: center; font-size: 75%;"
| (e) 
| (f) 
|- style="text-align: center; font-size: 75%;"
| colspan="2" style="padding:10px;"| '''Figure 20'''. Temperature and Flux magnitude graphs of the second example along a diameter of the circle for different meshes: mesh in Figure Figure [[#img-18|18]](a), a-Temperature, b-Flux; mesh in Figure [[#img-18|18]](b), c-Temperature, d-Flux; mesh in Figure [[#img-18|18]](c), e-Temperature, f-Flux
|}


Table [[#table-4|4]] shows some global error metrics for different meshes.   Figure [[#img-21|21]] shows the  error evolution in <math display="inline">L^2</math> norm for this example.

<div id='table-4'></div>
<div class="center" style="font-size: 75%;">'''Table 4'''. DEC errors in the second example</div>

{|  class="floating_tableSCP wikitable" style="text-align: center; margin: 1em auto;min-width:50%;font-size:85%;"
|- style="border-top: 2px solid;"
| style="border-left: 2px solid;border-right: 2px solid;" |    '''Mesh''' 
| style="border-left: 2px solid;border-right: 2px solid;" | '''Nodes''' 
| style="border-left: 2px solid;border-right: 2px solid;" | '''<math>    \sum (u-u_i)^2 \over nodes    </math>''' 
| style="border-left: 2px solid;border-right: 2px solid;" | '''<math>L^2</math> norm''' 
|- style="border-top: 2px solid;"
| style="border-left: 2px solid;border-right: 1px solid;" |     1 
| style="border-left: 1px solid;border-right: 1px solid;" | 17 
| style="border-left: 1px solid;border-right: 1px solid;" | 1.5818e-04 
| style="border-left: 1px solid;border-right: 2px solid;" | 2.3555e-06 
|- style="border-top: 1px solid;"
| style="border-left: 2px solid;border-right: 1px solid;" |     2 
| style="border-left: 1px solid;border-right: 1px solid;" | 51 
| style="border-left: 1px solid;border-right: 1px solid;" | 2.5395e-05 
| style="border-left: 1px solid;border-right: 2px solid;" | 1.5639e-07 
|- style="border-top: 1px solid;"
| style="border-left: 2px solid;border-right: 1px solid;" |     3 
| style="border-left: 1px solid;border-right: 1px solid;" | 201 
| style="border-left: 1px solid;border-right: 1px solid;" | 3.3643e-06 
| style="border-left: 1px solid;border-right: 2px solid;" | 1.0517e-08 
|- style="border-top: 1px solid;"
| style="border-left: 2px solid;border-right: 1px solid;" |     4 
| style="border-left: 1px solid;border-right: 1px solid;" | 713 
| style="border-left: 1px solid;border-right: 1px solid;" | 5.1563e-07 
| style="border-left: 1px solid;border-right: 2px solid;" | 8.3543e-10 
|- style="border-top: 1px solid;"
| style="border-left: 2px solid;border-right: 1px solid;" |     5 
| style="border-left: 1px solid;border-right: 1px solid;" | 2,455 
| style="border-left: 1px solid;border-right: 1px solid;" | 8.9235e-08 
| style="border-left: 1px solid;border-right: 2px solid;" | 7.6073e-11 
|- style="border-top: 1px solid;"
| style="border-left: 2px solid;border-right: 1px solid;" |     6 
| style="border-left: 1px solid;border-right: 1px solid;" | 8,180 
| style="border-left: 1px solid;border-right: 1px solid;" | 3.1731e-08 
| style="border-left: 1px solid;border-right: 2px solid;" | 2.9858e-11 
|- style="border-top: 1px solid;"
| style="border-left: 2px solid;border-right: 1px solid;" |     7 
| style="border-left: 1px solid;border-right: 1px solid;" | 20,016 
| style="border-left: 1px solid;border-right: 1px solid;" | 2.0217e-08 
| style="border-left: 1px solid;border-right: 2px solid;" | 2.6580e-11 
|- style="border-top: 1px solid;"
| style="border-left: 2px solid;border-right: 1px solid;" |     8 
| style="border-left: 1px solid;border-right: 1px solid;" | 42,306 
| style="border-left: 1px solid;border-right: 1px solid;" | 1.4421e-08 
| style="border-left: 1px solid;border-right: 2px solid;" | 2.7062e-11 
|- style="border-top: 1px solid;"
| style="border-left: 2px solid;border-right: 1px solid;" |     9 
| style="border-left: 1px solid;border-right: 1px solid;" | 82,722 
| style="border-left: 1px solid;border-right: 1px solid;" | 9.8533e-09 
| style="border-left: 1px solid;border-right: 2px solid;" | 2.6164e-11 
|- style="border-top: 1px solid;"
| style="border-left: 2px solid;border-right: 1px solid;" |     10 
| style="border-left: 1px solid;border-right: 1px solid;" | 156,274 
| style="border-left: 1px solid;border-right: 1px solid;" | 7.1352e-09 
| style="border-left: 1px solid;border-right: 2px solid;" | 2.5954e-11 
|- style="border-top: 1px solid;"
| style="border-left: 2px solid;border-right: 1px solid;" |     11 
| style="border-left: 1px solid;border-right: 1px solid;" | 420,013 
| style="border-left: 1px solid;border-right: 1px solid;" | 4.4277e-09 
| style="border-left: 1px solid;border-right: 2px solid;" | 2.6151e-11 
|- style="border-top: 1px solid;border-bottom: 2px solid;"
| style="border-left: 2px solid;border-right: 1px solid;" |     12 
| style="border-left: 1px solid;border-right: 1px solid;" | 935,016 
| style="border-left: 1px solid;border-right: 1px solid;" | 2.9635e-09 
| style="border-left: 1px solid;border-right: 2px solid;" | 2.6003e-11 
|}


<div id='img-21'></div>
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: auto;"
|-
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-DEC_Error_Graph_Circle.png|430px|DEC ]]
|- style="text-align: center; font-size: 75%;"
| colspan="1" style="padding-bottom:10px;"| '''Figure 21'''. DEC <math>L^2</math> error evolution in the second example
|}

===7.3 Third example: Heterogeneity and anisotropy===

Let us solve the Poisson equation in a circle of radius on the following domain (Figure [[#img-22|22]]) with various material properties. The geometry of the domain is defined by segments of ellipses passing through the given points which have centers at the origin <math display="inline">(0,0)</math>.  

<div id='img-23'></div>
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: auto;"
|-
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-Fig_20_HuevoGeometry.png|340px|Egg-like domain with different materials.]]
|- style="text-align: center; font-size: 75%;"
| colspan="1" style="padding:10px;"| '''Figure 22'''. Egg-like domain with different materials
|}


{|  class="floating_tableSCP wikitable" style="text-align: center; margin: 1em auto;min-width:40%;font-size:85%;"
|- style="border-top: 2px solid;"
| style="border-left: 2px solid;border-right: 2px solid;" |    Point 
| style="border-left: 2px solid;border-right: 2px solid;" | <math>x</math>
| style="border-left: 2px solid;border-right: 2px solid;" | <math>y</math>
| style="border-left: 2px solid;border-right: 2px solid;" | Point 
| style="border-left: 2px solid;border-right: 2px solid;" | <math>x</math>
| style="border-left: 2px solid;border-right: 2px solid;" | <math>y</math>
|- style="border-top: 2px solid;"
| style="border-left: 2px solid;border-right: 1px solid;" |    a 
| style="border-left: 1px solid;border-right: 1px solid;" | -5 
| style="border-left: 1px solid;border-right: 1px solid;" | 0 
| style="border-left: 1px solid;border-right: 1px solid;" | A 
| style="border-left: 1px solid;border-right: 1px solid;" | 0 
| style="border-left: 1px solid;border-right: 2px solid;" | -4 
|-
| style="border-left: 2px solid;border-right: 1px solid;" | b 
| style="border-left: 1px solid;border-right: 1px solid;" | -4 
| style="border-left: 1px solid;border-right: 1px solid;" | 0 
| style="border-left: 1px solid;border-right: 1px solid;" | B 
| style="border-left: 1px solid;border-right: 1px solid;" | 0 
| style="border-left: 1px solid;border-right: 2px solid;" | -3 
|-
| style="border-left: 2px solid;border-right: 1px solid;" | c 
| style="border-left: 1px solid;border-right: 1px solid;" | -3 
| style="border-left: 1px solid;border-right: 1px solid;" | 0 
| style="border-left: 1px solid;border-right: 1px solid;" | C 
| style="border-left: 1px solid;border-right: 1px solid;" | 0 
| style="border-left: 1px solid;border-right: 2px solid;" | -2 
|-
| style="border-left: 2px solid;border-right: 1px solid;" | d 
| style="border-left: 1px solid;border-right: 1px solid;" | -1 
| style="border-left: 1px solid;border-right: 1px solid;" | 0 
| style="border-left: 1px solid;border-right: 1px solid;" | D 
| style="border-left: 1px solid;border-right: 1px solid;" | 0 
| style="border-left: 1px solid;border-right: 2px solid;" | -1 
|-
| style="border-left: 2px solid;border-right: 1px solid;" | e 
| style="border-left: 1px solid;border-right: 1px solid;" | 1 
| style="border-left: 1px solid;border-right: 1px solid;" | 0 
| style="border-left: 1px solid;border-right: 1px solid;" | E 
| style="border-left: 1px solid;border-right: 1px solid;" | 0 
| style="border-left: 1px solid;border-right: 2px solid;" | 1 
|-
| style="border-left: 2px solid;border-right: 1px solid;" | f 
| style="border-left: 1px solid;border-right: 1px solid;" | 6 
| style="border-left: 1px solid;border-right: 1px solid;" | 0 
| style="border-left: 1px solid;border-right: 1px solid;" | F 
| style="border-left: 1px solid;border-right: 1px solid;" | 0 
| style="border-left: 1px solid;border-right: 2px solid;" | 2 
|-
| style="border-left: 2px solid;border-right: 1px solid;" | g 
| style="border-left: 1px solid;border-right: 1px solid;" | 7 
| style="border-left: 1px solid;border-right: 1px solid;" | 0 
| style="border-left: 1px solid;border-right: 1px solid;" | G 
| style="border-left: 1px solid;border-right: 1px solid;" | 0 
| style="border-left: 1px solid;border-right: 2px solid;" | 3 
|- style="border-bottom: 2px solid;"
| style="border-left: 2px solid;border-right: 1px solid;" | h 
| style="border-left: 1px solid;border-right: 1px solid;" | 8 
| style="border-left: 1px solid;border-right: 1px solid;" | 0 
| style="border-left: 1px solid;border-right: 1px solid;" | H 
| style="border-left: 1px solid;border-right: 1px solid;" | 0 
| style="border-left: 1px solid;border-right: 2px solid;" | 4 
|}


* The Dirichlet boundary condition is <math display="inline">u=10</math> and material properties (anisotropic heat diffusion constants, material angles and source terms) are given according to Figure [[#img-23|23]] and the table below.  

<div id='img-24'></div>
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: auto;"
|-
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-Fig_21_HuevoWithConditions.png|300px|Dirichlet condition.]]
|- style="text-align: center; font-size: 75%;"
| colspan="1" style="padding:10px;"| '''Figure 23'''. Dirichlet condition
|}


{|  class="floating_tableSCP wikitable" style="text-align: center; margin: 1em auto;min-width:40%;font-size:85%;"
|- style="border-top: 2px solid;"
| style="text-align: left;border-left: 2px solid;border-right: 2px solid;" |     
| style="border-left: 2px solid;border-right: 2px solid;" | <math>K_x</math>
| style="border-left: 2px solid;border-right: 2px solid;" | <math>K_y</math>
| style="border-left: 2px solid;border-right: 2px solid;" | Angle
| style="border-left: 2px solid;border-right: 2px solid;" | <math>q</math>
|- style="border-top: 2px solid;"
| style="text-align: left;border-left: 2px solid;border-right: 1px solid;" |     Domain mat1
| style="border-left: 1px solid;border-right: 1px solid;" | 5 
| style="border-left: 1px solid;border-right: 1px solid;" | 25 
| style="border-left: 1px solid;border-right: 1px solid;" | 30 
| style="border-left: 1px solid;border-right: 2px solid;" | 15
|- style="border-top: 1px solid;"
| style="text-align: left;border-left: 2px solid;border-right: 1px solid;" |     Domain mat2
| style="border-left: 1px solid;border-right: 1px solid;" | 25 
| style="border-left: 1px solid;border-right: 1px solid;" | 5 
| style="border-left: 1px solid;border-right: 1px solid;" | 0 
| style="border-left: 1px solid;border-right: 2px solid;" | 5
|- style="border-top: 1px solid;"
| style="text-align: left;border-left: 2px solid;border-right: 1px solid;" |     Domain mat3
| style="border-left: 1px solid;border-right: 1px solid;" | 50 
| style="border-left: 1px solid;border-right: 1px solid;" | 12
| style="border-left: 1px solid;border-right: 1px solid;" | 45
| style="border-left: 1px solid;border-right: 2px solid;" | 5
|- style="border-top: 1px solid;border-bottom: 2px solid;"
| style="text-align: left;border-left: 2px solid;border-right: 1px solid;" |     Domain mat4
| style="border-left: 1px solid;border-right: 1px solid;" | 10 
| style="border-left: 1px solid;border-right: 1px solid;" | 35 
| style="border-left: 1px solid;border-right: 1px solid;" | 0 
| style="border-left: 1px solid;border-right: 2px solid;" | 5
|}


The meshes used in this example are shown in Figure [[#img-24|24]]. 

<div id='img-24a'></div>
<div id='img-24b'></div>
<div id='img-24c'></div>
<div id='img-24d'></div>
<div id='img-24'></div>
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: auto;"
|-
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-Fig_22_HuevoMesh1.png|210px|]]
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-Fig_22_HuevoMesh2.png|210px|]]
|- style="text-align: center; font-size: 75%;"
| (a) 
| (b) 
|-
| style="padding:10px;"|[[Image:Esqueda_et_al_2020a-Fig_22_HuevoMesh3.png|210px|]]
| style="padding:10px;"|[[Image:Esqueda_et_al_2020a-Fig_22_HuevoMesh4.png|210px|]]
|- style="text-align: center; font-size: 75%;"
| (c) 
| (d) 
|- style="text-align: center; font-size: 75%;"
| colspan="3" style="padding:10px;"| '''Figure 24'''. Meshes for layered egg-like figure
|}


The numerical results for the maximum temperature value (<math display="inline">u(0,0)=10.2</math>) are exemplified in Table [[#table-5|5]] where  a comparison with the Finite Element Method with linear interpolation functions (FEML) is also shown.

<div id='table-5'></div>
<div style="text-align: center; font-size: 75%; ">'''Table 5'''. Maximum temperature and Flux magnitude values in the numerical simulations of the third example</div>  

{|  class="floating_tableSCP wikitable" style="text-align: center; margin: 1em auto;min-width:50%;font-size:85%;"
|- style="border-top: 2px solid;"
| rowspan='2' style="border-left: 2px solid;border-right: 2px solid;" | Mesh
| rowspan='2' style="border-left: 2px solid;border-right: 2px solid;" | Nodes
| rowspan='2' style="border-left: 2px solid;border-right: 2px solid;" | Elements
| colspan='2' style="border-right: 2px solid;border-left: 2px solid;border-right: 2px solid;" | Max. Temp. Value
| colspan='2' style="border-right: 2px solid;border-left: 2px solid;border-right: 2px solid;" | Max. Flux Magnitude
|-style="border-top: 2px solid;"
| style="border-left: 2px solid;border-right: 2px solid;" | DEC  
| style="border-left: 2px solid;border-right: 2px solid;" | FEML
| style="border-left: 2px solid;border-right: 2px solid;" | DEC  
| style="border-left: 2px solid;border-right: 2px solid;" | FEML  
|- style="border-top: 2px solid;"
| style="border-left: 2px solid;border-right: 1px solid;" |    Figure [[#img-24|24]](a) 
| style="border-left: 1px solid;border-right: 1px solid;" | 342 
| style="border-left: 1px solid;border-right: 1px solid;" | 616 
| style="border-left: 1px solid;border-right: 1px solid;" | 2.79221 
| style="border-left: 1px solid;border-right: 1px solid;" | 2.79854 
| style="border-left: 1px solid;border-right: 1px solid;" | 18.41066 
| style="border-left: 1px solid;border-right: 2px solid;" | 18.40573 
|-
| style="border-left: 2px solid;border-right: 1px solid;" | Figure [[#img-24|24]](b) 
| style="border-left: 1px solid;border-right: 1px solid;" | 1,259 
| style="border-left: 1px solid;border-right: 1px solid;" | 2,384 
| style="border-left: 1px solid;border-right: 1px solid;" | 2.83929 
| style="border-left: 1px solid;border-right: 1px solid;" | 2.84727 
| style="border-left: 1px solid;border-right: 1px solid;" | 18.93838 
| style="border-left: 1px solid;border-right: 2px solid;" | 18.91532 
|-
| style="border-left: 2px solid;border-right: 1px solid;" | Figure [[#img-24|24]](c) 
| style="border-left: 1px solid;border-right: 1px solid;" | 4,467 
| style="border-left: 1px solid;border-right: 1px solid;" | 8,668 
| style="border-left: 1px solid;border-right: 1px solid;" | 2.85608 
| style="border-left: 1px solid;border-right: 1px solid;" | 2.85717 
| style="border-left: 1px solid;border-right: 1px solid;" | 19.13297 
| style="border-left: 1px solid;border-right: 2px solid;" | 19.13193 
|-
| style="border-left: 2px solid;border-right: 1px solid;" | Figure [[#img-24|24]](d)  
| style="border-left: 1px solid;border-right: 1px solid;" | 14,250  
| style="border-left: 1px solid;border-right: 1px solid;" | 28,506 
| style="border-left: 1px solid;border-right: 1px solid;" | 2.85994 
| style="border-left: 1px solid;border-right: 1px solid;" | 2.86056 
| style="border-left: 1px solid;border-right: 1px solid;" | 19.20982 
| style="border-left: 1px solid;border-right: 2px solid;" | 19.20909 
|-
| style="border-left: 2px solid;border-right: 1px solid;" | 
| style="border-left: 1px solid;border-right: 1px solid;" | 20,493 
| style="border-left: 1px solid;border-right: 1px solid;" | 40,316 
| style="border-left: 1px solid;border-right: 1px solid;" | 2.86120 
| style="border-left: 1px solid;border-right: 1px solid;" | 2.86177 
| style="border-left: 1px solid;border-right: 1px solid;" | 19.23120 
| style="border-left: 1px solid;border-right: 2px solid;" | 19.23457 
|-
| style="border-left: 2px solid;border-right: 1px solid;" | 
| style="border-left: 1px solid;border-right: 1px solid;" | 60,380 
| style="border-left: 1px solid;border-right: 1px solid;" | 119,418 
| style="border-left: 1px solid;border-right: 1px solid;" | 2.86219 
| style="border-left: 1px solid;border-right: 1px solid;" | 2.86231 
| style="border-left: 1px solid;border-right: 1px solid;" | 19.26655 
| style="border-left: 1px solid;border-right: 2px solid;" | 19.26628 
|-
| style="border-left: 2px solid;border-right: 1px solid;" | 
| style="border-left: 1px solid;border-right: 1px solid;" | 142,702 
| style="border-left: 1px solid;border-right: 1px solid;" | 283,162 
| style="border-left: 1px solid;border-right: 1px solid;" | 2.86249 
| style="border-left: 1px solid;border-right: 1px solid;" | 2.86256 
| style="border-left: 1px solid;border-right: 1px solid;" | 19.28045 
| style="border-left: 1px solid;border-right: 2px solid;" | 19.28028 
|-
| style="border-left: 2px solid;border-right: 1px solid;" | 
| style="border-left: 1px solid;border-right: 1px solid;" | 291,363 
| style="border-left: 1px solid;border-right: 1px solid;" | 579,360 
| style="border-left: 1px solid;border-right: 1px solid;" | 2.86263 
| style="border-left: 1px solid;border-right: 1px solid;" | 2.86267 
| style="border-left: 1px solid;border-right: 1px solid;" | 19.28727 
| style="border-left: 1px solid;border-right: 2px solid;" | 19.28755 
|-
| style="border-left: 2px solid;border-right: 1px solid;" | 
| style="border-left: 1px solid;border-right: 1px solid;" | 495,607 
| style="border-left: 1px solid;border-right: 1px solid;" | 986,724 
| style="border-left: 1px solid;border-right: 1px solid;" | 2.86275 
| style="border-left: 1px solid;border-right: 1px solid;" | 2.86269 
| style="border-left: 1px solid;border-right: 1px solid;" | 19.29057 
| style="border-left: 1px solid;border-right: 2px solid;" | 19.29081 
|-
| style="border-left: 2px solid;border-right: 1px solid;" | 
| style="border-left: 1px solid;border-right: 1px solid;" | 1,064,447 
| style="border-left: 1px solid;border-right: 1px solid;" | 2,122,160 
| style="border-left: 1px solid;border-right: 1px solid;" | 2.86272 
| style="border-left: 1px solid;border-right: 1px solid;" | 2.86273 
| style="border-left: 1px solid;border-right: 1px solid;" | 19.29385 
| style="border-left: 1px solid;border-right: 2px solid;" | 19.29389 
|-
| style="border-left: 2px solid;border-right: 1px solid;" | 
| style="border-left: 1px solid;border-right: 1px solid;" | 2,106,077 
| style="border-left: 1px solid;border-right: 1px solid;" | 4,202,536 
| style="border-left: 1px solid;border-right: 1px solid;" | 2.86274 
| style="border-left: 1px solid;border-right: 1px solid;" | 2.86274 
| style="border-left: 1px solid;border-right: 1px solid;" | 19.29618 
| style="border-left: 1px solid;border-right: 2px solid;" | 19.29615 
|- style="border-bottom: 2px solid;"
| style="border-left: 2px solid;border-right: 1px solid;" | 
| style="border-left: 1px solid;border-right: 1px solid;" | 4,031,557 
| style="border-left: 1px solid;border-right: 1px solid;" | 8,049,644 
| style="border-left: 1px solid;border-right: 1px solid;" | 2.86275 
| style="border-left: 1px solid;border-right: 1px solid;" | 2.86275 
| style="border-left: 1px solid;border-right: 1px solid;" | 19.29763 
| style="border-left: 1px solid;border-right: 2px solid;" | 19.29765 
|}


The temperature distribution and Flux magnitude fields for the finest mesh are shown in Figure [[#img-25|25]]. 

<div id='img-25a'></div>
<div id='img-25b'></div>
<div id='img-25'></div>
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: auto;"
|-
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-Fig_23_HuevoContourTemp.png|300px|Contour Fill of Temperatures]]
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-Fig_23_HuevoContourFlux.png|300px|Contour Fill of Flux vectors on Elems]]
|- style="text-align: center; font-size: 75%;"
| (a) Contour fill of temperatures
| (b) Contour till of flux vectors on elems
|- style="text-align: center; font-size: 75%;"
| colspan="2" style="padding:10px;"| '''Figure 25'''. Temperature distribution and Flux magnitude fields for the finest mesh of the third example
|}


Figure [[#img-26|26]] shows the graphs of the temperature and flux magnitude values along a diameter of the circle for different meshes of Figure  [[#img-24|24]]. 

<div id='img-26a'></div>
<div id='img-26b'></div>
<div id='img-26c'></div>
<div id='img-26d'></div>
<div id='img-26e'></div>
<div id='img-26f'></div>
<div id='img-26'></div>
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: auto;"
|-
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-Fig_24_HuevoTempCrossSection01.png|300px|]]
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-Fig_24_HuevoFluxCrossSection01.png|300px|]]
|- style="text-align: center; font-size: 75%;"
| (a) 
| (b) 
|-
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-Fig_24_HuevoTempCrossSection02.png|300px|]]
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-Fig_24_HuevoFluxCrossSection02.png|300px|]]
|- style="text-align: center; font-size: 75%;"
| (c) 
| (d) 
|-
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-Fig_24_HuevoTempCrossSection03.png|300px|]]
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-Fig_24_HuevoFluxCrossSection03.png|300px|]]
|- style="text-align: center; font-size: 75%;"
| (e) 
| (f) 
|- style="text-align: center; font-size: 75%;"
| colspan="2" style="padding:10px;"| '''Figure 26'''. Temperature and Flux magnitude graphs of the third example along a cross-section of the domain
for different meshes: Mesh in Figure [[#img-24|24]](a), a-Temperature, b-Flux; Mesh in Figure [[#img-24|24]](b), c-Temperature, d-Flux; Mesh in Figure [[#img-24|24]](c), e-Temperature, f-Flux
|}

Table [[#table-6|6]] shows some global error metrics for different meshes.   Figure [[#img-27|27]] shows the error evolution in <math display="inline">L^2</math> norm for this example.

<div id='table-6'></div>
<div class="center" style="font-size: 75%;">'''Table 6'''. DEC errors in the third example</div>

{|  class="floating_tableSCP wikitable" style="text-align: center; margin: 1em auto;min-width:50%;font-size:85%;"
|- style="border-top: 2px solid;"
| style="border-left: 2px solid;border-right: 2px solid;" | '''Mesh''' 
| style="border-left: 2px solid;border-right: 2px solid;" | '''Nodes''' 
| style="border-left: 2px solid;border-right: 2px solid;" | '''<math>    \sum (u-u_i)^2 \over nodes    </math>''' 
| style="border-left: 2px solid;border-right: 2px solid;" | '''<math>L^2</math> norm''' 
|- style="border-top: 2px solid;"
| style="border-left: 2px solid;border-right: 1px solid;" |     1 
| style="border-left: 1px solid;border-right: 1px solid;" | 342 
| style="border-left: 1px solid;border-right: 1px solid;" | 9.3638e-04 
| style="border-left: 1px solid;border-right: 2px solid;" | 2.7786e-02 
|- style="border-top: 1px solid;"
| style="border-left: 2px solid;border-right: 1px solid;" |     2 
| style="border-left: 1px solid;border-right: 1px solid;" | 1,259 
| style="border-left: 1px solid;border-right: 1px solid;" | 1.5173e-04 
| style="border-left: 1px solid;border-right: 2px solid;" | 3.2731e-03 
|- style="border-top: 1px solid;"
| style="border-left: 2px solid;border-right: 1px solid;" |     3 
| style="border-left: 1px solid;border-right: 1px solid;" | 4,467 
| style="border-left: 1px solid;border-right: 1px solid;" | 2.2110e-05 
| style="border-left: 1px solid;border-right: 2px solid;" | 2.2937e-04 
|- style="border-top: 1px solid;"
| style="border-left: 2px solid;border-right: 1px solid;" |     4 
| style="border-left: 1px solid;border-right: 1px solid;" | 14,250 
| style="border-left: 1px solid;border-right: 1px solid;" | 3.7233e-06 
| style="border-left: 1px solid;border-right: 2px solid;" | 1.9151e-05 
|- style="border-top: 1px solid;"
| style="border-left: 2px solid;border-right: 1px solid;" |     5 
| style="border-left: 1px solid;border-right: 1px solid;" | 20,492 
| style="border-left: 1px solid;border-right: 1px solid;" | 2.2311e-06 
| style="border-left: 1px solid;border-right: 2px solid;" | 9.5930e-06 
|- style="border-top: 1px solid;"
| style="border-left: 2px solid;border-right: 1px solid;" |     6 
| style="border-left: 1px solid;border-right: 1px solid;" | 60,380 
| style="border-left: 1px solid;border-right: 1px solid;" | 4.3769e-07 
| style="border-left: 1px solid;border-right: 2px solid;" | 8.7330e-07 
|- style="border-top: 1px solid;"
| style="border-left: 2px solid;border-right: 1px solid;" |     7 
| style="border-left: 1px solid;border-right: 1px solid;" | 142,702 
| style="border-left: 1px solid;border-right: 1px solid;" | 1.1664e-07 
| style="border-left: 1px solid;border-right: 2px solid;" | 1.3926e-07 
|- style="border-top: 1px solid;"
| style="border-left: 2px solid;border-right: 1px solid;" |     8 
| style="border-left: 1px solid;border-right: 1px solid;" | 291,369 
| style="border-left: 1px solid;border-right: 1px solid;" | 3.9764e-08 
| style="border-left: 1px solid;border-right: 2px solid;" | 3.3314e-08 
|- style="border-top: 1px solid;"
| style="border-left: 2px solid;border-right: 1px solid;" |     9 
| style="border-left: 1px solid;border-right: 1px solid;" | 497,378 
| style="border-left: 1px solid;border-right: 1px solid;" | 1.6680e-08 
| style="border-left: 1px solid;border-right: 2px solid;" | 1.0275e-08 
|- style="border-top: 1px solid;"
| style="border-left: 2px solid;border-right: 1px solid;" |     10 
| style="border-left: 1px solid;border-right: 1px solid;" | 1,067,171 
| style="border-left: 1px solid;border-right: 1px solid;" | 3.9594e-09 
| style="border-left: 1px solid;border-right: 2px solid;" | 1.2190e-09 
|- style="border-top: 1px solid;border-bottom: 2px solid;"
| style="border-left: 2px solid;border-right: 1px solid;" |     11 
| style="border-left: 1px solid;border-right: 1px solid;" | 2,106,248 
| style="border-left: 1px solid;border-right: 1px solid;" | 9.6949e-10 
| style="border-left: 1px solid;border-right: 2px solid;" | 1.4415e-10 
|}


<div id='img-27'></div>
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: auto;"
|-
|style="padding:10px;"|[[Image:Esqueda_et_al_2020a-DEC_Error_Graph_Eliptic_Egg.png|410px|DEC ]]
|- style="text-align: center; font-size: 75%;"
| colspan="1" style="padding-bottom:10px;"| '''Figure 27'''. DEC <math display="inline">L^2</math> evolution in the third example.
|}


'''Remark'''. As can be seen from the previous examples, DEC behaves well on coarse meshes.  As expected, the results of DEC and FEML are identical for fine meshes.  We would also like to point out the the computational costs of DEC and FEML are very similar  since the local matrices are very similar in profile, or even identical in some cases.

==8. Conclusions==

DEC is a relatively recent discretization scheme for PDE's which takes into account the geometric and analytic  features of the operators and the domains involved.  The main contributions of this paper are the following:

<ol>

<li>We have made explicit the local formulation of DEC, i.e. on each triangle of the mesh. As is customary, the local pieces can be assembled, which facilitates the implementation of DEC by the interested reader.   Furthermore, the profiles of the assembled DEC matrices are equal to those of assembled FEML matrices. </li>

<li> Guided by the local formulation, we have deduced a natural way to approximate the flux/gradient vector of a discretized function, which coincides with that of FEML. Such approximate flux vector automatically gives  the discretized version of the anisotropic flux. </li>

<li> We have discretized the pullback operator on continuous 1-forms using Whitney interpolation forms and have found the discrete anisotropy operator for primal 1-forms.  </li>

<li>We have deduced the local DEC formulation of the 2D anisotropic Poisson equation, and have proved that the DEC and FEML diffusion terms are identical, while the source terms are not  &#8211; due to the different area-weight allocation for the nodes. </li>

<li>Local DEC allows a simple treatment of heterogeneous material properties assigned to subdomains (element by element),  which eliminates the need of dealing with it through ad hoc modifications of the global discrete Hodge star operator matrix.  </li>
</ol>


On the other hand we would like to point the following features:

*   The area weights assigned to the nodes of the mesh when solving the 2D anisotropic Poisson equation can even be negative (when a triangle has an inner angle greater that <math display="inline">120^\circ </math>), in stark contrast to the FEML formulation.

* The computational cost of DEC is similar to that of FEML. While the numerical results of DEC and FEML on fine meshes are virtually identical, the DEC solutions are better than those of FEML on coarse meshes. Furthermore, DEC  solutions display numerical convergence,  as shown by the  error measurements.


Our future work will include the DEC discretization of convective terms, and DEC on 2-dimensional simplicial surfaces in 3D.  Preliminary results on both problems are promising and competitive with FEML.

==Acknowledgements==

The second named author was partially supported by a CONACyT grant, and would like to thank  the International Centre for Numerical Methods in Engineering (CIMNE) and the University of Swansea for their hospitality. We gratefully acknowledge the support of NVIDIA Corporation with the donation of the Titan X Pascal GPU used for this research.

<span style="text-align: center; font-size: 75%;">

==References==

<div class="auto" style="text-align: left;width: auto; margin-left: auto; margin-right: auto;font-size: 85%;">

<div id="cite-1"></div>
[[#citeF-1|[1]]]  Bossavit A. Mixed finite elements and the complex of Whitney forms. In J. Whiteman (ed.), The Mathematics of Finite Elements and Applications VI, pp. 137–144, Academic Press, 1988.  

<div id="cite-2"></div>
[[#citeF-2|[2]]]  Botello S., Moreles M.Z., Oñate E. Módulo de aplicaciones del método de los elementos finitos para resolver la ecuación de Poisson: MEFIPOISS.   Aula CIMNE-CIMAT, Septiembre 2010.  

<div id="cite-3"></div>
[[#citeF-3|[3]]]  Cartan E. Sur certaines expressions différentielles et le problème de Pfaff. Annales Scientifiques de l'École Normale Supérieure, Série 3, Paris, Gauthier-Villars, Tome 16, pp. 239-332, 1899.   

<div id="cite-4"></div>
[[#citeF-4|[4]]]  Crane K., et al. Digital geometry processing with discrete exterior calculus. ACM SIGGRAPH 2013 Courses, pp. 1-126, July 2013.  

<div id="cite-5"></div>
[[#citeF-5|[5]]]  Dassios I., et al. A mathematical model for plasticity and damage: A discrete calculus formulation. Journal of Computational and Applied Mathematics, 312:27-38, 2017.  

<div id="cite-6"></div>
[[#citeF-6|[6]]]  Esqueda H., Herrera R., Botello S., Moreles M.A. A geometric description of discrete exterior calculus for general triangulations. Rev. Int. Métodos Numér. Cálc. Diseño Ing., 35(1), 2, 2019.   https://www.scipedia.com/public/Herrera_et_al_2018b  

<div id="cite-7"></div>
[[#citeF-7|[7]]]  Griebel, M., Rieger C.,  Schier A. Upwind schemes for scalar advection-dominated problems in the Discrete Exterior Calculus. In Bothe D., Reusken A. (eds), Transport Processes at Fluidic Interfaces, Birkhäuser, Cham,  145-175, 2017.  

<div id="cite-8"></div>
[[#citeF-8|[8]]]  Hirani A.N. Discrete exterior calculus. Thesis, California Institute of Technology, 2003.  

<div id="cite-9"></div>
[[#citeF-9|[9]]]  Hirani A.N., Nakshatrala K.B.,   Chaudhry J.H. Numerical method for Darcy flow derived using discrete exterior calculus. International Journal for Computational Methods in Engineering Science and Mechanics, 16(3):151-169, 2015.  

<div id="cite-10"></div>
[[#citeF-10|[10]]]  Mohamed M.S.,  Hirani A.N.,   Samtaney R. Discrete exterior calculus discretization of incompressible Navier-Stokes equations over surface simplicial meshes. Journal of Computational Physics, 312: 175-191, 2016.  

<div id="cite-11"></div>
[[#citeF-11|[11]]]  Oñate E.  2D solids. Linear triangular and rectangular elements. In Structural Analysis with the Finite Element Method.  Linear Statics, Volume 1: Basis and Solids, Chapter 4, pp. 117-157, CIMNE-Springer, Barcelona, 2009. 

<div id="cite-12"></div>
[[#citeF-12|[12]]]  Whitney H. Geometric integration theory. Princeton University Press, 1957.  

<div id="cite-13"></div>
[[#citeF-13|[13]]]    Zienkiewicz O.C.,  Taylor R.L., Zhu J.Z. Generalization of the finite element concepts. Galerkin-weighted residual and variational approaches. In The Finite Element Method Set (Sixth Edition), Butterworth-Heinemann, Oxford, Chapter 3, pp. 54-102, 2005.  

</div>