Peano 4
Loading...
Searching...
No Matches
Restriction.cpp
Go to the documentation of this file.
1#include "Restriction.h"
2#include "Interpolation.h"
3#include "Enumeration.h"
4
5#include "peano4/utils/Loop.h"
6#include "tarch/Assertions.h"
8
9namespace {
10 [[maybe_unused]] tarch::logging::Log _log( "toolbox::blockstructured" );
11}
12
14 [[maybe_unused]] const peano4::datamanagement::CellMarker& marker,
15 [[maybe_unused]] int numberOfDoFsPerAxisInPatch,
16 [[maybe_unused]] int unknowns,
17 [[maybe_unused]] const double* __restrict__ tangentialRestrictionMatrix1d,
18 [[maybe_unused]] double* fineGridValues,
19 [[maybe_unused]] double* coarseGridValues
20) {
21 assertionMsg( false, "@Han, we have to implement this one" );
22}
23
26 int numberOfDoFsPerAxisInPatch,
27 int overlap,
28 int unknowns,
29 const double* __restrict__ normalRestrictionMatrix1d,
30 const double* __restrict__ tangentialRestrictionMatrix1d,
31 double* fineGridValues,
32 double* coarseGridValues
33) {
34 restrictInnerHalfOfHaloLayer_AoS_tensor_product( marker, numberOfDoFsPerAxisInPatch, overlap, unknowns, normalRestrictionMatrix1d, tangentialRestrictionMatrix1d, fineGridValues, coarseGridValues, false );
35 restrictInnerHalfOfHaloLayer_AoS_tensor_product( marker, numberOfDoFsPerAxisInPatch, overlap, unknowns, normalRestrictionMatrix1d, tangentialRestrictionMatrix1d, fineGridValues, coarseGridValues, true );
36}
37
40 int numberOfDoFsPerAxisInPatch,
41 int overlap,
42 int unknowns,
43 const double* __restrict__ normalRestrictionMatrix1d,
44 const double* __restrict__ tangentialRestrictionMatrix1d,
45 double* fineGridValues,
46 double* coarseGridValues,
47 bool swapInsideOutside
48) {
49 logTraceInWith3Arguments( "restrictInnerHalfOfHaloLayer_AoS_tensor_product(...)", marker.toString(), numberOfDoFsPerAxisInPatch, overlap );
50
51 const int normal = marker.getSelectedFaceNumber() % Dimensions;
52
53 dfore(kCoarse, numberOfDoFsPerAxisInPatch, normal, 0) {
54 for (int iCoarse=0; iCoarse<overlap; iCoarse++) {
56 dest(normal) = ((marker.getSelectedFaceNumber() < Dimensions) xor swapInsideOutside) ? iCoarse + overlap: overlap - iCoarse - 1;
57
58 int destLinearised = serialiseVoxelIndexInOverlap(
59 dest,
60 numberOfDoFsPerAxisInPatch, overlap, normal
61 );
62
63 logDebug( "restrictInnerHalfOfHaloLayer_AoS_tensor_product(...)", "coarse volume/dof " << dest << " prior to update: " << coarseGridValues[destLinearised*5+0] );
64 assertion( coarseGridValues[destLinearised*5+0]==coarseGridValues[destLinearised*5+0] );
65
66 dfore(kFine, numberOfDoFsPerAxisInPatch, normal, 0) {
67 for (int iFine = 0; iFine < 2 * overlap; iFine++) {
69 src(normal) = ((marker.getSelectedFaceNumber() < Dimensions) xor swapInsideOutside) ? iFine : 2 * overlap - iFine - 1;
70 int srcLinearised = serialiseVoxelIndexInOverlap(
71 src,
72 numberOfDoFsPerAxisInPatch, overlap, normal
73 );
74
75 double weight = 1.0;
76 for (int d=0; d<Dimensions; d++) {
77 if (d==normal) {
78 int col = iFine;
79 int row = iCoarse;
80
81 assertion4(col>=0,row,col,src,dest);
82 assertion4(row>=0,row,col,src,dest);
83
84 assertion4(col<2*overlap,row,col,src,dest);
85 assertion4(row< overlap,row,col,src,dest);
86
87 int index = col + row * 2 * overlap;
88 weight *= normalRestrictionMatrix1d[ index ];
89 logDebug( "restrictInnerHalfOfHaloLayer_AoS_tensor_product(...)", "use normal matrix entry " << row << ", " << col << ", i.e. index " << index << ": " << normalRestrictionMatrix1d[ index ]);
90 } else {
91 assertion( marker.getRelativePositionWithinFatherCell()(d)>=0 );
92 assertion( marker.getRelativePositionWithinFatherCell()(d)<3 );
93
94 assertion( src(d)>=0 );
95 int col = src(d)+ marker.getRelativePositionWithinFatherCell()(d) * numberOfDoFsPerAxisInPatch;
96 int row = dest(d);
97 assertion5(col>=0,row,col,src,dest, marker.getRelativePositionWithinFatherCell());
98 assertion5(row>=0,row,col,src,dest, marker.getRelativePositionWithinFatherCell());
99 assertion5(col<3*numberOfDoFsPerAxisInPatch,row,col,src,dest, marker.getRelativePositionWithinFatherCell());
100 assertion5(row< numberOfDoFsPerAxisInPatch,row,col,src,dest, marker.getRelativePositionWithinFatherCell());
101 int index = col + row * 3 * numberOfDoFsPerAxisInPatch;
102 weight *= tangentialRestrictionMatrix1d[ index ];
103 logDebug( "restrictInnerHalfOfHaloLayer_AoS_tensor_product(...)", "use tangential matrix entry " << row << ", " << row << ", i.e. index " << index << ": " << tangentialRestrictionMatrix1d[ index ]);
104 }
105 }
106
107 assertion2(weight==weight, srcLinearised, destLinearised);
108
109 for (int unknown=0; unknown<unknowns; unknown++) {
110 assertion2( fineGridValues[srcLinearised*unknowns+unknown] == fineGridValues[srcLinearised*unknowns+unknown], srcLinearised, destLinearised );
111 coarseGridValues[destLinearised*unknowns+unknown] += weight * fineGridValues[srcLinearised*unknowns+unknown];
112 assertion2( coarseGridValues[destLinearised*unknowns+unknown] == coarseGridValues[destLinearised*unknowns+unknown], srcLinearised, destLinearised );
113 }
114
115 logDebug(
116 "restrictInnerHalfOfHaloLayer_AoS_tensor_product(...)",
117 "add dof " << src << " (" << srcLinearised << ") to "
118 << dest << " (" << destLinearised << ") with weight " << weight << " for marker " << marker.toString() <<
119 ": " << coarseGridValues[destLinearised*unknowns+0] << "<-" << fineGridValues[srcLinearised*unknowns+0] );
120 }
121 }
122 logDebug( "restrictInnerHalfOfHaloLayer_AoS_tensor_product(...)", "update coarse volume/dof " << dest << ": " << coarseGridValues[destLinearised*5+0] );
123 }
124 }
125
126 logTraceOut( "restrictInnerHalfOfHaloLayer_AoS_tensor_product(...)" );
127}
128
131 int numberOfDoFsPerAxisInPatch,
132 int overlap,
133 int unknowns,
134 double* fineGridValues,
135 double* coarseGridValues
136) {
137 restrictInnerHalfOfHaloLayer_AoS_inject( marker, numberOfDoFsPerAxisInPatch, overlap, unknowns, fineGridValues, coarseGridValues, false );
138 restrictInnerHalfOfHaloLayer_AoS_inject( marker, numberOfDoFsPerAxisInPatch, overlap, unknowns, fineGridValues, coarseGridValues, true );
139}
140
142 int numberOfDoFsPerAxisInSourcePatch,
143 int numberOfDoFsPerAxisInDestinationPatch,
144 int unknowns,
145 double* sourceValues,
146 double* destinationValues
147) {
148 // how many fine grid cells are mapped onto one destination cell
149 double sourceToDestinationRatio = static_cast<double>(numberOfDoFsPerAxisInSourcePatch) / static_cast<double>(numberOfDoFsPerAxisInDestinationPatch);
150 dfor( destinationVolume, numberOfDoFsPerAxisInDestinationPatch ) {
151 const int baseIndexDestination = peano4::utils::dLinearised(destinationVolume,numberOfDoFsPerAxisInDestinationPatch) * unknowns;
152
154 for (int d=0; d<Dimensions; d++) {
155 sourceVolume(d) = std::floor(
156 static_cast<double>(destinationVolume(d)) * sourceToDestinationRatio + 0.5 * sourceToDestinationRatio
157 );
158 }
159
160 const int baseIndexSource = peano4::utils::dLinearised(sourceVolume,numberOfDoFsPerAxisInSourcePatch) * unknowns;
161 for (int unknown=0; unknown<unknowns; unknown++) {
162 destinationValues[baseIndexDestination+unknown] = sourceValues[baseIndexSource+unknown];
163 }
164 }
165}
166
168 int numberOfDoFsPerAxisInSourcePatch,
169 int numberOfDoFsPerAxisInDestinationPatch,
170 int unknowns,
171 double* sourceValues,
172 double* destinationValues,
173 double weightOfInjectedValue
174) {
175 assertion1( weightOfInjectedValue>=0.0, weightOfInjectedValue );
176 assertion1( weightOfInjectedValue<=1.0, weightOfInjectedValue );
177
178 // how many fine grid cells are mapped onto one destination cell
179 double sourceToDestinationRatio = static_cast<double>(numberOfDoFsPerAxisInSourcePatch) / static_cast<double>(numberOfDoFsPerAxisInDestinationPatch);
180 dfor( destinationVolume, numberOfDoFsPerAxisInDestinationPatch ) {
181 const int baseIndexDestination = peano4::utils::dLinearised(destinationVolume,numberOfDoFsPerAxisInDestinationPatch) * unknowns;
182
184 for (int d=0; d<Dimensions; d++) {
185 sourceVolume(d) = std::floor(
186 static_cast<double>(destinationVolume(d)) * sourceToDestinationRatio + 0.5 * sourceToDestinationRatio
187 );
188 }
189
190 const int baseIndexSource = peano4::utils::dLinearised(sourceVolume,numberOfDoFsPerAxisInSourcePatch) * unknowns;
191 for (int unknown=0; unknown<unknowns; unknown++) {
192 destinationValues[baseIndexDestination+unknown] =
193 weightOfInjectedValue * sourceValues[baseIndexSource+unknown]
194 + (1.0-weightOfInjectedValue) * destinationValues[baseIndexDestination+unknown];
195 }
196 }
197}
198
201 int numberOfDoFsPerAxisInPatch,
202 int unknowns,
203 double* fineGridValues,
204 double* coarseGridValues
205) {
206 dfor(fineVolume,numberOfDoFsPerAxisInPatch) {
207 tarch::la::Vector<Dimensions, int> fineVolumeWithintCoarsePatch = (fineVolume + marker.getRelativePositionWithinFatherCell() * numberOfDoFsPerAxisInPatch);
209
210 bool restrictThisVolume = true;
211 for (int d = 0; d < Dimensions; d++) {
212 restrictThisVolume &= (fineVolumeWithintCoarsePatch(d) % 3) == 1;
213 coarseVolume(d) = fineVolumeWithintCoarsePatch(3) / 3;
214 }
215
216 if (restrictThisVolume) {
217 int coarseVolumeLinearised = peano4::utils::dLinearised( coarseVolume, numberOfDoFsPerAxisInPatch );
218 int fineVolumeLinearised = peano4::utils::dLinearised( fineVolume, numberOfDoFsPerAxisInPatch );
219 for (int j=0; j<unknowns; j++) {
220 assertion3(coarseGridValues[coarseVolumeLinearised*unknowns+j]==coarseGridValues[coarseVolumeLinearised*unknowns+j], coarseVolume, fineVolume, marker.toString());
221 assertion3(fineGridValues[fineVolumeLinearised*unknowns+j]==fineGridValues[fineVolumeLinearised*unknowns+j], coarseVolume, fineVolume, marker.toString());
222 coarseGridValues[coarseVolumeLinearised*unknowns+j] = fineGridValues[fineVolumeLinearised*unknowns+j];
223 }
224 }
225 }
226}
227
230 int numberOfDoFsPerAxisInPatch,
231 int overlap,
232 int unknowns,
233 double* fineGridValues,
234 double* coarseGridValues
235) {
236 restrictInnerHalfOfHaloLayer_AoS_averaging( marker, numberOfDoFsPerAxisInPatch, overlap, unknowns, fineGridValues, coarseGridValues, false );
237 restrictInnerHalfOfHaloLayer_AoS_averaging( marker, numberOfDoFsPerAxisInPatch, overlap, unknowns, fineGridValues, coarseGridValues, true );
238}
239
242 int numberOfDoFsPerAxisInPatch,
243 int unknowns,
244 double* fineGridValues,
245 double* coarseGridValues
246) {
247 double scaleFineGridVolume = std::pow(3.0,-static_cast<double>(Dimensions));
248 internal::projectCells_AoS(
249 marker,
250 numberOfDoFsPerAxisInPatch,
251 [&](
252 [[maybe_unused]] tarch::la::Vector<Dimensions, int> coarseVolume,
253 [[maybe_unused]] tarch::la::Vector<Dimensions, int> fineVolume,
254 [[maybe_unused]] tarch::la::Vector<Dimensions, double> coarseVolumeCentre,
255 [[maybe_unused]] tarch::la::Vector<Dimensions, double> fineVolumeCentre,
256 [[maybe_unused]] double coarseVolumeH,
257 [[maybe_unused]] double fineVolumeH
258 ) -> void {
259 int coarseVolumeLinearised = peano4::utils::dLinearised( coarseVolume, numberOfDoFsPerAxisInPatch );
260 int fineVolumeLinearised = peano4::utils::dLinearised( fineVolume, numberOfDoFsPerAxisInPatch );
261 logDebug( "restrictCell_AoS_averaging(...)", fineVolume << " -> " << coarseVolume );
262 for (int j=0; j<unknowns; j++) {
263 coarseGridValues[coarseVolumeLinearised*unknowns+j] += scaleFineGridVolume * fineGridValues[fineVolumeLinearised*unknowns+j];
264 }
265 }
266 );
267}
268
271 int numberOfDoFsPerAxisInPatch,
272 int overlap,
273 int unknowns,
274 double* value
275) {
276 const int normal = marker.getSelectedFaceNumber() % Dimensions;
277
278 // clear fine grid data structure
279 dfore(k,numberOfDoFsPerAxisInPatch,normal,0) {
280 for (int i=0; i<overlap*2; i++) {
282 targetCell(normal) = i;
283 int targetCellSerialised = serialiseVoxelIndexInOverlap(
284 targetCell,
285 numberOfDoFsPerAxisInPatch,
286 overlap,
287 normal
288 );
289 assertion(targetCellSerialised>=0);
290 for (int j=0; j<unknowns; j++) {
291 value[targetCellSerialised*unknowns+j] = 0.0;
292 }
293 }
294 }
295}
296
298 [[maybe_unused]] const peano4::datamanagement::CellMarker& marker,
299 [[maybe_unused]] int numberOfDoFsPerAxisInPatch,
300 [[maybe_unused]] int unknowns,
301 [[maybe_unused]] double* values
302) {
303 #if Dimensions==3
304 for (int i=0; i<numberOfDoFsPerAxisInPatch*numberOfDoFsPerAxisInPatch*numberOfDoFsPerAxisInPatch*unknowns; i++) {
305 #else
306 for (int i=0; i<numberOfDoFsPerAxisInPatch*numberOfDoFsPerAxisInPatch*unknowns; i++) {
307 #endif
308 values[i] = 0.0;
309 }
310}
311
314 int numberOfDoFsPerAxisInPatch,
315 int overlap,
316 int unknowns,
317 double* fineGridValues,
318 double* coarseGridValues,
319 bool swapInsideOutside
320) {
321 logTraceInWith4Arguments( "restrictInnerHalfOfHaloLayer_AoS_averaging(...)", marker.toString(), numberOfDoFsPerAxisInPatch, overlap, unknowns );
322
323 assertion1( overlap==1 or overlap%3==0, overlap );
324
325 bool mapFromInnerHalfOfHalo = not swapInsideOutside; // mapOuterCoarseGridHaloOntoInnerFineGridHalo
326 const int normal = marker.getSelectedFaceNumber() % Dimensions;
327 double scaleFineGridVolume = std::pow(3.0,-static_cast<double>(Dimensions-1)) / std::min(3.0,static_cast<double>(overlap));
328 const bool pickLeftHalfOfHalo = (marker.getSelectedFaceNumber() < Dimensions) xor mapFromInnerHalfOfHalo;
329
330 const double volumeHCoarse = marker.h()(0) / static_cast<double>(numberOfDoFsPerAxisInPatch) * 3.0;
331 const double volumeHFine = marker.h()(0) / static_cast<double>(numberOfDoFsPerAxisInPatch);
332
333 tarch::la::Vector<Dimensions,double> leftBottomCornerOfHaloFine = marker.x();
334 tarch::la::Vector<Dimensions,double> leftBottomCornerOfHaloCoarse = marker.x();
335 for (int d=0; d<Dimensions; d++) {
336 if (d==normal) {
337 leftBottomCornerOfHaloFine(d) -= static_cast<double>(overlap) * volumeHFine;
338 leftBottomCornerOfHaloCoarse(d) -= static_cast<double>(overlap) * volumeHCoarse;
339 }
340 else {
341 leftBottomCornerOfHaloFine(d) = marker.x()(d)-marker.h()(d)/2.0;
342 leftBottomCornerOfHaloCoarse(d) = leftBottomCornerOfHaloFine(d)-marker.getRelativePositionWithinFatherCell()(d)*marker.h()(d);
343 }
344 }
345
346 dfore(kFine,numberOfDoFsPerAxisInPatch,normal,0) {
347 for (int iFine=0; iFine<overlap; iFine++) {
348 tarch::la::Vector<Dimensions,int> fineVolume = kFine;
349 fineVolume(normal) += pickLeftHalfOfHalo ? iFine : iFine + overlap;
350
351 tarch::la::Vector<Dimensions,int> coarseVolume = fineVolume;
352 for (int d=0; d<Dimensions; d++) {
353 if (d!=normal) {
354 coarseVolume(d) += marker.getRelativePositionWithinFatherCell()(d) * numberOfDoFsPerAxisInPatch;
355 }
356 else {
357 coarseVolume(d) = pickLeftHalfOfHalo ? fineVolume(d) : fineVolume(d) + overlap * 3;
358 }
359 }
360 coarseVolume /= 3;
361
362 assertion3( coarseVolume(0)>=0, kFine, iFine, marker );
363 assertion3( coarseVolume(1)>=0, kFine, iFine, marker );
364 assertion3( coarseVolume(0)<numberOfDoFsPerAxisInPatch, kFine, iFine, marker );
365 assertion3( coarseVolume(1)<numberOfDoFsPerAxisInPatch, kFine, iFine, marker );
366
367 [[maybe_unused]] tarch::la::Vector<Dimensions, double> volumeXCoarse = leftBottomCornerOfHaloCoarse
368 + volumeHCoarse * tarch::la::convertScalar<double>(coarseVolume)
369 + 0.5 * tarch::la::Vector<Dimensions, double>(volumeHCoarse);
370 [[maybe_unused]] tarch::la::Vector<Dimensions, double> volumeXFine = leftBottomCornerOfHaloFine
371 + volumeHFine * tarch::la::convertScalar<double>(fineVolume)
372 + 0.5 * tarch::la::Vector<Dimensions, double>(volumeHFine);
373
374 int coarseVolumeLinearised = serialiseVoxelIndexInOverlap(
375 coarseVolume,
376 numberOfDoFsPerAxisInPatch, overlap, normal
377 );
378 int fineVolumeLinearised = serialiseVoxelIndexInOverlap(
379 fineVolume,
380 numberOfDoFsPerAxisInPatch, overlap, normal
381 );
382 for (int j=0; j<unknowns; j++) {
383 coarseGridValues[coarseVolumeLinearised*unknowns+j] += scaleFineGridVolume * fineGridValues[fineVolumeLinearised*unknowns+j];
384 }
385 }
386 }
387
388 logTraceOut( "restrictInnerHalfOfHaloLayer_AoS_averaging(...)" );
389}
390
393 int numberOfDoFsPerAxisInPatch,
394 int overlap,
395 int unknowns,
396 double* fineGridValues,
397 double* coarseGridValues,
398 bool swapInsideOutside
399) {
400 logTraceInWith6Arguments( "restrictInnerHalfOfHaloLayer_AoS_inject(...)", marker.toString(), numberOfDoFsPerAxisInPatch, overlap, unknowns, fineGridValues, coarseGridValues );
401
402 assertion1( overlap==1 or overlap%3==1, overlap );
403
404 const int normal = marker.getSelectedFaceNumber() % Dimensions;
405 const bool pickLeftHalfOfHaloOnFineGrid = (marker.getSelectedFaceNumber() >= Dimensions) xor swapInsideOutside;
406
407 dfore(kFine,numberOfDoFsPerAxisInPatch,normal,0) {
408 tarch::la::Vector<Dimensions,int> fineVolume = kFine;
409 tarch::la::Vector<Dimensions,int> fineVolumeAlongCoarseFace = kFine;
410
411 for (int iFine=0; iFine<overlap; iFine++) {
412 fineVolume(normal) = iFine;
413
414 bool restrictThisVolume = true;
415 for (int d=0; d<Dimensions; d++) {
416 if (d==normal) {
417 restrictThisVolume &= (overlap==1 or fineVolumeAlongCoarseFace(d)%3==1);
418 }
419 else {
420 fineVolumeAlongCoarseFace(d)+=marker.getRelativePositionWithinFatherCell()(d)*numberOfDoFsPerAxisInPatch;
421 restrictThisVolume &= (fineVolumeAlongCoarseFace(d)%3)==1;
422 }
423 }
424
425 if (restrictThisVolume) {
426 tarch::la::Vector<Dimensions,int> coarseVolume = fineVolumeAlongCoarseFace/3;
427
428 fineVolume(normal) = pickLeftHalfOfHaloOnFineGrid ? iFine : iFine + overlap;
429 coarseVolume(normal) = pickLeftHalfOfHaloOnFineGrid ? iFine : iFine + overlap;
430
431 int coarseVolumeLinearised = serialiseVoxelIndexInOverlap(
432 coarseVolume,
433 numberOfDoFsPerAxisInPatch, overlap, normal
434 );
435 int fineVolumeLinearised = serialiseVoxelIndexInOverlap(
436 fineVolume,
437 numberOfDoFsPerAxisInPatch, overlap, normal
438 );
439 for (int j=0; j<unknowns; j++) {
440 assertion3(coarseGridValues[coarseVolumeLinearised*unknowns+j]==coarseGridValues[coarseVolumeLinearised*unknowns+j], coarseVolume, fineVolume, marker.toString());
441 assertion3(fineGridValues[fineVolumeLinearised*unknowns+j]==fineGridValues[fineVolumeLinearised*unknowns+j], coarseVolume, fineVolume, marker.toString());
442 coarseGridValues[coarseVolumeLinearised*unknowns+j] = fineGridValues[fineVolumeLinearised*unknowns+j];
443 }
444 } // end of restrictThisVolume
445 } // for iFine loop, i.e. loop over halo layer depth
446 } // dfore loop, i.e. loop over submanifold
447
448 logTraceOut( "restrictInnerHalfOfHaloLayer_AoS_restrict(...)" );
449}
450
453 int numberOfDoFsPerAxisInPatch,
454 int overlap,
455 int unknowns,
456 bool clearInnerPart,
457 double* value
458) {
459 const int normal = marker.getSelectedFaceNumber() % Dimensions;
460 const bool left = (marker.getSelectedFaceNumber() < Dimensions) xor clearInnerPart;
461
462 // clear fine grid data structure
463 dfore(k,numberOfDoFsPerAxisInPatch,normal,0) {
464 for (int i=0; i<overlap; i++) {
466 targetCell(normal) = left ? i : i + overlap;
467 int targetCellSerialised = serialiseVoxelIndexInOverlap(
468 targetCell,
469 numberOfDoFsPerAxisInPatch,
470 overlap,
471 normal
472 );
473 assertion(targetCellSerialised>=0);
474 for (int j=0; j<unknowns; j++) {
475 value[targetCellSerialised*unknowns+j] = 0.0;
476 }
477 }
478 }
479}
480
483 int numberOfDoFsPerAxisInPatch,
484 std::function<void(
487 tarch::la::Vector<Dimensions,double> coarseVolumeCentre,
489 double coarseVolumeH,
490 double fineVolumeH
491 )> update
492) {
493 logTraceInWith2Arguments( "projectCells_AoS(...)", marker.toString(), numberOfDoFsPerAxisInPatch );
494
495 const double volumeHCoarse = marker.h()(0) / static_cast<double>(numberOfDoFsPerAxisInPatch) * 3.0;
496 const double volumeHFine = marker.h()(0) / static_cast<double>(numberOfDoFsPerAxisInPatch);
497
498 tarch::la::Vector<Dimensions,double> leftBottomOfFineCell = marker.getOffset();
499 tarch::la::Vector<Dimensions,double> leftBottomOfCoarseCell = marker.getOffset() - tarch::la::multiplyComponents( marker.h(), tarch::la::convertScalar<double>(marker.getRelativePositionWithinFatherCell()) );
500
501 dfor(kFine,numberOfDoFsPerAxisInPatch) {
502 tarch::la::Vector<Dimensions,int> kCoarse = (kFine + numberOfDoFsPerAxisInPatch * marker.getRelativePositionWithinFatherCell()) / 3;
503
504 tarch::la::Vector<Dimensions,double> volumeXCoarse = leftBottomOfCoarseCell
505 + volumeHCoarse * tarch::la::convertScalar<double>(kCoarse)
506 + 0.5 * tarch::la::Vector<Dimensions,double>(volumeHCoarse);
507 tarch::la::Vector<Dimensions,double> volumeXFine = leftBottomOfFineCell
508 + volumeHFine * tarch::la::convertScalar<double>(kFine)
509 + 0.5 * tarch::la::Vector<Dimensions,double>(volumeHFine);
510
511 update(
512 kCoarse,
513 kFine,
514 volumeXCoarse,
515 volumeXFine,
516 volumeHCoarse,
517 volumeHFine
518 );
519 }
520
521 logTraceOut( "projectCells_AoS(...)" );
522}
#define assertion2(expr, param0, param1)
#define assertion4(expr, param0, param1, param2, param3)
#define assertion3(expr, param0, param1, param2)
#define assertion1(expr, param)
#define assertionMsg(expr, message)
#define assertion(expr)
#define assertion5(expr, param0, param1, param2, param3, param4)
And from this we can write down f$ nabla phi_i nabla phi_i dx but since we are constructing matrix let s investigate the f$ our matrix elements will nabla phi_i dx f By this will be a sparse as these basis functions are chosen to not overlap with each other almost everywhere In other they have only local support We can read off the right hand side values
And from this we can write down f$ nabla phi_i nabla phi_i dx but since we are constructing matrix let s investigate the f$ j
we integrate over each cell and then take the sum across each of the cells We also consider the terms that enter the f$ k
#define logDebug(methodName, logMacroMessageStream)
Definition Log.h:50
#define logTraceOut(methodName)
Definition Log.h:379
#define logTraceInWith4Arguments(methodName, argument0, argument1, argument2, argument3)
Definition Log.h:373
#define logTraceInWith3Arguments(methodName, argument0, argument1, argument2)
Definition Log.h:372
#define logTraceInWith6Arguments(methodName, argument0, argument1, argument2, argument3, argument4, argument5)
Definition Log.h:375
#define logTraceInWith2Arguments(methodName, argument0, argument1)
Definition Log.h:371
#define dfore(counter, max, dim, value)
This is an exclusive d-dimensional for loop.
Definition Loop.h:937
#define dfor(counter, max)
d-dimensional Loop
Definition Loop.h:308
Log Device.
Definition Log.h:516
CPUGPUMethod int dLinearised(const tarch::la::Vector< Dimensions, int > &counter, int max)
Map d-dimensional vector onto integer index.
Definition Loop.cpp:106
Matrix< Rows, Cols, Scalar > multiplyComponents(const Matrix< Rows, X, Scalar > &lhs, const Matrix< X, Cols, Scalar > &rhs)
void clearHalfOfHaloLayerAoS(const peano4::datamanagement::FaceMarker &marker, int numberOfDoFsPerAxisInPatch, int overlap, int unknowns, bool clearInnerPart, double *values)
Clear half of a halo layer.
void projectCells_AoS(const peano4::datamanagement::CellMarker &fineGridCellMarker, int numberOfDoFsPerAxisInPatch, std::function< void(tarch::la::Vector< Dimensions, int > coarseVolume, tarch::la::Vector< Dimensions, int > fineVolume, tarch::la::Vector< Dimensions, double > coarseVolumeCentre, tarch::la::Vector< Dimensions, double > fineVolumeCentre, double coarseVolumeH, double fineVolumeH)> update)
Helper function.
void restrictHaloLayer_AoS_averaging(const peano4::datamanagement::FaceMarker &marker, int numberOfDoFsPerAxisInPatch, int overlap, int unknowns, double *fineGridValues, double *coarseGridValues)
Consult commend on interpolation that clarifies why we need two different halo layer restrictions,...
void restrictCellIntoOverlappingCell_inject_and_average(int numberOfDoFsPerAxisInSourcePatch, int numberOfDoFsPerAxisInDestinationPatch, int unknowns, double *sourceValues, double *destinationValues, double weightOfInjectedValue=0.5)
Flavour of restrictCellIntoOverlappingCell_inject() where we inject the solution but then take the av...
void restrictCellIntoOverlappingCell_inject(int numberOfDoFsPerAxisInSourcePatch, int numberOfDoFsPerAxisInDestinationPatch, int unknowns, double *sourceValues, double *destinationValues)
This routine should be used if a cell hosts two sets of unknowns.
void restrictHaloLayer_AoS_inject(const peano4::datamanagement::FaceMarker &marker, int numberOfDoFsPerAxisInPatch, int overlap, int unknowns, double *fineGridValues, double *coarseGridValues)
Restrict data by injection.
void clearCell(const peano4::datamanagement::CellMarker &marker, int numberOfDoFsPerAxisInPatch, int unknowns, double *values)
void clearHaloLayerAoS(const peano4::datamanagement::FaceMarker &marker, int numberOfDoFsPerAxisInPatch, int overlap, int unknowns, double *values)
Clear halo layer of face.
void restrictInnerHalfOfHaloLayer_AoS_inject(const peano4::datamanagement::FaceMarker &marker, int numberOfDoFsPerAxisInPatch, int overlap, int unknowns, double *fineGridValues, double *coarseGridValues, bool swapInsideOutside=false)
Restrict data by injection.
void restrictInnerHalfOfHaloLayer_AoS_averaging(const peano4::datamanagement::FaceMarker &marker, int numberOfDoFsPerAxisInPatch, int overlap, int unknowns, double *fineGridValues, double *coarseGridValues, bool swapInsideOutside=false)
Restrict with averaging.
void restrictHaloLayer_AoS_tensor_product(const peano4::datamanagement::FaceMarker &marker, int numberOfDoFsPerAxisInPatch, int overlap, int unknowns, double *fineGridValues, double *coarseGridValues)
void restrictCell_AoS_averaging(const peano4::datamanagement::CellMarker &marker, int numberOfDoFsPerAxisInPatch, int unknowns, double *fineGridValues, double *coarseGridValues)
This routine is used when we delete a cell due to dynamic AMR.
int serialiseVoxelIndexInOverlap(const tarch::la::Vector< Dimensions, int > &overlapCell, int numberOfDoFsPerAxisInPatch, int overlap, int normal)
The volumes or elements within an overlap are always enumerated lexicographically.
void restrictCell_AoS_inject(const peano4::datamanagement::CellMarker &marker, int numberOfDoFsPerAxisInPatch, int unknowns, double *fineGridValues, double *coarseGridValues)
void restrictCell_AoS_tensor_product(const peano4::datamanagement::CellMarker &marker, int numberOfDoFsPerAxisInPatch, int unknowns, double *fineGridValues, double *coarseGridValues)
void restrictInnerHalfOfHaloLayer_AoS_tensor_product(const peano4::datamanagement::FaceMarker &marker, int numberOfDoFsPerAxisInPatch, int overlap, int unknowns, double *fineGridValues, double *coarseGridValues, bool swapInsideOutside=false)
tarch::logging::Log _log("examples::unittests")
Provide information about selected face.
Definition FaceMarker.h:35
Simple vector class.
Definition Vector.h:134