Peano 4
Loading...
Searching...
No Matches
DiscontinuousGalerkinDiscretisation.template.cpp
Go to the documentation of this file.
1#include "{{CLASSNAME}}.h"
2#include "repositories/SolverRepository.h" //bad template practice; quick fix
3constexpr double pi = 3.14159265358979323846264338328;
4
5double f(
7){
8 double f = -8 * pi * pi;
9 for (int i=0; i<Dimensions; i++){
10 f *= std::sin( 2 * pi * x[i] );
11 }
12 return f;
13}
14
15{{NAMESPACE | join("::")}}::{{CLASSNAME}}::{{CLASSNAME}}() {
16 // @todo add your stuff here
17}
18
19
20{{NAMESPACE | join("::")}}::{{CLASSNAME}}::~{{CLASSNAME}}() {
21 // @todo add your stuff here
22}
23
24
25void {{NAMESPACE | join("::")}}::{{CLASSNAME}}::vertexInit(
28 {% if VERTEX_CARDINALITY==1 %}
29 double& value,
30 double& rhs
31 {% else %}
32 tarch::la::Vector< {{VERTEX_CARDINALITY}}, double >& value,
33 tarch::la::Vector< {{VERTEX_CARDINALITY}}, double >& rhs
34 {% endif %}
35) {
36 // @todo add your stuff here
37}
38
39/*
40
41todo: instead of harcoding
42 tarch::la::Vector<Dimensions, double> tempX = {x[0] - h[0]/2, x[1] - h[1]/2};
43make it:
44 tarch::la::Vector<Dimensions, double> tempX = {x[0] - h[0]/CARDINALITY, x[1] - h[1]/CARDINALITY};
45and do the same for the nested for loop below
46
47*/
48
49void {{NAMESPACE | join("::")}}::{{CLASSNAME}}::cellInit(
52 // for now, just use value, rhs. not sure what this will be in
53 // future.
54 tarch::la::Vector< {{CELL_CARDINALITY}}, double >& value,
55 tarch::la::Vector< {{CELL_CARDINALITY}}, double >& rhs
56) {
57 double cellSize = 1;
58 for (int i=0; i<Dimensions; i++)
59 cellSize *= h[i];
60
61 //set all entries to 0 by default. might not be necessary
62 for (int i=0; i<{{CELL_CARDINALITY}}; i++)
63 {
64 value[i] = 0;
65 rhs[i] = 0;
66 }
67
68 const int nPlusOne = {{POLYNOMIAL_DEGREE}}+1;
69 //proper loop for rhs construction
70 //we loop over point i in the mesh:
71 for (int i=0; i<{{CELL_DOFS}}; i++)
72 {
73 //we take the integral over the whole cell
74 double integral=0;
75 int p=0; int w=0;
76 for (; p<repositories::quadWeights.size(); p++, w++){
77 const std::vector<double> localQuadPoint = {{NAMESPACE | join("::")}}::repositories::quadPoints[p];
78 const double localQuadWeight = {{NAMESPACE | join("::")}}::repositories::quadWeights[w];
79 const tarch::la::Vector<Dimensions,double> quadPointInCell = {{NAMESPACE | join("::")}}::repositories::convertQuadPointToRealLocation(
80 p, x, h
81 );
82
83 integral += 0.25 * cellSize * localQuadWeight * repositories::evalPoly(i, localQuadPoint) * f(quadPointInCell);
84 }
85 /*
86 this i runs over the number of dofs in this cell. ie we don't have the
87 required multicplicity here.
88 we want to set all the u values with this integral. at time of writing
89 we have no policy for setting phi_x and phi_y. so we multiply by
90 number of dofs per cell to "push" this further.
91 */
92 rhs[i * {{CELL_DOF_UNKNOWNS}}] = integral;
93 }
94
95}
96
97void {{NAMESPACE | join("::")}}::{{CLASSNAME}}::faceInit(
100 // for now, just use value, rhs. not sure what this will be in
101 // future.
102 tarch::la::Vector< {{FACE_CARDINALITY}}, double >& value,
103 tarch::la::Vector< {{FACE_CARDINALITY}}, double >& rhs
104) {}
double f(const tarch::la::Vector< Dimensions, double > &x)
tuple w
Definition ModeCalc.py:195
int i
Definition makeIC.py:51
Simple vector class.
Definition Vector.h:134