From 3155a7ce9317a864bc4785d58804489728e322c1 Mon Sep 17 00:00:00 2001
From: Nicolas Marsic <nicolas.marsic@gmail.com>
Date: Wed, 27 Mar 2013 15:07:20 +0000
Subject: [PATCH] Fix LineEdgeBasis: went wrong way -- Add LineLagrangeBasis

---
 FunctionSpace/BasisGenerator.cpp    |  3 +-
 FunctionSpace/CMakeLists.txt        |  1 +
 FunctionSpace/LineEdgeBasis.cpp     | 18 ++++-----
 FunctionSpace/LineLagrangeBasis.cpp | 61 +++++++++++++++++++++++++++++
 FunctionSpace/LineLagrangeBasis.h   | 41 +++++++++++++++++++
 FunctionSpace/LineNedelecBasis.cpp  |  4 +-
 FunctionSpace/TriLagrangeBasis.cpp  |  4 +-
 7 files changed, 118 insertions(+), 14 deletions(-)
 create mode 100644 FunctionSpace/LineLagrangeBasis.cpp
 create mode 100644 FunctionSpace/LineLagrangeBasis.h

diff --git a/FunctionSpace/BasisGenerator.cpp b/FunctionSpace/BasisGenerator.cpp
index 7eeb73f0e4..1d8b774cb5 100644
--- a/FunctionSpace/BasisGenerator.cpp
+++ b/FunctionSpace/BasisGenerator.cpp
@@ -5,6 +5,7 @@
 #include "LineNodeBasis.h"
 #include "LineEdgeBasis.h"
 #include "LineNedelecBasis.h"
+#include "LineLagrangeBasis.h"
 
 #include "QuadNodeBasis.h"
 #include "QuadEdgeBasis.h"
@@ -68,7 +69,7 @@ BasisLocal* BasisGenerator::generateLagrange(unsigned int elementType,
 		basisType);
 
   switch(elementType){
-  case TYPE_LIN: throw Exception("Lagrange Basis on Lines not Implemented");
+  case TYPE_LIN: return new LineLagrangeBasis(order);
   case TYPE_TRI: return new TriLagrangeBasis(order);
   case TYPE_QUA: throw Exception("Lagrange Basis on Quads not Implemented");
   case TYPE_TET: return new TetLagrangeBasis(order);
diff --git a/FunctionSpace/CMakeLists.txt b/FunctionSpace/CMakeLists.txt
index 5e052f9b30..d4c99a1186 100644
--- a/FunctionSpace/CMakeLists.txt
+++ b/FunctionSpace/CMakeLists.txt
@@ -23,6 +23,7 @@ set(SRC
   LineNodeBasis.cpp
   LineEdgeBasis.cpp
   LineNedelecBasis.cpp
+  LineLagrangeBasis.cpp
 
 #  QuadNodeBasis.cpp
 #  QuadEdgeBasis.cpp
diff --git a/FunctionSpace/LineEdgeBasis.cpp b/FunctionSpace/LineEdgeBasis.cpp
index 408061a8c1..2fa03098b3 100644
--- a/FunctionSpace/LineEdgeBasis.cpp
+++ b/FunctionSpace/LineEdgeBasis.cpp
@@ -14,7 +14,7 @@ LineEdgeBasis::LineEdgeBasis(unsigned int order){
 
   // Set Basis Type //
   this->order = order;
-  
+
   type = 1;
   dim  = 1;
 
@@ -29,18 +29,18 @@ LineEdgeBasis::LineEdgeBasis(unsigned int order){
   Polynomial* intLegendre = new Polynomial[orderPlus];
 
   vector<Polynomial> first(3);
-  first[0] = Polynomial(+0.5, 0, 0, 0);
+  first[0] = Polynomial(-0.5, 0, 0, 0);
   first[1] = Polynomial( 0  , 0, 0, 0);
   first[2] = Polynomial( 0  , 0, 0, 0);
 
   vector<Polynomial> second(3);
-  second[0] = Polynomial(-0.5, 0, 0, 0);
+  second[0] = Polynomial(+0.5, 0, 0, 0);
   second[1] = Polynomial( 0  , 0, 0, 0);
   second[2] = Polynomial( 0  , 0, 0, 0);
 
   const Polynomial x[2] = {
-    Polynomial(+1, 1, 0, 0),
-    Polynomial(-1, 1, 0, 0)
+    Polynomial(-1, 1, 0, 0),
+    Polynomial(+1, 1, 0, 0)
   };
 
   // Legendre Polynomial //
@@ -52,19 +52,19 @@ LineEdgeBasis::LineEdgeBasis(unsigned int order){
   for(unsigned int s = 0; s < nRefSpace; s++)
     basis[s] = new vector<Polynomial>*[nFunction];
 
-  // Edge Based (Nedelec) // 
+  // Edge Based (Nedelec) //
   basis[0][0] = new vector<Polynomial>(first);
   basis[1][0] = new vector<Polynomial>(second);
 
   // Edge Based (High Order) //
   for(unsigned int s = 0; s < nRefSpace; s++){
     unsigned int i = 1;
-    
+
     for(unsigned int l = 1; l < orderPlus; l++){
-      basis[s][i] = 
+      basis[s][i] =
 	new vector<Polynomial>((intLegendre[l].compose
 				(x[(*(*edgeV[s])[0])[0]])).gradient());
-      
+
       i++;
     }
   }
diff --git a/FunctionSpace/LineLagrangeBasis.cpp b/FunctionSpace/LineLagrangeBasis.cpp
new file mode 100644
index 0000000000..c224d58021
--- /dev/null
+++ b/FunctionSpace/LineLagrangeBasis.cpp
@@ -0,0 +1,61 @@
+#include "Exception.h"
+#include "LineLagrangeBasis.h"
+
+LineLagrangeBasis::LineLagrangeBasis(unsigned int order){
+  // If order 0 (Nedelec): use order 1
+  if(order == 0)
+    order = 1;
+
+  // Set Basis Type //
+  this->order = order;
+
+  type = 0;
+  dim  = 1;
+
+  nVertex   = 2;
+  nEdge     = (order - 1);
+  nFace     = 0;
+  nCell     = 0;
+  nFunction = nVertex + nEdge + nFace + nCell;
+
+  // Init polynomialBasis //
+  lBasis = new polynomialBasis(getTag(order));
+
+  // Init Lagrange Point //
+  lPoint = new fullMatrix<double>(linePoint(order));
+}
+
+LineLagrangeBasis::~LineLagrangeBasis(void){
+  delete lBasis;
+  delete lPoint;
+}
+
+unsigned int LineLagrangeBasis::getTag(unsigned int order){
+  unsigned int tag = nodalBasis::getTag(TYPE_LIN, order, false);
+
+  if(tag)
+    return tag;
+
+  else
+    throw Exception
+      ("Can't instanciate an order %d Lagrangian Basis for a Line",
+       order);
+}
+
+fullMatrix<double> LineLagrangeBasis::
+linePoint(unsigned int order){
+  fullMatrix<double> line(order + 1, 1);
+  line(0 ,0) = 0;
+
+  if(order > 0){
+    line(0, 0) = -1.;
+    line(1, 0) =  1.;
+
+    double dd = 2. / order;
+
+    for(unsigned int i = 2; i < order + 1; i++)
+      line(i, 0) = -1. + dd * (i - 1);
+  }
+
+  return line;
+}
diff --git a/FunctionSpace/LineLagrangeBasis.h b/FunctionSpace/LineLagrangeBasis.h
new file mode 100644
index 0000000000..a1b7d77faa
--- /dev/null
+++ b/FunctionSpace/LineLagrangeBasis.h
@@ -0,0 +1,41 @@
+#ifndef _LINELAGRANGEBASIS_H_
+#define _LINELAGRANGEBASIS_H_
+
+#include "BasisLagrange.h"
+
+/**
+   @class LineLagrangeBasis
+   @brief Lagrange Basis for Lines
+
+   This class can instantiate a @em Lagrange @em Basis
+   for a Line and for a given Order.@n
+
+   It uses
+   <a href="http://geuz.org/gmsh/">gmsh</a> Basis.
+ */
+
+class LineLagrangeBasis: public BasisLagrange{
+ public:
+  //! @param order A natural number
+  //!
+  //! Returns a new LineLagrangeBasis
+  //! of the given Order
+  LineLagrangeBasis(unsigned int order);
+
+  //! Deletes this Basis
+  //!
+  virtual ~LineLagrangeBasis(void);
+
+ private:
+  //! @param order A natural number
+  //! @return Returns the @em tag of a @em Line of
+  //! the given order
+  static unsigned int getTag(unsigned int order);
+
+  //! @param order A natural number
+  //! @return Returns Lagrangian Points on a Line
+  //! for the given Order
+  static fullMatrix<double> linePoint(unsigned int order);
+};
+
+#endif
diff --git a/FunctionSpace/LineNedelecBasis.cpp b/FunctionSpace/LineNedelecBasis.cpp
index 1cab63473b..3a9e9cff20 100644
--- a/FunctionSpace/LineNedelecBasis.cpp
+++ b/FunctionSpace/LineNedelecBasis.cpp
@@ -23,12 +23,12 @@ LineNedelecBasis::LineNedelecBasis(void){
 
   // Alloc Temporary Space //
   vector<Polynomial> first(3);
-  first[0] = Polynomial(+0.5, 0, 0, 0);
+  first[0] = Polynomial(-0.5, 0, 0, 0);
   first[1] = Polynomial( 0  , 0, 0, 0);
   first[2] = Polynomial( 0  , 0, 0, 0);
 
   vector<Polynomial> second(3);
-  second[0] = Polynomial(-0.5, 0, 0, 0);
+  second[0] = Polynomial(+0.5, 0, 0, 0);
   second[1] = Polynomial( 0  , 0, 0, 0);
   second[2] = Polynomial( 0  , 0, 0, 0);
 
diff --git a/FunctionSpace/TriLagrangeBasis.cpp b/FunctionSpace/TriLagrangeBasis.cpp
index 4d4c03925c..0812ff9b9d 100644
--- a/FunctionSpace/TriLagrangeBasis.cpp
+++ b/FunctionSpace/TriLagrangeBasis.cpp
@@ -14,8 +14,8 @@ TriLagrangeBasis::TriLagrangeBasis(unsigned int order){
 
   nVertex   = 3;
   nEdge     = 3 * (order - 1);
-  nFace     = 0;
-  nCell     =     (order - 1) * (order - 2) / 2;
+  nFace     =     (order - 1) * (order - 2) / 2;
+  nCell     = 0;
   nFunction = nVertex + nEdge + nFace + nCell;
 
   // Init polynomialBasis //
-- 
GitLab