New HoleShape parameter in MakeHole and updated Bosch process

examples/boschProcess/boschProcessRayTracing.py

@@ -0,0 +1,146 @@
+#####################################################################
+# 3D Bosch process simulation using ray tracing
+# Replicates the results from:
+#   Ertl and Selberherr https://doi.org/10.1016/j.mee.2009.05.011
+# Execute: 
+#   python boschProcessRayTracing.py -D 3 configBoschRayTracing.txt
+#####################################################################
+
+from argparse import ArgumentParser
+
+# parse config file name and simulation dimension
+parser = ArgumentParser(
+    prog="boschProcess",
+    description="Run a Bosch process on a trench geometry.",
+)
+parser.add_argument("-D", "-DIM", dest="dim", type=int, default=2)
+parser.add_argument("filename")
+args = parser.parse_args()
+
+# switch between 2D and 3D mode
+if args.dim == 2:
+    print("Running 2D simulation.")
+    import viennaps2d as vps
+else:
+    print("Running 3D simulation.")
+    import viennaps3d as vps
+
+params = vps.ReadConfigFile(args.filename)
+
+# print only error output surfaces during the process
+vps.Logger.setLogLevel(vps.LogLevel.ERROR)
+
+# Map the shape string to the corresponding vps.HoleShape enum
+shape_map = {
+    "Full": vps.HoleShape.Full,
+    "Half": vps.HoleShape.Half,
+    "Quarter": vps.HoleShape.Quarter,
+}
+
+hole_shape_str = params.get("holeShape", "Full").strip()
+
+# geometry setup, all units in um
+geometry = vps.Domain()
+vps.MakeHole(
+    domain=geometry,
+    gridDelta=params["gridDelta"],
+    xExtent=params["xExtent"],
+    yExtent=params["yExtent"],
+    holeRadius=params["holeRadius"],
+    holeDepth=params["maskHeight"],
+    taperingAngle=params["taperAngle"],
+    holeShape=shape_map[hole_shape_str],
+    periodicBoundary=False,
+    makeMask=True,
+    material=vps.Material.Si,
+).apply()
+
+depoModel = vps.MultiParticleProcess()
+depoModel.addNeutralParticle(params["stickingDep"])
+depoModel.addIonParticle(params["ionSourceExponent"])
+
+etchModel = vps.MultiParticleProcess()
+materialStickingEtch = {vps.Material.Si: params["stickingEtchPoly"], 
+                         vps.Material.Mask: params["stickingEtchMask"],
+                         vps.Material.Polymer: params["stickingEtchSubs"]}
+etchModel.addNeutralParticle(materialStickingEtch, defaultStickingProbability=1.0)
+etchModel.addIonParticle(params["ionSourceExponent"])
+
+# Conversion factors
+A2c = 1e-8 # 1 angstrom = 1e-8 cm
+c2u = 1e4  # 1 cm = 1e4 micrometers
+
+# Rate calculations for the depo model
+alphaDep   = params["alphaDep"] * (A2c ** 3)      # [cm^3/atom]
+ionDepRate = alphaDep * params["Flux_ionD"] * c2u # [um/s]
+betaDep    = params["betaDep"] * (A2c ** 3)       # [cm^3/atom]
+neuDepRate = betaDep * params["Flux_neuD"] * c2u  # [um/s]
+# Print the result
+print(f"Deposition cycle:")
+print(f"Poly deposit rate (ion, neutral): {ionDepRate:.3e} μm/s, {neuDepRate:.3e} μm/s")
+
+# Rate calculations for the etch model
+alphaPoly   = params["alphaPoly"] * (A2c ** 3)     # [cm^3/atom]
+alphaSubs   = params["alphaSubs"] * (A2c ** 3)     # [cm^3/atom]
+alphaMask   = params["alphaMask"] * (A2c ** 3)     # [cm^3/atom]
+ionEtchRatePoly = alphaPoly * params["Flux_ionE"] * c2u # [um/s]
+ionEtchRateSubs = alphaSubs * params["Flux_ionE"] * c2u # [um/s]
+ionEtchRateMask = alphaMask * params["Flux_ionE"] * c2u # [um/s]
+betaPoly   = params["betaPoly"] * (A2c ** 3)       # [cm^3/atom]
+betaSubs   = params["betaSubs"] * (A2c ** 3)       # [cm^3/atom]
+betaMask   = params["betaMask"] * (A2c ** 3)       # [cm^3/atom]
+neuEtchRatePoly = betaPoly * params["Flux_neuE"] * c2u # [um/s]
+neuEtchRateSubs = betaSubs * params["Flux_neuE"] * c2u # [um/s]
+neuEtchRateMask = betaMask * params["Flux_neuE"] * c2u # [um/s]
+# Print the result
+print(f"Etching cycle:")
+print(f"Mask etching rate (ion, neutral): {ionEtchRateMask:.3e} μm/s, {neuEtchRateMask:.3e} μm/s")
+print(f"Poly etching rate (ion, neutral): {ionEtchRatePoly:.3e} μm/s, {neuEtchRatePoly:.3e} μm/s")
+print(f"Subs etching rate (ion, neutral): {ionEtchRateSubs:.3e} μm/s, {neuEtchRateSubs:.3e} μm/s")
+
+def rateFunctionDep(fluxes, material):
+    rate = fluxes[1] * ionDepRate + fluxes[0] * neuDepRate
+    return rate
+
+# Custom rate function for the etch model
+def rateFunctionEtch(fluxes, material):
+    rate = 0.0
+    if material == vps.Material.Mask:
+        rate = fluxes[1] * ionEtchRateMask + fluxes[0] * neuEtchRateMask
+    if material == vps.Material.Polymer:
+        rate = fluxes[1] * ionEtchRatePoly + fluxes[0] * neuEtchRatePoly
+    if material == vps.Material.Si:
+        rate = fluxes[1] * ionEtchRateSubs + fluxes[0] * neuEtchRateSubs
+    return rate
+
+depoModel.setRateFunction(rateFunctionDep)
+etchModel.setRateFunction(rateFunctionEtch)
+
+depoProcess = vps.Process(geometry, depoModel, params["depTime"])
+depoProcess.setTimeStepRatio(0.2)
+etchProcess = vps.Process(geometry, etchModel, params["etchTime"])
+etchProcess.setTimeStepRatio(0.2)
+
+geometry.saveSurfaceMesh("initial.vtp")
+geometry.saveVolumeMesh("initial")
+
+geometry.duplicateTopLevelSet(vps.Material.Polymer)
+
+numCycles = int(params["numCycles"])
+print(f"---------------------------")
+print(f"Starting {numCycles} cycles")
+for i in range(numCycles):
+    print(f"Cycle {i+1}")
+    # geometry.saveLevelSetMesh(filename=f"run_{2*i}", width=3)
+    geometry.saveSurfaceMesh(f"run_{2*i}.vtp")
+    # geometry.saveVolumeMesh(f"run_{2*i}")
+    print(f"  - Deposition -")
+    depoProcess.apply()
+    # geometry.saveLevelSetMesh(filename=f"run_{2*i+1}", width=3)
+    geometry.saveSurfaceMesh(f"run_{2*i+1}.vtp")
+    # geometry.saveVolumeMesh(f"run_{2*i+1}")
+    print(f"  - Etching -")
+    etchProcess.apply()
+
+geometry.saveSurfaceMesh(f"final.vtp")
+geometry.saveVolumeMesh(f"final")

examples/boschProcess/configBoschRayTracing.txt

@@ -0,0 +1,44 @@
+# Geometry
+holeShape = Quarter
+# all length units are in micrometers (um)
+gridDelta=0.025  # um
+xExtent=3.5      # um
+yExtent=3.5      # um
+
+holeRadius=1.0   # um
+maskHeight=0.6   # um
+taperAngle=0.0   # degree
+
+# Model parameters
+ionSourceExponent = 1000
+numCycles = 10
+
+# Deposition step parameters
+# Sticking
+stickingDep = 0.1
+# Fluxes
+Flux_ionD = 3.125e15    # atoms/cm^2.s
+Flux_neuD = 2.0e18      # atoms/cm^2.s
+# Effecrtive deposition rates
+alphaDep = 10.          # A^3/atom
+betaDep  = 0.5          # A^3/atom
+# Time
+depTime = 5             # s
+
+# Etching step parameters
+# Sticking
+stickingEtchPoly = 0.1
+stickingEtchMask = 0.2
+stickingEtchSubs = 0.2
+# Fluxes
+Flux_ionE = 4.375e15    # atoms/cm^2.s
+Flux_neuE = 1.0e19      # atoms/cm^2.s
+# Effective etch rates
+alphaPoly = -125        # A^3/atom
+alphaSubs = -270        # A^3/atom
+alphaMask = -13.5       # A^3/atom
+betaPoly  = -0.03       # A^3/atom
+betaSubs  = -0.9        # A^3/atom
+betaMask  = -0.045      # A^3/atom
+# Time
+etchTime = 11           # s

include/viennaps/geometries/psMakeHole.hpp

@@ -8,6 +8,8 @@
 
 namespace viennaps {
 
+enum class HoleShape { Full, Half, Quarter };
+
 using namespace viennacore;
 
 /// Generates new a hole geometry in the z direction, which, in 2D mode,
@@ -40,19 +42,33 @@ template <class NumericType, int D> class MakeHole {
   const bool periodicBoundary_;
   const Material material_;
 
+  const HoleShape shape_;
+
 public:
   MakeHole(psDomainType domain, NumericType gridDelta, NumericType xExtent,
            NumericType yExtent, NumericType holeRadius, NumericType holeDepth,
            NumericType taperAngle = 0., NumericType baseHeight = 0.,
            bool periodicBoundary = false, bool makeMask = false,
-           Material material = Material::None)
+           Material material = Material::None,
+           HoleShape shape = HoleShape::Full)
       : domain_(domain), gridDelta_(gridDelta), xExtent_(xExtent),
         yExtent_(yExtent), holeRadius_(holeRadius), holeDepth_(holeDepth),
-        taperAngle_(taperAngle), baseHeight_(baseHeight),
-        periodicBoundary_(periodicBoundary), makeMask_(makeMask),
-        material_(material) {}
+        shape_(shape), taperAngle_(taperAngle), baseHeight_(baseHeight),
+        periodicBoundary_(periodicBoundary && (shape != HoleShape::Half &&
+                                               shape != HoleShape::Quarter)),
+        makeMask_(makeMask), material_(material) {
+    if (periodicBoundary &&
+        (shape == HoleShape::Half || shape == HoleShape::Quarter)) {
+      Logger::getInstance()
+          .addWarning("MakeHole: 'Half' or 'Quarter' shapes do not support "
+                      "periodic boundaries! "
+                      "Defaulting to reflective boundaries!")
+          .print();
+    }
+  }
 
   void apply() {
+
     if constexpr (D != 3) {
       Logger::getInstance()
           .addWarning("MakeHole: Hole geometry can only be created in 3D! "
@@ -65,13 +81,14 @@ template <class NumericType, int D> class MakeHole {
 
       return;
     }
+
     domain_->clear();
     double bounds[2 * D];
-    bounds[0] = -xExtent_ / 2.;
+    bounds[0] = (shape_ != HoleShape::Full) ? 0. : -xExtent_ / 2.;
     bounds[1] = xExtent_ / 2.;
 
     if constexpr (D == 3) {
-      bounds[2] = -yExtent_ / 2.;
+      bounds[2] = (shape_ == HoleShape::Quarter) ? 0. : -yExtent_ / 2.;
       bounds[3] = yExtent_ / 2.;
       bounds[4] = baseHeight_ - gridDelta_;
       bounds[5] = baseHeight_ + holeDepth_ + gridDelta_;

python/pyWrap.cpp

@@ -1007,20 +1007,41 @@ PYBIND11_MODULE(VIENNAPS_MODULE_NAME, module) {
       .def("apply", &MakeTrench<T, D>::apply, "Create a trench geometry.");
 
   // Hole
+  pybind11::enum_<HoleShape>(module, "HoleShape")
+      .value("Full", HoleShape::Full)
+      .value("Half", HoleShape::Half)
+      .value("Quarter", HoleShape::Quarter)
+      .export_values();
+
   pybind11::class_<MakeHole<T, D>>(module, "MakeHole")
       .def(pybind11::init<DomainType, const T, const T, const T, const T,
                           const T, const T, const T, const bool, const bool,
-                          const Material>(),
+                          const Material, HoleShape>(),
            pybind11::arg("domain"), pybind11::arg("gridDelta"),
            pybind11::arg("xExtent"), pybind11::arg("yExtent"),
            pybind11::arg("holeRadius"), pybind11::arg("holeDepth"),
            pybind11::arg("taperingAngle") = 0.,
            pybind11::arg("baseHeight") = 0.,
            pybind11::arg("periodicBoundary") = false,
            pybind11::arg("makeMask") = false,
-           pybind11::arg("material") = Material::None)
+           pybind11::arg("material") = Material::None,
+           pybind11::arg("holeShape") = HoleShape::Full) // New argument
       .def("apply", &MakeHole<T, D>::apply, "Create a hole geometry.");
 
+  //  pybind11::class_<MakeHole<T, D>>(module, "MakeHole")
+  //      .def(pybind11::init<DomainType, const T, const T, const T, const T,
+  //                          const T, const T, const T, const bool, const bool,
+  //                          const Material>(),
+  //           pybind11::arg("domain"), pybind11::arg("gridDelta"),
+  //           pybind11::arg("xExtent"), pybind11::arg("yExtent"),
+  //           pybind11::arg("holeRadius"), pybind11::arg("holeDepth"),
+  //           pybind11::arg("taperingAngle") = 0.,
+  //           pybind11::arg("baseHeight") = 0.,
+  //           pybind11::arg("periodicBoundary") = false,
+  //           pybind11::arg("makeMask") = false,
+  //           pybind11::arg("material") = Material::None)
+  //      .def("apply", &MakeHole<T, D>::apply, "Create a hole geometry.");
+
   // Fin
   pybind11::class_<MakeFin<T, D>>(module, "MakeFin")
       .def(pybind11::init<DomainType, const T, const T, const T, const T,

tests/hole/hole.cpp

@@ -10,8 +10,33 @@ using namespace viennaps;
 template <class NumericType, int D> void RunTest() {
   auto domain = SmartPointer<Domain<NumericType, D>>::New();
 
+  // Test with HoleShape::Full
   MakeHole<NumericType, D>(domain, 1., 10., 10., 2.5, 5., 10., 1., false, true,
-                           Material::Si)
+                           Material::Si, HoleShape::Full)
+      .apply();
+
+  VC_TEST_ASSERT(domain->getLevelSets().size() == 2);
+  VC_TEST_ASSERT(domain->getMaterialMap());
+  VC_TEST_ASSERT(domain->getMaterialMap()->size() == 2);
+
+  LSTEST_ASSERT_VALID_LS(domain->getLevelSets().back(), NumericType, D);
+
+  // Test with HoleShape::Quarter
+  domain->clear(); // Reset the domain for a new test
+  MakeHole<NumericType, D>(domain, 1., 10., 10., 2.5, 5., 10., 1., false, true,
+                           Material::Si, HoleShape::Quarter)
+      .apply();
+
+  VC_TEST_ASSERT(domain->getLevelSets().size() == 2);
+  VC_TEST_ASSERT(domain->getMaterialMap());
+  VC_TEST_ASSERT(domain->getMaterialMap()->size() == 2);
+
+  LSTEST_ASSERT_VALID_LS(domain->getLevelSets().back(), NumericType, D);
+
+  // Test with HoleShape::Half
+  domain->clear(); // Reset the domain for a new test
+  MakeHole<NumericType, D>(domain, 1., 10., 10., 2.5, 5., 10., 1., false, true,
+                           Material::Si, HoleShape::Half)
       .apply();
 
   VC_TEST_ASSERT(domain->getLevelSets().size() == 2);