Add solid core to tubenetwork mesh

src/scaffoldmaker/meshtypes/meshtype_3d_tubenetwork1.py

@@ -1,15 +1,83 @@
 """
-Generates a 3-D Hermite bifurcating tube network.
+Generates a 3-D Hermite bifurcating tube network (with optional solid core).
 """
 from scaffoldmaker.meshtypes.meshtype_1d_network_layout1 import MeshType_1d_network_layout1
 from scaffoldmaker.meshtypes.scaffold_base import Scaffold_base
 from scaffoldmaker.scaffoldpackage import ScaffoldPackage
 from scaffoldmaker.utils.tubenetworkmesh import TubeNetworkMeshBuilder, TubeNetworkMeshGenerateData
 
 
+def calculateElementsCountAcrossMinor(cls, options):
+    # Calculate elementsCountAcrosMinor
+    elementsCountAcrossMinor = int(((options["Elements count around"] - 4) / 4 -
+                                    (options['Number of elements across major'] / 2)) * 2 + 6)
+
+    return elementsCountAcrossMinor
+
+def setParametersToDefault(cls, options):
+    options["Elements count around"] = 8
+    options["Number of elements across major"] = 4
+    options["Number of elements across transition"] = 1
+
+    return options
+
+def checkCoreParameters(cls, options, dependentChanges=False):
+    # Check elements count around are ok
+    if options["Elements count around"] < 8:
+        dependentChanges = True
+        options = setParametersToDefault(cls, options)
+    if options["Elements count around"] % 4:
+        options["Elements count around"] += 4 - options["Elements count around"] % 4
+
+    # Check elements count across major are ok
+    if options['Number of elements across major'] < 4:
+        options['Number of elements across major'] = 4
+    if options['Number of elements across major'] % 2:
+        options['Number of elements across major'] += 1
+
+    # Check elements count across transition
+    if options['Number of elements across transition'] < 1:
+        options['Number of elements across transition'] = 1
+
+    # Calculate elementsCountAcrossMinor based on the current set of elementsCountAround and elementsCountAcrossMajor
+    elementsCountAcrossMinor = calculateElementsCountAcrossMinor(cls, options)
+
+    # Rcrit check
+    Rcrit = max(
+        min(options['Number of elements across major'] - 4, elementsCountAcrossMinor - 4) // 2 + 1, 0)
+    if Rcrit < options['Number of elements across transition']:
+        dependentChanges = True
+        options['Number of elements across transition'] = Rcrit
+
+    # Number of elements around sanity check
+    eM = options["Number of elements across major"]
+    em = elementsCountAcrossMinor
+    eC = (eM - 1) * (em - 1) - ((eM - 3) * (em - 3))
+    if options["Elements count around"] != eC:
+        dependentChanges = True
+        # Reset parameter values
+        setParametersToDefault(cls, options)
+
+    annotationElementsCountsAround = options["Annotation elements counts around"]
+    if len(annotationElementsCountsAround) == 0:
+        options["Annotation elements count around"] = [0]
+    else:
+        for i in range(len(annotationElementsCountsAround)):
+            if annotationElementsCountsAround[i] <= 0:
+                annotationElementsCountsAround[i] = 0
+            elif annotationElementsCountsAround[i] < 4:
+                annotationElementsCountsAround[i] = 4
+
+    return options, dependentChanges
+
+
 class MeshType_3d_tubenetwork1(Scaffold_base):
     """
     Generates a 3-D hermite tube network, with linear or cubic basis through wall.
+    Number of elements generated are primarily controlled by "elements count around" and "elements count across major".
+    The elements count around parameter determines the number of elements around a circular rim.
+    The elements count across major determines the number of elements across the major axis of an ellipse (y-axis in the scaffold).
+    The number of elements across minor axis (z-axis) is calculated internally based on the two elements count parameters.
     """
 
     @classmethod
@@ -31,7 +99,11 @@ def getDefaultOptions(cls, parameterSetName="Default"):
             "Annotation elements counts around": [0],
             "Target element density along longest segment": 4.0,
             "Use linear through wall": False,
-            "Show trim surfaces": False
+            "Show trim surfaces": False,
+            "Core": False,
+            'Number of elements across major': 4,
+            'Number of elements across transition': 1,
+            "Annotation elements counts across major": [0]
         }
         return options
 
@@ -44,7 +116,11 @@ def getOrderedOptionNames():
             "Annotation elements counts around",
             "Target element density along longest segment",
             "Use linear through wall",
-            "Show trim surfaces"
+            "Show trim surfaces",
+            "Core",
+            'Number of elements across major',
+            'Number of elements across transition',
+            "Annotation elements counts across major"
         ]
 
     @classmethod
@@ -82,23 +158,33 @@ def getOptionScaffoldPackage(cls, optionName, scaffoldType, parameterSetName=Non
     def checkOptions(cls, options):
         if not options["Network layout"].getScaffoldType() in cls.getOptionValidScaffoldTypes("Network layout"):
             options["Network layout"] = cls.getOptionScaffoldPackage("Network layout", MeshType_1d_network_layout1)
-        elementsCountAround = options["Elements count around"]
-        if options["Elements count around"] < 4:
-            options["Elements count around"] = 4
+        dependentChanges = False
+
+        # Parameters specific to the core
+        if options["Core"] == True:
+            options, dependentChanges = checkCoreParameters(cls, options, dependentChanges)
+
+        # When core option is false
+        else:
+            if options["Elements count around"] < 4:
+                options["Elements count around"] = 4
+            annotationElementsCountsAround = options["Annotation elements counts around"]
+            if len(annotationElementsCountsAround) == 0:
+                options["Annotation elements count around"] = [0]
+            else:
+                for i in range(len(annotationElementsCountsAround)):
+                    if annotationElementsCountsAround[i] <= 0:
+                        annotationElementsCountsAround[i] = 0
+                    elif annotationElementsCountsAround[i] < 4:
+                        annotationElementsCountsAround[i] = 4
+
+        # Common parameters
         if options["Elements count through wall"] < 1:
             options["Elements count through wall"] = 1
-        annotationElementsCountsAround = options["Annotation elements counts around"]
-        if len(annotationElementsCountsAround) == 0:
-            options["Annotation elements count around"] = [0]
-        else:
-            for i in range(len(annotationElementsCountsAround)):
-                if annotationElementsCountsAround[i] <= 0:
-                    annotationElementsCountsAround[i] = 0
-                elif annotationElementsCountsAround[i] < 4:
-                    annotationElementsCountsAround[i] = 4
+
         if options["Target element density along longest segment"] < 1.0:
             options["Target element density along longest segment"] = 1.0
-        dependentChanges = False
+
         return dependentChanges
 
     @classmethod
@@ -112,15 +198,21 @@ def generateBaseMesh(cls, region, options):
         layoutRegion = region.createRegion()
         networkLayout = options["Network layout"]
         networkLayout.generate(layoutRegion)  # ask scaffold to generate to get user-edited parameters
+        layoutAnnotationGroups = networkLayout.getAnnotationGroups()
         networkMesh = networkLayout.getConstructionObject()
 
         tubeNetworkMeshBuilder = TubeNetworkMeshBuilder(
             networkMesh,
             targetElementDensityAlongLongestSegment=options["Target element density along longest segment"],
             defaultElementsCountAround=options["Elements count around"],
             elementsCountThroughWall=options["Elements count through wall"],
-            layoutAnnotationGroups=networkLayout.getAnnotationGroups(),
-            annotationElementsCountsAround=options["Annotation elements counts around"])
+            layoutAnnotationGroups=layoutAnnotationGroups,
+            annotationElementsCountsAround=options["Annotation elements counts around"],
+            defaultElementsCountAcrossMajor=options['Number of elements across major'],
+            elementsCountTransition=options['Number of elements across transition'],
+            annotationElementsCountsAcrossMajor=options["Annotation elements counts across major"],
+            isCore=options["Core"])
+
         tubeNetworkMeshBuilder.build()
         generateData = TubeNetworkMeshGenerateData(
             region, 3,