Enable intermediate velocity calculation for RK2/RK3

.github/workflows/build.yml

@@ -54,9 +54,15 @@ jobs:
       - name: 🦥 Cache Dependencies
         uses: actions/cache@v4
         with:
-          key: test-${{ matrix.os }}-${{ matrix.config }}
+          key: test-${{ matrix.os }}-${{ matrix.config }}-${{ hashFiles('CMakeLists.txt', 'cmake/**', 'pyproject.toml') }}
           path: build
 
+      - name: 🧹 Clean FetchContent deps (Windows)
+        if: ${{ matrix.os == 'windows-latest' }}
+        shell: pwsh
+        run: |
+          if (Test-Path build/_deps) { Remove-Item -Recurse -Force build/_deps }
+
       - name: 🏗️ Compile
         run: cmake -DVIENNAPS_BUILD_TESTS=ON -B build && cmake --build build --config ${{ matrix.config }}
 

.github/workflows/python.yml

@@ -51,9 +51,15 @@ jobs:
       - name: 🦥 Cache Dependencies
         uses: actions/cache@v4
         with:
-          key: python-${{ matrix.os }}
+          key: python-${{ matrix.os }}-${{ hashFiles('CMakeLists.txt', 'cmake/**', 'pyproject.toml') }}
           path: build
 
+      - name: 🧹 Clean FetchContent deps (Windows)
+        if: ${{ matrix.os == 'windows-latest' }}
+        shell: pwsh
+        run: |
+          if (Test-Path build/_deps) { Remove-Item -Recurse -Force build/_deps }
+
       - name: 🛠️ Disable IPO (Alpine)
         if: ${{ matrix.os == 'ubuntu-latest' }}
         run: |
@@ -164,4 +170,4 @@ jobs:
           merge-multiple: true
 
       - name: 🚀 Publish Wheels
-        uses: pypa/gh-action-pypi-publish@release/v1
+        uses: pypa/gh-action-pypi-publish@release/v1

CMakeLists.txt

@@ -105,7 +105,7 @@ include("cmake/cpm.cmake")
 
 CPMAddPackage(
   NAME ViennaCore
-  VERSION 1.9.2
+  VERSION 1.9.4
   GIT_REPOSITORY "https://github.com/ViennaTools/ViennaCore"
   EXCLUDE_FROM_ALL ${VIENNAPS_BUILD_PYTHON}
   OPTIONS "VIENNACORE_USE_GPU ${VIENNAPS_USE_GPU}")
@@ -124,7 +124,7 @@ CPMAddPackage(
 
 CPMAddPackage(
   NAME ViennaLS
-  VERSION 5.4.0
+  VERSION 5.5.0
   GIT_REPOSITORY "https://github.com/ViennaTools/ViennaLS"
   EXCLUDE_FROM_ALL ${VIENNAPS_BUILD_PYTHON})
 

examples/holeEtching/config.txt

@@ -30,5 +30,7 @@ A_Si=7         # Si yield constant
 etchStopDepth=-10   # maximum etching depth
 spatialScheme=EO_1
 temporalScheme=RK3
+fluxEngine=CD
 
 raysPerPoint=1000
+outputFile=final.vtp

examples/holeEtching/holeEtching.cpp

@@ -8,9 +8,10 @@ using namespace viennaps;
 
 int main(int argc, char *argv[]) {
   using NumericType = double;
-  constexpr int D = 3;
+  constexpr int D = 2;
 
-  omp_set_num_threads(16);
+  Logger::setLogLevel(LogLevel::ERROR);
+  omp_set_num_threads(8);
 
   // Parse the parameters
   util::Parameters params;
@@ -30,54 +31,77 @@ int main(int argc, char *argv[]) {
   units::Length::setUnit(params.get<std::string>("lengthUnit"));
   units::Time::setUnit(params.get<std::string>("timeUnit"));
 
-  // geometry setup
-  auto geometry = Domain<NumericType, D>::New(
-      params.get("gridDelta"), params.get("xExtent"), params.get("yExtent"));
-  MakeHole<NumericType, D>(geometry, params.get("holeRadius"),
-                           0.0, // holeDepth
-                           0.0, // holeTaperAngle
-                           params.get("maskHeight"), params.get("taperAngle"),
-                           HoleShape::HALF)
-      .apply();
-
-  // use pre-defined model SF6O2 etching model
-  auto modelParams = SF6O2Etching<NumericType, D>::defaultParameters();
-  modelParams.ionFlux = params.get("ionFlux");
-  modelParams.etchantFlux = params.get("etchantFlux");
-  modelParams.passivationFlux = params.get("oxygenFlux");
-  modelParams.Ions.meanEnergy = params.get("meanEnergy");
-  modelParams.Ions.sigmaEnergy = params.get("sigmaEnergy");
-  modelParams.Ions.exponent = params.get("ionExponent");
-  modelParams.Passivation.A_ie = params.get("A_O");
-  modelParams.Substrate.A_ie = params.get("A_Si");
-  modelParams.etchStopDepth = params.get("etchStopDepth");
-  auto model = SmartPointer<SF6O2Etching<NumericType, D>>::New(modelParams);
-
-  CoverageParameters coverageParams;
-  coverageParams.tolerance = 1e-4;
-
-  RayTracingParameters rayTracingParams;
-  rayTracingParams.raysPerPoint = params.get<unsigned>("raysPerPoint");
-
-  AdvectionParameters advectionParams;
-  advectionParams.spatialScheme =
-      util::convertSpatialScheme(params.get<std::string>("spatialScheme"));
-  advectionParams.temporalScheme =
-      util::convertTemporalScheme(params.get<std::string>("temporalScheme"));
-
-  // process setup
-  Process<NumericType, D> process(geometry, model);
-  process.setProcessDuration(params.get("processTime"));
-  process.setParameters(coverageParams);
-  process.setParameters(rayTracingParams);
-  process.setParameters(advectionParams);
-
-  // print initial surface
-  geometry->saveSurfaceMesh("initial.vtp");
-
-  // run the process
-  process.apply();
-
-  // print final surface
-  geometry->saveSurfaceMesh("final.vtp", true, 0.01, true);
+  auto runSimulation = [&](bool intermediateVelocities, std::string suffix) {
+    // geometry setup
+    auto geometry = Domain<NumericType, D>::New(
+        params.get("gridDelta"), params.get("xExtent"), params.get("yExtent"));
+    MakeHole<NumericType, D>(geometry, params.get("holeRadius"),
+                             0.0, // holeDepth
+                             0.0, // holeTaperAngle
+                             params.get("maskHeight"), params.get("taperAngle"),
+                             HoleShape::QUARTER)
+        .apply();
+
+    // use pre-defined model SF6O2 etching model
+    auto modelParams = SF6O2Etching<NumericType, D>::defaultParameters();
+    modelParams.ionFlux = params.get("ionFlux");
+    modelParams.etchantFlux = params.get("etchantFlux");
+    modelParams.passivationFlux = params.get("oxygenFlux");
+    modelParams.Ions.meanEnergy = params.get("meanEnergy");
+    modelParams.Ions.sigmaEnergy = params.get("sigmaEnergy");
+    modelParams.Ions.exponent = params.get("ionExponent");
+    modelParams.Passivation.A_ie = params.get("A_O");
+    modelParams.Substrate.A_ie = params.get("A_Si");
+    modelParams.etchStopDepth = params.get("etchStopDepth");
+    auto model = SmartPointer<SF6O2Etching<NumericType, D>>::New(modelParams);
+
+    CoverageParameters coverageParams;
+    coverageParams.tolerance = 1e-4;
+
+    RayTracingParameters rayTracingParams;
+    rayTracingParams.raysPerPoint = params.get<unsigned>("raysPerPoint");
+
+    AdvectionParameters advectionParams;
+    advectionParams.spatialScheme =
+        util::convertSpatialScheme(params.get<std::string>("spatialScheme"));
+    advectionParams.temporalScheme =
+        util::convertTemporalScheme(params.get<std::string>("temporalScheme"));
+    advectionParams.calculateIntermediateVelocities = intermediateVelocities;
+
+    // process setup
+    const std::string fluxEngineStr = params.get<std::string>("fluxEngine");
+    const auto fluxEngine = util::convertFluxEngineType(fluxEngineStr);
+
+    Process<NumericType, D> process(geometry, model);
+    process.setProcessDuration(params.get("processTime"));
+    process.setParameters(coverageParams);
+    process.setParameters(rayTracingParams);
+    process.setParameters(advectionParams);
+    process.setFluxEngineType(fluxEngine);
+
+    // print initial surface
+    if (suffix == "_noIntermediate")
+      geometry->saveSurfaceMesh("initial.vtp");
+
+    // run the process
+    process.apply();
+
+    // print final surface
+    std::string outputFile = params.get<std::string>("outputFile");
+    auto pos = outputFile.find_last_of('.');
+    if (pos != std::string::npos) {
+      outputFile.insert(pos, suffix);
+    } else {
+      outputFile += suffix;
+    }
+    geometry->saveSurfaceMesh(outputFile, true, 0.01, true);
+  };
+
+  std::cout << "Running simulation without intermediate velocity calculation..."
+            << std::endl;
+  runSimulation(false, "_noIntermediate");
+
+  std::cout << "Running simulation with intermediate velocity calculation..."
+            << std::endl;
+  runSimulation(true, "_intermediate");
 }

examples/holeEtching/holeEtching.py

@@ -13,62 +13,90 @@
 else:
     print("Running 3D simulation.")
 ps.setDimension(args.dim)
+ps.setNumThreads(16)
 
 params = ps.readConfigFile(args.filename)
 
 ps.Length.setUnit(params["lengthUnit"])
 ps.Time.setUnit(params["timeUnit"])
 
-# geometry setup, all units in um
-geometry = ps.Domain(
-    gridDelta=params["gridDelta"],
-    xExtent=params["xExtent"],
-    yExtent=params["yExtent"],
-)
-ps.MakeHole(
-    domain=geometry,
-    holeRadius=params["holeRadius"],
-    holeDepth=0.0,
-    maskHeight=params["maskHeight"],
-    maskTaperAngle=params["taperAngle"],
-    holeShape=ps.HoleShape.HALF,
-).apply()
-
-# use pre-defined model SF6O2 etching model
-modelParams = ps.SF6O2Etching.defaultParameters()
-modelParams.ionFlux = params["ionFlux"]
-modelParams.etchantFlux = params["etchantFlux"]
-modelParams.passivationFlux = params["oxygenFlux"]
-modelParams.Ions.meanEnergy = params["meanEnergy"]
-modelParams.Ions.sigmaEnergy = params["sigmaEnergy"]
-modelParams.Ions.exponent = params["ionExponent"]
-modelParams.Passivation.A_ie = params["A_O"]
-modelParams.Substrate.A_ie = params["A_Si"]
-modelParams.etchStopDepth = params["etchStopDepth"]
-model = ps.SF6O2Etching(modelParams)
-
-covParams = ps.CoverageParameters()
-covParams.tolerance = 1e-4
-
-rayParams = ps.RayTracingParameters()
-rayParams.raysPerPoint = int(params["raysPerPoint"])
-rayParams.smoothingNeighbors = 2
-
-advParams = ps.AdvectionParameters()
-advParams.spatialScheme = ps.util.convertSpatialScheme(params["spatialScheme"])
-
-# process setup
-process = ps.Process(geometry, model)
-process.setProcessDuration(params["processTime"])  # seconds
-process.setParameters(covParams)
-process.setParameters(rayParams)
-process.setParameters(advParams)
-
-# print initial surface
-geometry.saveSurfaceMesh(filename="initial.vtp")
-
-# run the process
-process.apply()
-
-# print final surface
-geometry.saveSurfaceMesh(filename="final.vtp", addInterfaces=True)
+def run_simulation(intermediate_velocities, suffix):
+    # geometry setup, all units in um
+    geometry = ps.Domain(
+        gridDelta=params["gridDelta"],
+        xExtent=params["xExtent"],
+        yExtent=params["yExtent"],
+    )
+
+    hole_shape = ps.HoleShape.QUARTER if args.dim == 3 else ps.HoleShape.HALF
+
+    ps.MakeHole(
+        domain=geometry,
+        holeRadius=params["holeRadius"],
+        holeDepth=0.0,
+        maskHeight=params["maskHeight"],
+        maskTaperAngle=params["taperAngle"],
+        holeShape=hole_shape,
+    ).apply()
+
+    # use pre-defined model SF6O2 etching model
+    modelParams = ps.SF6O2Etching.defaultParameters()
+    modelParams.ionFlux = params["ionFlux"]
+    modelParams.etchantFlux = params["etchantFlux"]
+    modelParams.passivationFlux = params["oxygenFlux"]
+    modelParams.Ions.meanEnergy = params["meanEnergy"]
+    modelParams.Ions.sigmaEnergy = params["sigmaEnergy"]
+    modelParams.Ions.exponent = params["ionExponent"]
+    modelParams.Passivation.A_ie = params["A_O"]
+    modelParams.Substrate.A_ie = params["A_Si"]
+    modelParams.etchStopDepth = params["etchStopDepth"]
+    model = ps.SF6O2Etching(modelParams)
+
+    covParams = ps.CoverageParameters()
+    covParams.tolerance = 1e-4
+
+    rayParams = ps.RayTracingParameters()
+    rayParams.raysPerPoint = int(params["raysPerPoint"])
+
+    advParams = ps.AdvectionParameters()
+    advParams.spatialScheme = ps.util.convertSpatialScheme(
+        params["spatialScheme"]
+    )
+    advParams.temporalScheme = ps.util.convertTemporalScheme(
+        params["temporalScheme"]
+    )
+    advParams.calculateIntermediateVelocities = intermediate_velocities
+
+    # process setup
+    process = ps.Process(geometry, model)
+    process.setProcessDuration(params["processTime"])  # seconds
+    process.setParameters(covParams)
+    process.setParameters(rayParams)
+    process.setParameters(advParams)
+
+    fluxEngineStr = params["fluxEngine"]
+    fluxEngine = ps.util.convertFluxEngineType(fluxEngineStr)
+    process.setFluxEngineType(fluxEngine)
+
+    # print initial surface
+    if suffix == "_noIntermediate":
+        geometry.saveSurfaceMesh(filename="initial.vtp")
+
+    # run the process
+    process.apply()
+
+    # print final surface
+    output_file = params["outputFile"]
+    if "." in output_file:
+        parts = output_file.rsplit(".", 1)
+        output_file = f"{parts[0]}{suffix}.{parts[1]}"
+    else:
+        output_file += suffix
+
+    geometry.saveSurfaceMesh(filename=output_file, addInterfaces=True, boolMaterials=True)
+
+print("Running simulation without intermediate velocity calculation...")
+run_simulation(False, "_noIntermediate")
+
+print("Running simulation with intermediate velocity calculation...")
+run_simulation(True, "_intermediate")

examples/simpleEtching/CMakeLists.txt

@@ -0,0 +1,3 @@
+project(simpleEtching LANGUAGES CXX)
+
+viennaps_add_example(${PROJECT_NAME} "${PROJECT_NAME}.cpp")