Skip to content
Snippets Groups Projects
Commit e66e9e63 authored by Nicolas Marsic's avatar Nicolas Marsic
Browse files

LineEdgeBasis: Convergence Test -- Seems OK

parent 410245bd
No related branches found
No related tags found
No related merge requests found
......@@ -3,6 +3,8 @@
#include "Exception.h"
#include "LineNodeBasis.h"
#include "LineEdgeBasis.h"
#include "LineNedelecBasis.h"
#include "QuadNodeBasis.h"
#include "QuadEdgeBasis.h"
......@@ -43,7 +45,10 @@ Basis* BasisGenerator::linGen(int basisType,
int order){
switch(basisType){
case 0: return new LineNodeBasis(order);
case 1: throw Exception("1-form not implemented on Lines");
case 1:
if (order == 0) return new LineNedelecBasis();
else return new LineEdgeBasis(order);
case 2: throw Exception("2-form not implemented on Lines");
case 3: throw Exception("3-form not implemented on Lines");
......
......@@ -21,6 +21,8 @@ set(SRC
TriLagrangeBasis.cpp
LineNodeBasis.cpp
LineEdgeBasis.cpp
LineNedelecBasis.cpp
QuadNodeBasis.cpp
QuadEdgeBasis.cpp
......
#include "LineEdgeBasis.h"
#include "Legendre.h"
using namespace std;
LineEdgeBasis::LineEdgeBasis(int order){
// Set Basis Type //
this->order = order;
type = 1;
dim = 1;
nVertex = 0;
nEdge = (order + 1);
nFace = 0;
nCell = 0;
nEdgeClosure = 2;
nFaceClosure = 0;
size = nVertex + nEdge + nFace + nCell;
// Alloc Temporary Space //
const int orderPlus = order + 1;
Polynomial* intLegendre = new Polynomial[orderPlus];
const Polynomial plusOneHalf(+0.5, 0, 0, 0);
const Polynomial minusOneHalf(-0.5, 0, 0, 0);
const Polynomial zero(0, 0, 0, 0);
const Polynomial x[2] = {
Polynomial(+1, 1, 0, 0),
Polynomial(-1, 1, 0, 0)
};
// Legendre Polynomial //
Legendre::integrated(intLegendre, orderPlus);
// Permutations //
const int permutation[2] = {0, 1};
// Basis //
node = new vector<vector<Polynomial>*>(nVertex);
edge = new vector<vector<vector<Polynomial>*>*>(nEdgeClosure);
face = new vector<vector<vector<Polynomial>*>*>(nFaceClosure);
cell = new vector<vector<Polynomial>*>(nCell);
(*edge)[0] = new vector<vector<Polynomial>*>(nEdge);
(*edge)[1] = new vector<vector<Polynomial>*>(nEdge);
// Edge Based (Nedelec) //
(*(*edge)[0])[0] = new vector<Polynomial>(3);
(*(*edge)[0])[0]->at(0) = plusOneHalf;
(*(*edge)[0])[0]->at(1) = zero;
(*(*edge)[0])[0]->at(2) = zero;
(*(*edge)[1])[0] = new vector<Polynomial>(3);
(*(*edge)[1])[0]->at(0) = minusOneHalf;
(*(*edge)[1])[0]->at(1) = zero;
(*(*edge)[1])[0]->at(2) = zero;
// Edge Based (High Order) //
for(int c = 0; c < 2; c++){
unsigned int i = 0;
for(int l = 1; l < orderPlus; l++){
(*(*edge)[c])[i + 1] =
new vector<Polynomial>((intLegendre[l].compose(x[permutation[c]])).gradient());
i++;
}
}
// Free Temporary Space //
delete[] intLegendre;
}
LineEdgeBasis::~LineEdgeBasis(void){
// Vertex Based //
for(int i = 0; i < nVertex; i++)
delete (*node)[i];
delete node;
// Edge Based //
for(int c = 0; c < 2; c++){
for(int i = 0; i < nEdge; i++)
delete (*(*edge)[c])[i];
delete (*edge)[c];
}
delete edge;
// Face Based //
delete face;
// Cell Based //
delete cell;
}
#ifndef _LINEEDGEBASIS_H_
#define _LINEEDGEBASIS_H_
#include "BasisVector.h"
/**
@class LineEdgeBasis
@brief An Edge Basis for Lines
This class can instantiate an Edge-Based Basis
(high or low order) for Lines.@n
It uses an @em adaptation of
<a href="http://www.hpfem.jku.at/publications/szthesis.pdf">Zaglmayr's</a>
Basis for @em high @em order Polynomial%s generation.@n
This Basis is a restriction of a Quad Basis to @f$y = 0@f$.@n
It also uses the following mapping: @f$x = \frac{u + 1}{2}@f$.
*/
class LineEdgeBasis: public BasisVector{
public:
//! @param order The order of the Basis
//!
//! Returns a new Edge-Basis for Lines of the given order
LineEdgeBasis(int order);
//! Deletes this Basis
//!
virtual ~LineEdgeBasis(void);
};
#endif
#include "LineNedelecBasis.h"
#include "Legendre.h"
using namespace std;
LineNedelecBasis::LineNedelecBasis(void){
// Set Basis Type //
order = 1;
type = 1;
dim = 1;
nVertex = 0;
nEdge = 1;
nFace = 0;
nCell = 0;
nEdgeClosure = 2;
nFaceClosure = 0;
size = 1;
// Alloc Temporary Space //
const Polynomial plusOneHalf(+0.5, 0, 0, 0);
const Polynomial minusOneHalf(-0.5, 0, 0, 0);
const Polynomial zero(0, 0, 0, 0);
// Basis //
node = new vector<vector<Polynomial>*>(nVertex);
edge = new vector<vector<vector<Polynomial>*>*>(nEdgeClosure);
face = new vector<vector<vector<Polynomial>*>*>(nFaceClosure);
cell = new vector<vector<Polynomial>*>(nCell);
(*edge)[0] = new vector<vector<Polynomial>*>(nEdge);
(*edge)[1] = new vector<vector<Polynomial>*>(nEdge);
// Nedelec //
(*(*edge)[0])[0] = new vector<Polynomial>(3);
(*(*edge)[0])[0]->at(0) = plusOneHalf;
(*(*edge)[0])[0]->at(1) = zero;
(*(*edge)[0])[0]->at(2) = zero;
(*(*edge)[1])[0] = new vector<Polynomial>(3);
(*(*edge)[1])[0]->at(0) = minusOneHalf;
(*(*edge)[1])[0]->at(1) = zero;
(*(*edge)[1])[0]->at(2) = zero;
}
LineNedelecBasis::~LineNedelecBasis(void){
// Vertex Based //
for(int i = 0; i < nVertex; i++)
delete (*node)[i];
delete node;
// Edge Based //
for(int c = 0; c < 2; c++){
for(int i = 0; i < nEdge; i++)
delete (*(*edge)[c])[i];
delete (*edge)[c];
}
delete edge;
// Face Based //
delete face;
// Cell Based //
delete cell;
}
#ifndef _LINENEDELECBASIS_H_
#define _LINENEDELECBASIS_H_
#include "BasisVector.h"
/**
@class LineNedelecBasis
@brief Nedelec Basis for Lines
This class can instantiate a Nedelec Basis
for Lines.@n
*/
class LineNedelecBasis: public BasisVector{
public:
//! Returns a new Nedelec Basis for Lines
//!
LineNedelecBasis(void);
//! Deletes this Basis
//!
virtual ~LineNedelecBasis(void);
};
#endif
......@@ -20,10 +20,20 @@ LineNodeBasis::LineNodeBasis(int order){
size = nVertex + nEdge + nFace + nCell;
// Legendre Polynomial //
// Alloc Temporary Space //
Polynomial* intLegendre = new Polynomial[order];
const Polynomial x[2] = {
Polynomial(+1, 1, 0, 0),
Polynomial(-1, 1, 0, 0)
};
// Legendre Polynomial //
Legendre::integrated(intLegendre, order);
// Permutations //
const int permutation[2] = {0, 1};
// Basis //
node = new vector<Polynomial*>(nVertex);
edge = new vector<vector<Polynomial*>*>(nEdgeClosure);
......@@ -45,12 +55,6 @@ LineNodeBasis::LineNodeBasis(int order){
// Edge Based //
const int permutation[2] = {0, 1};
const Polynomial x[2] = {
Polynomial(+1, 1, 0, 0),
Polynomial(-1, 1, 0, 0)
};
for(int c = 0; c < 2; c++){
unsigned int i = 0;
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment