Peano
Loading...
Searching...
No Matches
FD4_patchwise_functors.cpp
Go to the documentation of this file.
1#include "../FD4.h"
2
5
8
11 int numberOfGridCellsPerPatchPerAxis,
12 int haloSize,
13 int unknowns,
14 int auxiliaryVariables,
15 double KOSigma,
16 bool evaluateFlux,
17 bool evaluateDifferentialSource, //for ncp
18 bool evaluateAlgebraicSource, //real source
19 bool copyOldTimeStepAndScaleWithTimeStepSize,
20 DifferentialSourceTermVariant variant,
21 Flux flux,
22 NonconservativeProduct DifferentialSource,
23 Source AlgebraicSource
24) {
25 static tarch::logging::Log _log( "exahype2::fd" );
26
27 assertionEquals(haloSize,3);
28
29 bool evaluateKODissipation = KOSigma>0.0;
30
31 exahype2::enumerator::AoSLexicographicEnumerator QInEnumerator (1,numberOfGridCellsPerPatchPerAxis,haloSize,unknowns,auxiliaryVariables);
32 exahype2::enumerator::AoSLexicographicEnumerator QOutEnumerator(1,numberOfGridCellsPerPatchPerAxis,0, unknowns,0);
33
34 logDebug( "timeStep_patchwise_heap_functors(...)", "size of temp flux field=" << QOutEnumerator.size() );
35 double* tempFluxX = evaluateFlux ? tarch::allocateMemory<double>( QOutEnumerator.size(),tarch::MemoryLocation::Heap ) : nullptr;
36 double* tempFluxY = evaluateFlux ? tarch::allocateMemory<double>( QOutEnumerator.size(),tarch::MemoryLocation::Heap ) : nullptr;
37 double* tempFluxZ = (evaluateFlux and Dimensions==3) ? tarch::allocateMemory<double>( QOutEnumerator.size(),tarch::MemoryLocation::Heap ) : nullptr;
38 double* tempDifferentialSourceX = evaluateDifferentialSource ? tarch::allocateMemory<double>( QOutEnumerator.size(),tarch::MemoryLocation::Heap ) : nullptr;
39 double* tempDifferentialSourceY = evaluateDifferentialSource ? tarch::allocateMemory<double>( QOutEnumerator.size(),tarch::MemoryLocation::Heap ) : nullptr;
40 double* tempDifferentialSourceZ = (evaluateDifferentialSource and Dimensions==3) ? tarch::allocateMemory<double>( QOutEnumerator.size(),tarch::MemoryLocation::Heap ) : nullptr;
41 double* tempKODissipationX = evaluateKODissipation ? tarch::allocateMemory<double>( QOutEnumerator.size(),tarch::MemoryLocation::Heap ) : nullptr;
42 double* tempKODissipationY = evaluateKODissipation ? tarch::allocateMemory<double>( QOutEnumerator.size(),tarch::MemoryLocation::Heap ) : nullptr;
43 double* tempKODissipationZ = (evaluateKODissipation and Dimensions==3) ? tarch::allocateMemory<double>( QOutEnumerator.size(),tarch::MemoryLocation::Heap ) : nullptr;
44
45 for (int patchIndex=0; patchIndex<patchData.numberOfCells; patchIndex++) {
46 const double timeStepSize = copyOldTimeStepAndScaleWithTimeStepSize ? patchData.dt[patchIndex] : 1.0;
47
48 // STEP 1: copy solution from old array to new array, or set all to zero
49 if ( Dimensions==3 and copyOldTimeStepAndScaleWithTimeStepSize ) {
50 for (int x = 0; x < numberOfGridCellsPerPatchPerAxis; x++)
51 for (int y = 0; y < numberOfGridCellsPerPatchPerAxis; y++)
52 for (int z = 0; z < numberOfGridCellsPerPatchPerAxis; z++)
53 for (int unknown=0; unknown<unknowns; unknown++) {
55 patchData.QIn[patchIndex],
56 QInEnumerator,
57 patchIndex,
58 gridCellIndex3d(x,y,z),
59 unknown,
60 patchData.QOut[patchIndex],
61 QOutEnumerator
62 );
63 }
64 }
65 else if (Dimensions==3 and not copyOldTimeStepAndScaleWithTimeStepSize){
66 for (int x = 0; x < numberOfGridCellsPerPatchPerAxis; x++)
67 for (int y = 0; y < numberOfGridCellsPerPatchPerAxis; y++)
68 for (int z = 0; z < numberOfGridCellsPerPatchPerAxis; z++)
69 {
70 for (int unknown=0; unknown<unknowns; unknown++) {
72 patchIndex,
73 gridCellIndex3d(x,y,z),
74 unknown,
75 patchData.QOut[patchIndex],
76 QOutEnumerator
77 );
78 }
79 }
80 }
81 else if (Dimensions==2 and copyOldTimeStepAndScaleWithTimeStepSize) {
82 for (int x = 0; x < numberOfGridCellsPerPatchPerAxis; x++)
83 for (int y = 0; y < numberOfGridCellsPerPatchPerAxis; y++)
84 for (int unknown=0; unknown<unknowns; unknown++) {
86 patchData.QIn[patchIndex],
87 QInEnumerator,
88 patchIndex,
89 gridCellIndex2d(x,y),
90 unknown,
91 patchData.QOut[patchIndex],
92 QOutEnumerator
93 );
94 }
95 }
96 else {
97 for (int x = 0; x < numberOfGridCellsPerPatchPerAxis; x++)
98 for (int y = 0; y < numberOfGridCellsPerPatchPerAxis; y++){
99 for (int unknown=0; unknown<unknowns; unknown++) {
101 patchIndex,
102 gridCellIndex2d(x,y),
103 unknown,
104 patchData.QOut[patchIndex],
105 QOutEnumerator
106 );
107 }
108 }
109 }
110
111 // STEP 1.5: Calculate and update auxiliary variables
112 if (Dimensions==2) {
113
114 for (int x = 0; x < numberOfGridCellsPerPatchPerAxis; x++)
115 for (int y = 0; y < numberOfGridCellsPerPatchPerAxis; y++) {
116 internal::computeAuxiliaryVariables_LoopBody(
117 patchData.QIn[patchIndex],
118 QInEnumerator,
119 patchData.cellCentre[patchIndex],
120 patchData.cellSize[patchIndex],
121 patchIndex,
122 gridCellIndex2d(x,y),
123 0, // normal
124 QOutEnumerator
125 );
126 }
127
128 for (int x = 0; x < numberOfGridCellsPerPatchPerAxis; x++)
129 for (int y = 0; y < numberOfGridCellsPerPatchPerAxis; y++) {
130 internal::computeAuxiliaryVariables_LoopBody(
131 patchData.QIn[patchIndex],
132 QInEnumerator,
133 patchData.cellCentre[patchIndex],
134 patchData.cellSize[patchIndex],
135 patchIndex,
136 gridCellIndex2d(x,y),
137 1, // normal
138 QOutEnumerator
139 );
140 }
141
142 } else if (Dimensions==3) {
143
144 for (int x = 0; x < numberOfGridCellsPerPatchPerAxis; x++)
145 for (int y = 0; y < numberOfGridCellsPerPatchPerAxis; y++)
146 for (int z = 0; z < numberOfGridCellsPerPatchPerAxis; z++) {
147 internal::computeAuxiliaryVariables_LoopBody(
148 patchData.QIn[patchIndex],
149 QInEnumerator,
150 patchData.cellCentre[patchIndex],
151 patchData.cellSize[patchIndex],
152 patchIndex,
153 gridCellIndex3d(x,y,z),
154 0, // normal
155 QOutEnumerator
156 );
157 }
158
159 for (int x = 0; x < numberOfGridCellsPerPatchPerAxis; x++)
160 for (int y = 0; y < numberOfGridCellsPerPatchPerAxis; y++)
161 for (int z = 0; z < numberOfGridCellsPerPatchPerAxis; z++) {
162 internal::computeAuxiliaryVariables_LoopBody(
163 patchData.QIn[patchIndex],
164 QInEnumerator,
165 patchData.cellCentre[patchIndex],
166 patchData.cellSize[patchIndex],
167 patchIndex,
168 gridCellIndex3d(x,y,z),
169 1, // normal
170 QOutEnumerator
171 );
172 }
173
174 for (int x = 0; x < numberOfGridCellsPerPatchPerAxis; x++)
175 for (int y = 0; y < numberOfGridCellsPerPatchPerAxis; y++)
176 for (int z = 0; z < numberOfGridCellsPerPatchPerAxis; z++) {
177 internal::computeAuxiliaryVariables_LoopBody(
178 patchData.QIn[patchIndex],
179 QInEnumerator,
180 patchData.cellCentre[patchIndex],
181 patchData.cellSize[patchIndex],
182 patchIndex,
183 gridCellIndex3d(x,y,z),
184 2, // normal
185 QOutEnumerator
186 );
187 }
188
189 }
190
191
192
193 // Step 2: Add source terms
194 if (Dimensions==2 and evaluateAlgebraicSource) {
195 for (int x = 0; x < numberOfGridCellsPerPatchPerAxis; x++)
196 for (int y = 0; y < numberOfGridCellsPerPatchPerAxis; y++) {
198 patchData.QIn[patchIndex],
199 QInEnumerator,
200 AlgebraicSource,
201 patchData.cellCentre[patchIndex],
202 patchData.cellSize[patchIndex],
203 patchIndex,
204 gridCellIndex2d(x,y),
205 patchData.t[patchIndex],
206 timeStepSize,
207 patchData.QOut[patchIndex],
208 QOutEnumerator
209 );
210 }
211 }
212 else if (Dimensions==3 and evaluateAlgebraicSource) {
213 for (int x = 0; x < numberOfGridCellsPerPatchPerAxis; x++)
214 for (int y = 0; y < numberOfGridCellsPerPatchPerAxis; y++)
215 for (int z = 0; z < numberOfGridCellsPerPatchPerAxis; z++) {
217 patchData.QIn[patchIndex],
218 QInEnumerator,
219 AlgebraicSource,
220 patchData.cellCentre[patchIndex],
221 patchData.cellSize[patchIndex],
222 patchIndex,
223 gridCellIndex3d(x,y,z),
224 patchData.t[patchIndex],
225 timeStepSize,
226 patchData.QOut[patchIndex],
227 QOutEnumerator
228 );
229 }
230 }
231
232 // STEP 3: Evaluate flux and update solution accordingly
233 if (evaluateFlux and Dimensions==2) {
234 for (int x = 0; x < numberOfGridCellsPerPatchPerAxis; x++)
235 for (int y = 0; y < numberOfGridCellsPerPatchPerAxis; y++) {
236 internal::computeFlux_LoopBody(
237 patchData.QIn[patchIndex],
238 QInEnumerator,
239 flux,
240 patchData.cellCentre[patchIndex],
241 patchData.cellSize[patchIndex],
242 patchIndex,
243 gridCellIndex2d(x,y),
244 patchData.t[patchIndex],
245 timeStepSize,
246 0, // normal
247 tempFluxX,
248 QOutEnumerator
249 );
250 }
251 for (int x = 0; x < numberOfGridCellsPerPatchPerAxis; x++)
252 for (int y = 0; y < numberOfGridCellsPerPatchPerAxis; y++) {
253 internal::computeFlux_LoopBody(
254 patchData.QIn[patchIndex],
255 QInEnumerator,
256 flux,
257 patchData.cellCentre[patchIndex],
258 patchData.cellSize[patchIndex],
259 patchIndex,
260 gridCellIndex2d(x,y),
261 patchData.t[patchIndex],
262 timeStepSize,
263 1, // normal
264 tempFluxY,
265 QOutEnumerator
266 );
267 }
268 for (int x = 0; x < numberOfGridCellsPerPatchPerAxis; x++)
269 for (int y = 0; y < numberOfGridCellsPerPatchPerAxis; y++)
270 for (int unknown=0; unknown<unknowns; unknown++) {
271 // @todo 2d missing
272 internal::updateSolutionWithFlux_LoopBody(
273 tempFluxX, tempFluxY, tempFluxZ,
274 QOutEnumerator,
275 patchData.cellCentre[patchIndex],
276 patchData.cellSize[patchIndex],
277 patchIndex,
278 gridCellIndex2d(x,y),
279 unknown,
280 timeStepSize,
281 *(patchData.QOut + patchIndex),
282 QOutEnumerator
283 );
284 }
285 }
286 else if (evaluateFlux and Dimensions==3) {
287 for (int x = 0; x < numberOfGridCellsPerPatchPerAxis; x++)
288 for (int y = 0; y < numberOfGridCellsPerPatchPerAxis; y++)
289 for (int z = 0; z < numberOfGridCellsPerPatchPerAxis; z++) {
290 internal::computeFlux_LoopBody(
291 patchData.QIn[patchIndex],
292 QInEnumerator,
293 flux,
294 patchData.cellCentre[patchIndex],
295 patchData.cellSize[patchIndex],
296 patchIndex,
297 gridCellIndex3d(x,y,z),
298 patchData.t[patchIndex],
299 timeStepSize,
300 0, // normal
301 tempFluxX,
302 QOutEnumerator
303 );
304 }
305 for (int x = 0; x < numberOfGridCellsPerPatchPerAxis; x++)
306 for (int y = 0; y < numberOfGridCellsPerPatchPerAxis; y++)
307 for (int z = 0; z < numberOfGridCellsPerPatchPerAxis; z++) {
308 internal::computeFlux_LoopBody(
309 patchData.QIn[patchIndex],
310 QInEnumerator,
311 flux,
312 patchData.cellCentre[patchIndex],
313 patchData.cellSize[patchIndex],
314 patchIndex,
315 gridCellIndex3d(x,y,z),
316 patchData.t[patchIndex],
317 timeStepSize,
318 1, // normal
319 tempFluxY,
320 QOutEnumerator
321 );
322 }
323 for (int x = 0; x < numberOfGridCellsPerPatchPerAxis; x++)
324 for (int y = 0; y < numberOfGridCellsPerPatchPerAxis; y++)
325 for (int z = 0; z < numberOfGridCellsPerPatchPerAxis; z++) {
326 internal::computeFlux_LoopBody(
327 patchData.QIn[patchIndex],
328 QInEnumerator,
329 flux,
330 patchData.cellCentre[patchIndex],
331 patchData.cellSize[patchIndex],
332 patchIndex,
333 gridCellIndex3d(x,y,z),
334 patchData.t[patchIndex],
335 timeStepSize,
336 2, // normal
337 tempFluxZ,
338 QOutEnumerator
339 );
340 }
341
342 for (int x = 0; x < numberOfGridCellsPerPatchPerAxis; x++)
343 for (int y = 0; y < numberOfGridCellsPerPatchPerAxis; y++)
344 for (int z = 0; z < numberOfGridCellsPerPatchPerAxis; z++)
345 for (int unknown=0; unknown<unknowns; unknown++) {
346 internal::updateSolutionWithFlux_LoopBody(
347 tempFluxX, tempFluxY, tempFluxZ,
348 QOutEnumerator,
349 patchData.cellCentre[patchIndex],
350 patchData.cellSize[patchIndex],
351 patchIndex,
352 gridCellIndex3d(x,y,z),
353 unknown,
354 timeStepSize,
355 *(patchData.QOut + patchIndex),
356 QOutEnumerator
357 );
358 }
359 }
360
361 // STEP 4: Evaluate differential source and update solution accordingly
362 if (evaluateDifferentialSource and Dimensions==2) {
363 for (int x = 0; x < numberOfGridCellsPerPatchPerAxis; x++)
364 for (int y = 0; y < numberOfGridCellsPerPatchPerAxis; y++) {
365 internal::computeDifferentialSourceTerm_LoopBody(
366 patchData.QIn[patchIndex],
367 QInEnumerator,
368 DifferentialSource,
369 patchData.cellCentre[patchIndex],
370 patchData.cellSize[patchIndex],
371 patchIndex,
372 gridCellIndex2d(x,y),
373 patchData.t[patchIndex],
374 timeStepSize,
375 0, // normal
376 tempDifferentialSourceX,
377 QOutEnumerator,
378 variant
379 );
380 }
381 for (int x = 0; x < numberOfGridCellsPerPatchPerAxis; x++)
382 for (int y = 0; y < numberOfGridCellsPerPatchPerAxis; y++) {
383 internal::computeDifferentialSourceTerm_LoopBody(
384 patchData.QIn[patchIndex],
385 QInEnumerator,
386 DifferentialSource,
387 patchData.cellCentre[patchIndex],
388 patchData.cellSize[patchIndex],
389 patchIndex,
390 gridCellIndex2d(x,y),
391 patchData.t[patchIndex],
392 timeStepSize,
393 1, // normal
394 tempDifferentialSourceY,
395 QOutEnumerator,
396 variant
397 );
398 }
399 for (int x = 0; x < numberOfGridCellsPerPatchPerAxis; x++)
400 for (int y = 0; y < numberOfGridCellsPerPatchPerAxis; y++)
401 for (int unknown=0; unknown<unknowns; unknown++) {
402 internal::updateSolutionWithDifferentialSourceTerm_LoopBody(
403 tempDifferentialSourceX, tempDifferentialSourceY, tempDifferentialSourceZ,
404 QOutEnumerator,
405 patchData.cellCentre[patchIndex],
406 patchData.cellSize[patchIndex],
407 patchIndex,
408 gridCellIndex2d(x,y),
409 unknown,
410 timeStepSize,
411 *(patchData.QOut + patchIndex),
412 QOutEnumerator
413 );
414 }
415 }
416 else if (evaluateDifferentialSource and Dimensions==3) {
417 for (int x = 0; x < numberOfGridCellsPerPatchPerAxis; x++)
418 for (int y = 0; y < numberOfGridCellsPerPatchPerAxis; y++)
419 for (int z = 0; z < numberOfGridCellsPerPatchPerAxis; z++) {
420 internal::computeDifferentialSourceTerm_LoopBody(
421 patchData.QIn[patchIndex],
422 QInEnumerator,
423 DifferentialSource,
424 patchData.cellCentre[patchIndex],
425 patchData.cellSize[patchIndex],
426 patchIndex,
427 gridCellIndex3d(x,y,z),
428 patchData.t[patchIndex],
429 timeStepSize,
430 0, // normal
431 tempDifferentialSourceX,
432 QOutEnumerator,
433 variant
434 );
435 }
436 for (int x = 0; x < numberOfGridCellsPerPatchPerAxis; x++)
437 for (int y = 0; y < numberOfGridCellsPerPatchPerAxis; y++)
438 for (int z = 0; z < numberOfGridCellsPerPatchPerAxis; z++) {
439 internal::computeDifferentialSourceTerm_LoopBody(
440 patchData.QIn[patchIndex],
441 QInEnumerator,
442 DifferentialSource,
443 patchData.cellCentre[patchIndex],
444 patchData.cellSize[patchIndex],
445 patchIndex,
446 gridCellIndex3d(x,y,z),
447 patchData.t[patchIndex],
448 timeStepSize,
449 1, // normal
450 tempDifferentialSourceY,
451 QOutEnumerator,
452 variant
453 );
454 }
455 for (int x = 0; x < numberOfGridCellsPerPatchPerAxis; x++)
456 for (int y = 0; y < numberOfGridCellsPerPatchPerAxis; y++)
457 for (int z = 0; z < numberOfGridCellsPerPatchPerAxis; z++) {
458 internal::computeDifferentialSourceTerm_LoopBody(
459 patchData.QIn[patchIndex],
460 QInEnumerator,
461 DifferentialSource,
462 patchData.cellCentre[patchIndex],
463 patchData.cellSize[patchIndex],
464 patchIndex,
465 gridCellIndex3d(x,y,z),
466 patchData.t[patchIndex],
467 timeStepSize,
468 2, // normal
469 tempDifferentialSourceZ,
470 QOutEnumerator,
471 variant
472 );
473 }
474 for (int x = 0; x < numberOfGridCellsPerPatchPerAxis; x++)
475 for (int y = 0; y < numberOfGridCellsPerPatchPerAxis; y++)
476 for (int z = 0; z < numberOfGridCellsPerPatchPerAxis; z++)
477 for (int unknown=0; unknown<unknowns; unknown++) {
478 internal::updateSolutionWithDifferentialSourceTerm_LoopBody(
479 tempDifferentialSourceX, tempDifferentialSourceY, tempDifferentialSourceZ,
480 QOutEnumerator,
481 patchData.cellCentre[patchIndex],
482 patchData.cellSize[patchIndex],
483 patchIndex,
484 gridCellIndex3d(x,y,z),
485 unknown,
486 timeStepSize,
487 *(patchData.QOut + patchIndex),
488 QOutEnumerator
489 );
490 }
491 }
492
493 // STEP 5: compute the Kreiss Oliger dissipation and add it to the solution
494 if (evaluateKODissipation and Dimensions==2) {
495 for (int x = 0; x < numberOfGridCellsPerPatchPerAxis; x++)
496 for (int y = 0; y < numberOfGridCellsPerPatchPerAxis; y++) {
497 internal::computeKreissOligerDissipationTerm_LoopBody(
498 patchData.QIn[patchIndex],
499 QInEnumerator,
500 patchData.cellCentre[patchIndex],
501 patchData.cellSize[patchIndex],
502 patchIndex,
503 gridCellIndex2d(x,y),
504 patchData.t[patchIndex],
505 timeStepSize,
506 0, // normal
507 tempKODissipationX,
508 QOutEnumerator
509 );
510 }
511 for (int x = 0; x < numberOfGridCellsPerPatchPerAxis; x++)
512 for (int y = 0; y < numberOfGridCellsPerPatchPerAxis; y++) {
513 internal::computeKreissOligerDissipationTerm_LoopBody(
514 patchData.QIn[patchIndex],
515 QInEnumerator,
516 patchData.cellCentre[patchIndex],
517 patchData.cellSize[patchIndex],
518 patchIndex,
519 gridCellIndex2d(x,y),
520 patchData.t[patchIndex],
521 timeStepSize,
522 1, // normal
523 tempKODissipationY,
524 QOutEnumerator
525 );
526 }
527 for (int x = 0; x < numberOfGridCellsPerPatchPerAxis; x++)
528 for (int y = 0; y < numberOfGridCellsPerPatchPerAxis; y++)
529 for (int unknown=0; unknown<unknowns; unknown++) {
530 internal::updateSolutionWithKODissipationTerm_LoopBody(
531 KOSigma, tempKODissipationX, tempKODissipationY, tempKODissipationZ,
532 QOutEnumerator,
533 patchData.cellCentre[patchIndex],
534 patchData.cellSize[patchIndex],
535 patchIndex,
536 gridCellIndex2d(x,y),
537 unknown,
538 timeStepSize,
539 *(patchData.QOut + patchIndex),
540 QOutEnumerator
541 );
542 }
543 }
544 else if (evaluateKODissipation and Dimensions==3) {
545 for (int x = 0; x < numberOfGridCellsPerPatchPerAxis; x++)
546 for (int y = 0; y < numberOfGridCellsPerPatchPerAxis; y++)
547 for (int z = 0; z < numberOfGridCellsPerPatchPerAxis; z++) {
548 internal::computeKreissOligerDissipationTerm_LoopBody(
549 patchData.QIn[patchIndex],
550 QInEnumerator,
551 patchData.cellCentre[patchIndex],
552 patchData.cellSize[patchIndex],
553 patchIndex,
554 gridCellIndex3d(x,y,z),
555 patchData.t[patchIndex],
556 timeStepSize,
557 0, // normal
558 tempKODissipationX,
559 QOutEnumerator
560 );
561 }
562 for (int x = 0; x < numberOfGridCellsPerPatchPerAxis; x++)
563 for (int y = 0; y < numberOfGridCellsPerPatchPerAxis; y++)
564 for (int z = 0; z < numberOfGridCellsPerPatchPerAxis; z++) {
565 internal::computeKreissOligerDissipationTerm_LoopBody(
566 patchData.QIn[patchIndex],
567 QInEnumerator,
568 patchData.cellCentre[patchIndex],
569 patchData.cellSize[patchIndex],
570 patchIndex,
571 gridCellIndex3d(x,y,z),
572 patchData.t[patchIndex],
573 timeStepSize,
574 1, // normal
575 tempKODissipationY,
576 QOutEnumerator
577 );
578 }
579 for (int x = 0; x < numberOfGridCellsPerPatchPerAxis; x++)
580 for (int y = 0; y < numberOfGridCellsPerPatchPerAxis; y++)
581 for (int z = 0; z < numberOfGridCellsPerPatchPerAxis; z++) {
582 internal::computeKreissOligerDissipationTerm_LoopBody(
583 patchData.QIn[patchIndex],
584 QInEnumerator,
585 patchData.cellCentre[patchIndex],
586 patchData.cellSize[patchIndex],
587 patchIndex,
588 gridCellIndex3d(x,y,z),
589 patchData.t[patchIndex],
590 timeStepSize,
591 2, // normal
592 tempKODissipationZ,
593 QOutEnumerator
594 );
595 }
596 for (int x = 0; x < numberOfGridCellsPerPatchPerAxis; x++)
597 for (int y = 0; y < numberOfGridCellsPerPatchPerAxis; y++)
598 for (int z = 0; z < numberOfGridCellsPerPatchPerAxis; z++)
599 for (int unknown=0; unknown<unknowns; unknown++) {
600 internal::updateSolutionWithKODissipationTerm_LoopBody(
601 KOSigma, tempKODissipationX, tempKODissipationY, tempKODissipationZ,
602 QOutEnumerator,
603 patchData.cellCentre[patchIndex],
604 patchData.cellSize[patchIndex],
605 patchIndex,
606 gridCellIndex3d(x,y,z),
607 unknown,
608 timeStepSize,
609 *(patchData.QOut + patchIndex),
610 QOutEnumerator
611 );
612 }
613 }
614 }
615
616 if (tempFluxX!=nullptr) tarch::freeMemory(tempFluxX, tarch::MemoryLocation::Heap);
617 if (tempFluxY!=nullptr) tarch::freeMemory(tempFluxY, tarch::MemoryLocation::Heap);
618 if (tempFluxZ!=nullptr) tarch::freeMemory(tempFluxZ, tarch::MemoryLocation::Heap);
619 if (tempDifferentialSourceX!=nullptr) tarch::freeMemory(tempDifferentialSourceX, tarch::MemoryLocation::Heap);
620 if (tempDifferentialSourceY!=nullptr) tarch::freeMemory(tempDifferentialSourceY, tarch::MemoryLocation::Heap);
621 if (tempDifferentialSourceZ!=nullptr) tarch::freeMemory(tempDifferentialSourceZ, tarch::MemoryLocation::Heap);
622 if (tempKODissipationX!=nullptr) tarch::freeMemory(tempKODissipationX, tarch::MemoryLocation::Heap);
623 if (tempKODissipationY!=nullptr) tarch::freeMemory(tempKODissipationY, tarch::MemoryLocation::Heap);
624 if (tempKODissipationZ!=nullptr) tarch::freeMemory(tempKODissipationZ, tarch::MemoryLocation::Heap);
625}
#define assertionEquals(lhs, rhs)
#define logDebug(methodName, logMacroMessageStream)
Definition Log.h:50
Log Device.
Definition Log.h:516
tarch::logging::Log _log("exahype2::fv")
void timeStep_patchwise_heap_functors(::exahype2::CellData< double, double > &patchData, int numberOfGridCellsPerPatchPerAxis, int haloSize, int unknowns, int auxiliaryVariables, double KOSigma, bool evaluateFlux, bool evaluateDifferentialSource, bool evaluateAlgebraicSource, bool copyOldTimeStepAndScaleWithTimeStepSize, DifferentialSourceTermVariant variant, Flux flux, NonconservativeProduct DifferentialSource, Source AlgebraicSource)
Fourth-order Finite Differences.
static void copySolution_LoopBody(const double *__restrict__ QIn, const QInEnumeratorType &QInEnumerator, int patchIndex, const tarch::la::Vector< Dimensions, int > &volumeIndex, int unknown, double *__restrict__ QOut, const QOutEnumeratorType &QOutEnumerator) InlineMethod
static void addAlgebraicSourceTerm_LoopBody(const double *__restrict__ QIn, const QInEnumeratorType &QInEnumerator, exahype2::fd::Source AlgebraicSource, const tarch::la::Vector< Dimensions, double > &patchCentre, const tarch::la::Vector< Dimensions, double > &patchSize, int patchIndex, const tarch::la::Vector< Dimensions, int > &volumeIndex, double t, double dt, double *__restrict__ QOut, const QOutEnumeratorType &QOutEnumerator) InlineMethod
Copy previous solution over into new time step and add source term.
static void clearSolution_LoopBody(int patchIndex, const tarch::la::Vector< Dimensions, int > &volumeIndex, int unknown, double *__restrict__ QOut, const QOutEnumeratorType &QOutEnumerator) InlineMethod
The complicated way to write =0.
void freeMemory(void *data, MemoryLocation location, int device=accelerator::Device::HostDevice)
Free memory.
@ Heap
Create data on the heap of the local device.
Representation of a number of cells which contains all information that's required to process the sto...
Definition CellData.h:77
outType ** QOut
Out values.
Definition CellData.h:116
inType ** QIn
QIn may not be const, as some kernels delete it straightaway once the input data has been handled.
Definition CellData.h:82
const int numberOfCells
As we store data as SoA, we have to know how big the actual arrays are.
Definition CellData.h:99
tarch::la::Vector< Dimensions, double > * cellCentre
Definition CellData.h:83
tarch::la::Vector< Dimensions, double > * cellSize
Definition CellData.h:84