Skip to content

Commit 62b32c6

Browse files
tobre1FilipovicLadokenyastyle
authored
fix: intermediate velocity calculation in higher order time schemes (#211)
* Added removeMaterial info when removing top LS layer * Added WENO 5th order advection scheme * Reverted back to install script, with ViennaLS default version 5.2.1 * Add option to enable adaptive time stepping during advection * Add adptive time step threshold * Update ViennaLS version * Implement VTK rendering * Add python bindings * Use caching for switch rendering modes * Screenshot feature * renaming integration -> discretization * Domain offsets * Add more features to render window * Add material colors and scalar bar * Hex to RGB util * Fix sphere dist examples * Add comments * Format * Update README * Update advection handler to use RK3 time integration from ViennaLS. Tested with the multiParticleProcess test, holeEtching example shows how to select a temporal scheme. * format (newline at the end of CMakeLists.txt) * Fixed format. Binding errors persist due to ViennaLS version missmatch.. This PR must follow the ViennaLS PR. * Adapt ViennaPS to match temporal scheme from the latest ViennaLS commit. * Allow for calculating intermediate velocities during higher order time integration * Introduce velocity recalculation during intermediate time integration stpes when using higher order schemes (e.g., RK2, RK3). Added python bindings and a simpleEthing example which uses these schemes in C++ and Python. * adapted examples, tested against latest ViennaLS build for int-velocities * fix for removeStrayPoints * revert /python/viennaps/__init__.py * format * remove GIT_TAG * minor fixes * Set CMake to use ViennaLS5.5.0 for intermediate velocity calculation for RK2/3 * fix holeEtching config * format * format * avoid caching FetchContent deps on Windows * Revert to caching windows-latest * Clean FetchContent deps (Windows) * Add intermediate test and update install_ViennaPS.py * fix _core * format - empty line at end of intermediate.cpp * test seg fault * Update intermediate test to use 1st order LF * add advection steps and flux calculation num debug output * Reset velocity update callback in AdvectionHandler & improve intermediate test * remove test folder * refactor: small changes * refactor: reduce intermediate test complexity * chore: bump version --------- Co-authored-by: filipovic <filipovic@iue.tuwien.ac.at> Co-authored-by: Roman Kostal <roman.kstl@gmail.com>
1 parent ffe8d91 commit 62b32c6

25 files changed

+1381
-1089
lines changed

‎CMakeLists.txt‎

Lines changed: 3 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -2,7 +2,7 @@ cmake_minimum_required(VERSION 3.20 FATAL_ERROR)
22
project(
33
ViennaPS
44
LANGUAGES CXX C
5-
VERSION 4.2.1)
5+
VERSION 4.2.2)
66

77
# --------------------------------------------------------------------------------------------------------
88
# Include guards
@@ -119,14 +119,14 @@ CPMAddPackage(
119119

120120
CPMAddPackage(
121121
NAME ViennaRay
122-
VERSION 3.11.0
122+
VERSION 3.11.1
123123
GIT_REPOSITORY "https://github.com/ViennaTools/ViennaRay"
124124
EXCLUDE_FROM_ALL ${VIENNAPS_BUILD_PYTHON}
125125
OPTIONS "VIENNARAY_USE_GPU ${VIENNAPS_USE_GPU}")
126126

127127
CPMAddPackage(
128128
NAME ViennaLS
129-
VERSION 5.5.0
129+
VERSION 5.5.1
130130
GIT_REPOSITORY "https://github.com/ViennaTools/ViennaLS"
131131
EXCLUDE_FROM_ALL ${VIENNAPS_BUILD_PYTHON}
132132
OPTIONS "VIENNALS_PRECOMPILE_HEADERS ${VIENNAPS_PRECOMPILE_HEADERS}")

‎README.md‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -130,7 +130,7 @@ We recommend using [CPM.cmake](https://github.com/cpm-cmake/CPM.cmake) to consum
130130

131131
* Installation with CPM
132132
```cmake
133-
CPMAddPackage("gh:viennatools/viennaps@4.2.1")
133+
CPMAddPackage("gh:viennatools/viennaps@4.2.2")
134134
```
135135

136136
* With a local installation

‎include/viennaps/process/psAdvectionHandler.hpp‎

Lines changed: 12 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -11,7 +11,8 @@ namespace viennaps {
1111
VIENNAPS_TEMPLATE_ND(NumericType, D) class AdvectionHandler {
1212
viennals::Advect<NumericType, D> advectionKernel_;
1313
viennacore::Timer<> timer_;
14-
unsigned lsVelOutputCounter = 0;
14+
unsigned lsVelOutputCounter_ = 0;
15+
unsigned totalAdvectionSteps_ = 0;
1516

1617
public:
1718
ProcessResult initialize(ProcessContext<NumericType, D> &context) {
@@ -39,11 +40,10 @@ VIENNAPS_TEMPLATE_ND(NumericType, D) class AdvectionHandler {
3940

4041
context.resetTime();
4142

42-
advectionKernel_.setTemporalScheme(context.advectionParams.temporalScheme);
43-
4443
advectionKernel_.setSingleStep(true);
45-
advectionKernel_.setVelocityField(context.translationField);
4644
advectionKernel_.setSpatialScheme(context.advectionParams.spatialScheme);
45+
advectionKernel_.setTemporalScheme(context.advectionParams.temporalScheme);
46+
advectionKernel_.setVelocityField(context.translationField);
4747
advectionKernel_.setTimeStepRatio(context.advectionParams.timeStepRatio);
4848
advectionKernel_.setSaveAdvectionVelocities(
4949
context.advectionParams.velocityOutput);
@@ -56,6 +56,8 @@ VIENNAPS_TEMPLATE_ND(NumericType, D) class AdvectionHandler {
5656
context.advectionParams.adaptiveTimeStepping,
5757
context.advectionParams.adaptiveTimeStepSubdivisions);
5858

59+
advectionKernel_.setVelocityUpdateCallback(nullptr);
60+
5961
// normals vectors are only necessary for analytical velocity fields
6062
if (translationMethod > 0)
6163
advectionKernel_.setCalculateNormalVectors(false);
@@ -65,6 +67,8 @@ VIENNAPS_TEMPLATE_ND(NumericType, D) class AdvectionHandler {
6567
advectionKernel_.insertNextLevelSet(dom);
6668
}
6769

70+
totalAdvectionSteps_ = 0;
71+
6872
return ProcessResult::SUCCESS;
6973
}
7074

@@ -78,6 +82,8 @@ VIENNAPS_TEMPLATE_ND(NumericType, D) class AdvectionHandler {
7882
advectionKernel_.setVelocityUpdateCallback(callback);
7983
}
8084

85+
auto getTotalAdvectionSteps() const { return totalAdvectionSteps_; }
86+
8187
void disableSingleStep() { advectionKernel_.setSingleStep(false); }
8288

8389
void prepareAdvection(const ProcessContext<NumericType, D> &context) {
@@ -99,13 +105,14 @@ VIENNAPS_TEMPLATE_ND(NumericType, D) class AdvectionHandler {
99105
advectionKernel_.apply();
100106
timer_.finish();
101107

108+
++totalAdvectionSteps_;
102109
if (context.advectionParams.velocityOutput) {
103110
auto mesh = viennals::Mesh<NumericType>::New();
104111
viennals::ToMesh<NumericType, D>(context.domain->getSurface(), mesh)
105112
.apply();
106113
viennals::VTKWriter<NumericType>(
107114
mesh,
108-
"ls_velocities_" + std::to_string(lsVelOutputCounter++) + ".vtp")
115+
"ls_velocities_" + std::to_string(lsVelOutputCounter_++) + ".vtp")
109116
.apply();
110117
}
111118

‎include/viennaps/process/psCPUDiskEngine.hpp‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -141,6 +141,7 @@ class CPUDiskEngine final : public FluxEngine<NumericType, D> {
141141
}
142142
rayTracer_.setParticleType(particle);
143143
rayTracer_.apply();
144+
++this->fluxCalculationsCount_;
144145

145146
auto info = rayTracer_.getRayTraceInfo();
146147

‎include/viennaps/process/psCPUTriangleEngine.hpp‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -230,6 +230,7 @@ class CPUTriangleEngine final : public FluxEngine<NumericType, D> {
230230
}
231231
rayTracer_.setParticleType(particle);
232232
rayTracer_.apply();
233+
++this->fluxCalculationsCount_;
233234

234235
auto info = rayTracer_.getRayTraceInfo();
235236

‎include/viennaps/process/psFluxEngine.hpp‎

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -9,6 +9,7 @@ namespace viennaps {
99
VIENNAPS_TEMPLATE_ND(NumericType, D) class FluxEngine {
1010
protected:
1111
viennacore::Timer<> timer_;
12+
unsigned fluxCalculationsCount_ = 0;
1213

1314
public:
1415
virtual ~FluxEngine() = default;
@@ -28,6 +29,7 @@ VIENNAPS_TEMPLATE_ND(NumericType, D) class FluxEngine {
2829

2930
auto &getTimer() const { return timer_; }
3031
void resetTimer() { timer_.reset(); }
32+
auto getFluxCalculationsCount() const { return fluxCalculationsCount_; }
3133
};
3234

3335
} // namespace viennaps

‎include/viennaps/process/psFluxProcessStrategy.hpp‎

Lines changed: 45 additions & 40 deletions
Original file line numberDiff line numberDiff line change
@@ -137,8 +137,7 @@ class FluxProcessStrategy final : public ProcessStrategy<NumericType, D> {
137137
translator_ = SmartPointer<TranslatorType>::New();
138138
meshGenerator_.setTranslator(translator_);
139139
context.translationField->setTranslator(translator_);
140-
}
141-
if (translationMethod == 2) {
140+
} else if (translationMethod == 2) {
142141
if (!kdTree_)
143142
kdTree_ = SmartPointer<KDTree<NumericType, Vec3D<NumericType>>>::New();
144143
context.translationField->setKdTree(kdTree_);
@@ -243,7 +242,19 @@ class FluxProcessStrategy final : public ProcessStrategy<NumericType, D> {
243242
context.model->finalize(context.domain, context.processTime);
244243

245244
processTimer.finish();
246-
logProcessingTimes(context, processTimer);
245+
246+
if (Logger::hasTiming()) {
247+
logProcessingTimes(context, processTimer);
248+
}
249+
250+
if (Logger::hasDebug()) {
251+
auto numAdvectionSteps = advectionHandler_.getTotalAdvectionSteps();
252+
VIENNACORE_LOG_DEBUG("Total advection steps: " +
253+
std::to_string(numAdvectionSteps));
254+
auto numFluxCalculations = fluxEngine_->getFluxCalculationsCount();
255+
VIENNACORE_LOG_DEBUG("Total flux calculations: " +
256+
std::to_string(numFluxCalculations));
257+
}
247258

248259
if (static_cast<int>(context.domain->getMetaDataLevel()) > 1) {
249260
context.domain->addMetaData("ProcessTime", context.processTime);
@@ -268,7 +279,7 @@ class FluxProcessStrategy final : public ProcessStrategy<NumericType, D> {
268279
PROCESS_CHECK(fluxEngine_->updateSurface(context));
269280

270281
// Calculate fluxes
271-
auto fluxes = SmartPointer<viennals::PointData<NumericType>>::New();
282+
auto fluxes = viennals::PointData<NumericType>::New();
272283
PROCESS_CHECK(fluxEngine_->calculateFluxes(context, fluxes));
273284

274285
// Update coverages if needed
@@ -294,7 +305,9 @@ class FluxProcessStrategy final : public ProcessStrategy<NumericType, D> {
294305
advectionHandler_.copyCoveragesToLevelSet(context, translator_));
295306
}
296307

297-
outputIntermediateResults(context, velocities, fluxes);
308+
if (Logger::hasIntermediate()) {
309+
outputIntermediateResults(context, velocities, fluxes);
310+
}
298311

299312
// Perform advection, updates processTime, reduces level set to width 1
300313
PROCESS_CHECK(advectionHandler_.performAdvection(context));
@@ -324,35 +337,6 @@ class FluxProcessStrategy final : public ProcessStrategy<NumericType, D> {
324337
return ProcessResult::SUCCESS;
325338
}
326339

327-
void outputIntermediateResults(
328-
ProcessContext<NumericType, D> &context,
329-
SmartPointer<std::vector<NumericType>> &velocities,
330-
const SmartPointer<viennals::PointData<NumericType>> &fluxes) {
331-
if (Logger::hasIntermediate()) {
332-
auto const name = context.getProcessName();
333-
auto surfaceModel = context.model->getSurfaceModel();
334-
context.diskMesh->getCellData().insertNextScalarData(*velocities,
335-
"velocities");
336-
if (context.flags.useCoverages) {
337-
mergeScalarData(context.diskMesh->getCellData(),
338-
surfaceModel->getCoverages());
339-
}
340-
if (auto surfaceData = surfaceModel->getSurfaceData())
341-
mergeScalarData(context.diskMesh->getCellData(), surfaceData);
342-
mergeScalarData(context.diskMesh->getCellData(), fluxes);
343-
viennals::VTKWriter<NumericType>(
344-
context.diskMesh, context.intermediateOutputPath + name + "_" +
345-
std::to_string(context.currentIteration) +
346-
".vtp")
347-
.apply();
348-
if (context.domain->getCellSet()) {
349-
context.domain->getCellSet()->writeVTU(
350-
context.intermediateOutputPath + name + "_cellSet_" +
351-
std::to_string(context.currentIteration) + ".vtu");
352-
}
353-
}
354-
}
355-
356340
bool applyPreAdvectionCallback(ProcessContext<NumericType, D> &context) {
357341
callbackTimer_.start();
358342
bool result = context.model->getAdvectionCallback()->applyPreAdvect(
@@ -398,9 +382,7 @@ class FluxProcessStrategy final : public ProcessStrategy<NumericType, D> {
398382
void updateState(ProcessContext<NumericType, D> &context) {
399383
meshGenerator_.apply();
400384

401-
auto const translationMethod =
402-
context.translationField->getTranslationMethod();
403-
if (translationMethod == 2) {
385+
if (context.translationField->getTranslationMethod() == 2) {
404386
kdTree_->setPoints(context.diskMesh->getNodes());
405387
kdTree_->build();
406388
}
@@ -471,11 +453,34 @@ class FluxProcessStrategy final : public ProcessStrategy<NumericType, D> {
471453
}
472454
}
473455

456+
void outputIntermediateResults(
457+
const ProcessContext<NumericType, D> &context,
458+
const SmartPointer<std::vector<NumericType>> &velocities,
459+
const SmartPointer<viennals::PointData<NumericType>> &fluxes) {
460+
auto const name = context.getProcessName();
461+
auto surfaceModel = context.model->getSurfaceModel();
462+
context.diskMesh->getCellData().insertNextScalarData(*velocities,
463+
"velocities");
464+
if (context.flags.useCoverages) {
465+
mergeScalarData(context.diskMesh->getCellData(),
466+
surfaceModel->getCoverages());
467+
}
468+
if (auto surfaceData = surfaceModel->getSurfaceData())
469+
mergeScalarData(context.diskMesh->getCellData(), surfaceData);
470+
mergeScalarData(context.diskMesh->getCellData(), fluxes);
471+
viennals::VTKWriter<NumericType>(
472+
context.diskMesh, context.intermediateOutputPath + name + "_" +
473+
std::to_string(context.currentIteration) + ".vtp")
474+
.apply();
475+
if (context.domain->getCellSet()) {
476+
context.domain->getCellSet()->writeVTU(
477+
context.intermediateOutputPath + name + "_cellSet_" +
478+
std::to_string(context.currentIteration) + ".vtu");
479+
}
480+
}
481+
474482
void logProcessingTimes(const ProcessContext<NumericType, D> &context,
475483
const viennacore::Timer<> &processTimer) {
476-
if (!Logger::hasTiming())
477-
return;
478-
479484
Logger::getInstance()
480485
.addTiming("\nProcess " + context.getProcessName(),
481486
processTimer.currentDuration * 1e-9)

‎include/viennaps/process/psGPUDiskEngine.hpp‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -166,6 +166,7 @@ class GPUDiskEngine final : public FluxEngine<NumericType, D> {
166166
rayTracer_.normalizeResults();
167167
downloadResultsToPointData(*fluxes,
168168
context.rayTracingParams.smoothingNeighbors);
169+
++this->fluxCalculationsCount_;
169170

170171
// output
171172
if (Logger::hasIntermediate()) {

‎include/viennaps/process/psGPULineEngine.hpp‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -217,6 +217,7 @@ class GPULineEngine final : public FluxEngine<NumericType, D> {
217217
downloadResultsToPointData(*fluxes, context.diskMesh,
218218
context.rayTracingParams.smoothingNeighbors,
219219
context.domain->getGridDelta());
220+
++this->fluxCalculationsCount_;
220221

221222
// output
222223
if (Logger::hasIntermediate()) {

‎include/viennaps/process/psGPUTriangleEngine.hpp‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -209,6 +209,7 @@ class GPUTriangleEngine final : public FluxEngine<NumericType, D> {
209209

210210
// run the ray tracer
211211
rayTracer_.apply(); // device detach point here
212+
++this->fluxCalculationsCount_;
212213

213214
// Prepare post-processing
214215
PostProcessingType postProcessing(

0 commit comments

Comments
 (0)