Extended material map to cell set and translation field (#45)

* Material mapping in cell set

* Material mapping in velocity field and implementation in stack redeposition

* Check material utility function

* Remove top LS function

Examples/OxideRegrowth/OxideRegrowth.cpp

@@ -62,11 +62,10 @@ int main(int argc, char **argv) {
   // process in the cell set. The byproducts are then distributed by solving a
   // convection-diffusion equation on the cell set.
   auto model = psSmartPointer<OxideRegrowthModel<NumericType, D>>::New(
-      params.numLayers + 1, params.nitrideEtchRate / 60,
-      params.oxideEtchRate / 60, params.redepositionRate,
-      params.redepositionThreshold, params.redepositionTimeInt,
-      params.diffusionCoefficient, params.sink, params.scallopVelocity,
-      params.centerVelocity,
+      params.nitrideEtchRate / 60, params.oxideEtchRate / 60,
+      params.redepositionRate, params.redepositionThreshold,
+      params.redepositionTimeInt, params.diffusionCoefficient, params.sink,
+      params.scallopVelocity, params.centerVelocity,
       params.substrateHeight + params.numLayers * params.layerHeight,
       params.trenchWidth);
 

Examples/TrenchDeposition/Parameters.hpp

@@ -17,7 +17,7 @@ template <typename T> struct Parameters {
   T taperAngle = 0.;
 
   // Process
-  T processTime = 5.;
+  T processTime = 0.5;
   T stickingProbability = 0.1;
   T sourcePower = 1.;
 

include/CellSet/csDenseCellSet.hpp

@@ -14,6 +14,7 @@
 #include <rayUtil.hpp>
 
 #include <psLogger.hpp>
+#include <psMaterials.hpp>
 #include <psSmartPointer.hpp>
 #include <psVTKWriter.hpp>
 
@@ -26,11 +27,13 @@ template <class T, int D> class csDenseCellSet {
   using gridType = psSmartPointer<lsMesh<T>>;
   using levelSetsType =
       psSmartPointer<std::vector<psSmartPointer<lsDomain<T, D>>>>;
+  using materialMapType = psSmartPointer<psMaterialMap>;
 
   levelSetsType levelSets = nullptr;
   gridType cellGrid = nullptr;
   psSmartPointer<lsDomain<T, D>> surface = nullptr;
   psSmartPointer<csBVH<T, D>> BVH = nullptr;
+  materialMapType materialMap = nullptr;
   std::vector<std::array<int, 2 * D>> cellNeighbors; // -x, x, -y, y, -z, z
   T gridDelta;
   size_t numberOfCells;
@@ -50,8 +53,11 @@ template <class T, int D> class csDenseCellSet {
     fromLevelSets(passedLevelSets, passedDepth);
   }
 
-  void fromLevelSets(levelSetsType passedLevelSets, T passedDepth = 0.) {
+  void fromLevelSets(levelSetsType passedLevelSets,
+                     materialMapType passedMaterialMap = nullptr,
+                     T passedDepth = 0.) {
     levelSets = passedLevelSets;
+    materialMap = passedMaterialMap;
 
     if (cellGrid == nullptr)
       cellGrid = psSmartPointer<lsMesh<T>>::New();
@@ -79,8 +85,9 @@ template <class T, int D> class csDenseCellSet {
       levelSetsInOrder.push_back(plane);
     for (auto ls : *levelSets)
       levelSetsInOrder.push_back(ls);
-    if (cellSetAboveSurface)
+    if (cellSetAboveSurface) {
       levelSetsInOrder.push_back(plane);
+    }
 
     calculateMinMaxIndex(levelSetsInOrder);
     lsToVoxelMesh<T, D>(levelSetsInOrder, cellGrid).apply();
@@ -97,7 +104,7 @@ template <class T, int D> class csDenseCellSet {
     psVTKWriter<T>(cellGrid, "cellSet_debug_init.vtu").apply();
 #endif
 
-    if (!cellSetAboveSurface)
+    if (!cellSetAboveSurface || materialMap)
       adjustMaterialIds();
 
     // create filling fractions as default scalar cell data
@@ -177,7 +184,7 @@ template <class T, int D> class csDenseCellSet {
     return getFillingFractions()->at(idx);
   }
 
-  int getIndex(std::array<T, 3> &point) { return findIndex(point); }
+  int getIndex(const std::array<T, 3> &point) { return findIndex(point); }
 
   std::vector<T> *getScalarData(std::string name) {
     return cellGrid->getCellData().getScalarData(name);
@@ -358,6 +365,7 @@ template <class T, int D> class csDenseCellSet {
 
     // move iterator for lowest material id and then adjust others if they are
     // needed
+    const materialMapType matMapPtr = materialMap;
     unsigned cellIdx = 0;
     for (; iterators.front().getIndices() < maxIndex;
          iterators.front().next()) {
@@ -391,7 +399,12 @@ template <class T, int D> class csDenseCellSet {
           }
 
           if (isVoxel) {
-            materialIds->at(cellIdx++) = materialId;
+            if (matMapPtr) {
+              auto material = matMapPtr->getMaterialAtIdx(materialId);
+              materialIds->at(cellIdx++) = static_cast<int>(material);
+            } else {
+              materialIds->at(cellIdx++) = materialId;
+            }
           }
 
           // jump out of material for loop
@@ -544,8 +557,16 @@ template <class T, int D> class csDenseCellSet {
 
 #pragma omp parallel for
     for (size_t i = 0; i < matIds->size(); i++) {
-      if (matIds->at(i) > 0) {
-        matIds->at(i) -= 1;
+      int materialId = static_cast<int>(matIds->at(i));
+
+      if (!cellSetAboveSurface && materialId > 0) {
+        materialId -= 1;
+      }
+
+      assert(materialId >= 0);
+      if (materialMap) {
+        matIds->at(i) =
+            static_cast<int>(materialMap->getMaterialAtIdx(materialId));
       }
     }
   }

include/Models/FluorocarbonEtching.hpp

@@ -63,8 +63,7 @@ class FluorocarbonSurfaceModel : public psSurfaceModel<NumericType> {
         break;
       }
 
-      auto matId =
-          psMaterialMap::mapToMaterial(static_cast<int>(materialIds[i]));
+      auto matId = psMaterialMap::mapToMaterial(materialIds[i]);
       assert(matId == psMaterial::Mask || matId == psMaterial::Polymer ||
              matId == psMaterial::Si || matId == psMaterial::SiO2 ||
              matId == psMaterial::Si3N4 && "Unexptected material");

include/Models/SF6O2Etching.hpp

@@ -65,8 +65,7 @@ class SF6O2SurfaceModel : public psSurfaceModel<NumericType> {
         break;
       }
 
-      if (psMaterialMap::mapToMaterial(static_cast<int>(materialIds[i])) ==
-          psMaterial::Si) {
+      if (psMaterialMap::isMaterial(materialIds[i], psMaterial::Si)) {
 
         etchRate[i] =
             -(1 / rho_Si) *

include/Models/SimpleDeposition.hpp

@@ -22,23 +22,6 @@ class SimpleDepositionSurfaceModel : public psSurfaceModel<NumericType> {
   }
 };
 
-template <class T>
-class SimpleDepositionVelocityField : public psVelocityField<T> {
-public:
-  T getScalarVelocity(const std::array<T, 3> & /*coordinate*/, int /*material*/,
-                      const std::array<T, 3> & /*normalVector*/,
-                      unsigned long pointID) override {
-    return velocities->at(pointID);
-  }
-
-  void setVelocities(psSmartPointer<std::vector<T>> passedVelocities) override {
-    velocities = passedVelocities;
-  }
-
-private:
-  psSmartPointer<std::vector<T>> velocities = nullptr;
-};
-
 template <typename NumericType, int D>
 class SimpleDepositionParticle
     : public rayParticle<SimpleDepositionParticle<NumericType, D>,
@@ -55,7 +38,7 @@ class SimpleDepositionParticle
                         rayTracingData<NumericType> &localData,
                         const rayTracingData<NumericType> *globalData,
                         rayRNG &Rng) override final {
-    localData.getVectorData(0)[primID] += rayWeight * stickingProbability;
+    localData.getVectorData(0)[primID] += rayWeight;
   }
   std::pair<NumericType, rayTriple<NumericType>>
   surfaceReflection(NumericType rayWeight, const rayTriple<NumericType> &rayDir,
@@ -96,8 +79,7 @@ class SimpleDeposition : public psProcessModel<NumericType, D> {
         psSmartPointer<SimpleDepositionSurfaceModel<NumericType, D>>::New();
 
     // velocity field
-    auto velField =
-        psSmartPointer<SimpleDepositionVelocityField<NumericType>>::New();
+    auto velField = psSmartPointer<psDefaultVelocityField<NumericType>>::New();
 
     this->setSurfaceModel(surfModel);
     this->setVelocityField(velField);