Peano 4
Loading...
Searching...
No Matches
ApplySplit1DRiemannToPatchTest.cpp
Go to the documentation of this file.
1// This file is part of the ExaHyPE2 project. For conditions of distribution and
2// use, please see the copyright notice at www.peano-framework.org
4
6#include "tarch/la/Vector.h"
9
11 TestCase("exahype2::fv::rusanov::tests::ApplySplit1DRiemannToPatchTest") {}
12
14
16 double* numericalFluxL = ::tarch::allocateMemory<double>(5, tarch::MemoryLocation::Heap);
17 double* numericalFluxR = ::tarch::allocateMemory<double>(5, tarch::MemoryLocation::Heap);
18 numericalFluxL[0] = 1.0;
19 numericalFluxL[1] = 2.0;
20 numericalFluxL[2] = 3.0;
21 numericalFluxL[3] = 4.0;
22 numericalFluxL[4] = 5.0;
23 numericalFluxR[0] = 6.0;
24 numericalFluxR[1] = 7.0;
25 numericalFluxR[2] = 8.0;
26 numericalFluxR[3] = 9.0;
27 numericalFluxR[4] = 10.0;
28
29 constexpr int NumberOfVolumesPerAxisInPatch = 10;
30 constexpr int NumberOfUnknowns = 5;
31 constexpr int NumberOfAuxiliaryVariables = 0;
32 constexpr double dt = 0.1;
33 const tarch::la::Vector<Dimensions, double> volumeH = {0.01, 0.01};
34
35 double* Qout = new double[500];
36 double* Qcmp = new double[500];
37 for (unsigned int i = 0; i < 500; ++i) {
38 Qout[i] = 0.0;
39 Qcmp[i] = 1.0 * i;
40 }
41
42 for (int y = 0; y < NumberOfVolumesPerAxisInPatch; y++) {
43 for (int x = 0; x <= NumberOfVolumesPerAxisInPatch; x++) {
44 const int leftVoxelInImage = x - 1 + y * NumberOfVolumesPerAxisInPatch;
45 const int rightVoxelInImage = x + y * NumberOfVolumesPerAxisInPatch;
46 for (int unknown = 0; unknown < NumberOfUnknowns; unknown++) {
47 if (x > 0) {
48 Qout[leftVoxelInImage * (NumberOfUnknowns + NumberOfAuxiliaryVariables) + unknown] -= dt / volumeH(0) * numericalFluxL[unknown];
49 Qout[leftVoxelInImage * (NumberOfUnknowns + NumberOfAuxiliaryVariables) + unknown] += Qcmp[leftVoxelInImage * (NumberOfUnknowns + NumberOfAuxiliaryVariables) + unknown];
50 }
51 if (x < NumberOfVolumesPerAxisInPatch) {
52 Qout[rightVoxelInImage * (NumberOfUnknowns + NumberOfAuxiliaryVariables) + unknown] += dt / volumeH(0) * numericalFluxR[unknown];
53 Qout[rightVoxelInImage * (NumberOfUnknowns + NumberOfAuxiliaryVariables) + unknown] += Qcmp
54 [rightVoxelInImage * (NumberOfUnknowns + NumberOfAuxiliaryVariables) + unknown];
55 }
56 }
57 }
58 }
59
60 for (int y = 0; y <= NumberOfVolumesPerAxisInPatch; y++) {
61 for (int x = 0; x < NumberOfVolumesPerAxisInPatch; x++) {
62 const int lowerVoxelInImage = x + (y - 1) * NumberOfVolumesPerAxisInPatch;
63 const int upperVoxelInImage = x + y * NumberOfVolumesPerAxisInPatch;
64 for (int unknown = 0; unknown < NumberOfUnknowns; unknown++) {
65 if (y > 0) {
66 Qout[lowerVoxelInImage * (NumberOfUnknowns + NumberOfAuxiliaryVariables) + unknown] -= dt / volumeH(0) * numericalFluxL[unknown];
67 Qout[lowerVoxelInImage * (NumberOfUnknowns + NumberOfAuxiliaryVariables) + unknown] += Qcmp
68 [lowerVoxelInImage * (NumberOfUnknowns + NumberOfAuxiliaryVariables) + unknown];
69 }
70 if (y < NumberOfVolumesPerAxisInPatch) {
71 Qout[upperVoxelInImage * (NumberOfUnknowns + NumberOfAuxiliaryVariables) + unknown] += dt / volumeH(0) * numericalFluxR[unknown];
72 Qout[upperVoxelInImage * (NumberOfUnknowns + NumberOfAuxiliaryVariables) + unknown] += Qcmp
73 [upperVoxelInImage * (NumberOfUnknowns + NumberOfAuxiliaryVariables) + unknown];
74 }
75 }
76 }
77 }
78
79 double* Qouttest = new double[500];
80 for (unsigned int i = 0; i < 500; ++i) {
81 Qouttest[i] = 0.0;
82 }
83
84 for (int shift = 0; shift < 2; shift++) {
85 for (int x = shift; x <= NumberOfVolumesPerAxisInPatch; x += 2) {
86 for (int y = 0; y < NumberOfVolumesPerAxisInPatch; y++) {
87 const int leftVoxelInImage = x - 1 + y * NumberOfVolumesPerAxisInPatch;
88 const int rightVoxelInImage = x + y * NumberOfVolumesPerAxisInPatch;
89 for (int unknown = 0; unknown < NumberOfUnknowns; unknown++) {
90 if (x > 0) {
91 Qouttest[leftVoxelInImage * (NumberOfUnknowns + NumberOfAuxiliaryVariables) + unknown] -= dt / volumeH(0) * numericalFluxL[unknown];
92 Qouttest[leftVoxelInImage * (NumberOfUnknowns + NumberOfAuxiliaryVariables) + unknown] += Qcmp
93 [leftVoxelInImage * (NumberOfUnknowns + NumberOfAuxiliaryVariables) + unknown];
94 }
95 if (x < NumberOfVolumesPerAxisInPatch) {
96 Qouttest[rightVoxelInImage * (NumberOfUnknowns + NumberOfAuxiliaryVariables) + unknown] += dt / volumeH(0) * numericalFluxR[unknown];
97 Qouttest[rightVoxelInImage * (NumberOfUnknowns + NumberOfAuxiliaryVariables) + unknown] += Qcmp
98 [rightVoxelInImage * (NumberOfUnknowns + NumberOfAuxiliaryVariables) + unknown];
99 }
100 }
101 }
102 }
103 }
104
105 for (int shift = 0; shift < 2; shift++) {
106 for (int y = shift; y <= NumberOfVolumesPerAxisInPatch; y += 2) {
107 for (int x = 0; x < NumberOfVolumesPerAxisInPatch; x++) {
108 const int lowerVoxelInImage = x + (y - 1) * NumberOfVolumesPerAxisInPatch;
109 const int upperVoxelInImage = x + y * NumberOfVolumesPerAxisInPatch;
110 for (int unknown = 0; unknown < NumberOfUnknowns; unknown++) {
111 if (y > 0) {
112 Qouttest[lowerVoxelInImage * (NumberOfUnknowns + NumberOfAuxiliaryVariables) + unknown] -= dt / volumeH(0) * numericalFluxL[unknown];
113 Qouttest[lowerVoxelInImage * (NumberOfUnknowns + NumberOfAuxiliaryVariables) + unknown] += Qcmp
114 [lowerVoxelInImage * (NumberOfUnknowns + NumberOfAuxiliaryVariables) + unknown];
115 }
116 if (y < NumberOfVolumesPerAxisInPatch) {
117 Qouttest[upperVoxelInImage * (NumberOfUnknowns + NumberOfAuxiliaryVariables) + unknown] += dt / volumeH(0) * numericalFluxR[unknown];
118 Qouttest[upperVoxelInImage * (NumberOfUnknowns + NumberOfAuxiliaryVariables) + unknown] += Qcmp
119 [upperVoxelInImage * (NumberOfUnknowns + NumberOfAuxiliaryVariables) + unknown];
120 }
121 }
122 }
123 }
124 }
125
126 double* Qouttest2 = new double[500];
127 for (unsigned int i = 0; i < 500; ++i) {
128 Qouttest2[i] = 0.0;
129 }
130
131 for (int shift = 0; shift < 2; shift++) {
132 for (int x = shift; x <= NumberOfVolumesPerAxisInPatch; x += 2) {
133 for (int y = 0; y < NumberOfVolumesPerAxisInPatch; y++) {
134 const int leftVoxelInImage = x - 1 + y * NumberOfVolumesPerAxisInPatch;
135 const int rightVoxelInImage = x + y * NumberOfVolumesPerAxisInPatch;
136 for (int unknown = 0; unknown < NumberOfUnknowns; unknown++) {
137 if (x > 0) {
138 Qouttest2[leftVoxelInImage * (NumberOfUnknowns + NumberOfAuxiliaryVariables) + unknown] -= dt / volumeH(0) * numericalFluxL[unknown];
139 Qouttest2[leftVoxelInImage * (NumberOfUnknowns + NumberOfAuxiliaryVariables) + unknown] += Qcmp
140 [leftVoxelInImage * (NumberOfUnknowns + NumberOfAuxiliaryVariables) + unknown];
141 }
142 if (x < NumberOfVolumesPerAxisInPatch) {
143 Qouttest2[rightVoxelInImage * (NumberOfUnknowns + NumberOfAuxiliaryVariables) + unknown] += dt / volumeH(0) * numericalFluxR[unknown];
144 Qouttest2[rightVoxelInImage * (NumberOfUnknowns + NumberOfAuxiliaryVariables) + unknown] += Qcmp
145 [rightVoxelInImage * (NumberOfUnknowns + NumberOfAuxiliaryVariables) + unknown];
146 }
147 }
148 }
149 }
150 }
151
152 // Iterate over other normal
153 for (int shift = 0; shift < 2; shift++) {
154 for (int y = shift; y <= NumberOfVolumesPerAxisInPatch; y += 2) {
155 for (int x = 0; x < NumberOfVolumesPerAxisInPatch; x++) {
156 const int lowerVoxelInImage = x + (y - 1) * NumberOfVolumesPerAxisInPatch;
157 const int upperVoxelInImage = x + y * NumberOfVolumesPerAxisInPatch;
158 for (int unknown = 0; unknown < NumberOfUnknowns; unknown++) {
159 if (y > 0) {
160 Qouttest2[lowerVoxelInImage * (NumberOfUnknowns + NumberOfAuxiliaryVariables) + unknown] -= dt / volumeH(0) * numericalFluxL[unknown];
161 Qouttest2[lowerVoxelInImage * (NumberOfUnknowns + NumberOfAuxiliaryVariables) + unknown] += Qcmp
162 [lowerVoxelInImage * (NumberOfUnknowns + NumberOfAuxiliaryVariables) + unknown];
163 }
164 if (y < NumberOfVolumesPerAxisInPatch) {
165 Qouttest2[upperVoxelInImage * (NumberOfUnknowns + NumberOfAuxiliaryVariables) + unknown] += dt / volumeH(0) * numericalFluxR[unknown];
166 Qouttest2[upperVoxelInImage * (NumberOfUnknowns + NumberOfAuxiliaryVariables) + unknown] += Qcmp
167 [upperVoxelInImage * (NumberOfUnknowns + NumberOfAuxiliaryVariables) + unknown];
168 }
169 }
170 }
171 }
172 }
173
174 // Elementwise comparison
175 for (unsigned int i = 0; i < 500; ++i) {
176 validateNumericalEquals(Qout[i], Qouttest[i]);
177 validateNumericalEquals(Qouttest[i], Qouttest2[i]);
178 }
179
180 delete[] Qout;
181 delete[] Qouttest;
184}
#define testMethod(name)
Run a test method and check for errors.
Definition TestMacros.h:24
#define validateNumericalEquals(actualValue, validValue)
Definition TestMacros.h:428
virtual void run() override
This routine is triggered by the TestCaseCollection.
void freeMemory(void *data, MemoryLocation location)
@ Heap
Create data on the heap of the local device.
Simple vector class.
Definition Vector.h:134