Skip to content

Commit

Permalink
Add discontinuity on thorax, abdomen core boundaries
Browse files Browse the repository at this point in the history
  • Loading branch information
rchristie committed Oct 27, 2024
1 parent 33bf002 commit 47b22e5
Show file tree
Hide file tree
Showing 3 changed files with 95 additions and 20 deletions.
8 changes: 8 additions & 0 deletions src/scaffoldmaker/meshtypes/meshtype_3d_wholebody2.py
Original file line number Diff line number Diff line change
Expand Up @@ -775,9 +775,12 @@ def generateBaseMesh(cls, region, options):

annotationAlongCounts = []
annotationAroundCounts = []
defaultCoreBoundaryScalingMode = 1
annotationCoreBoundaryScalingMode = []
for layoutAnnotationGroup in layoutAnnotationGroups:
alongCount = 0
aroundCount = 0
coreBoundaryScalingMode = 0
name = layoutAnnotationGroup.getName()
if "head" in name:
alongCount = elementsCountAlongHead
Expand All @@ -788,9 +791,11 @@ def generateBaseMesh(cls, region, options):
elif "thorax" in name:
alongCount = elementsCountAlongThorax
aroundCount = elementsCountAroundTorso
coreBoundaryScalingMode = 2
elif "abdomen" in name:
alongCount = elementsCountAlongAbdomen
aroundCount = elementsCountAroundTorso
coreBoundaryScalingMode = 2
elif "arm to hand" in name:
alongCount = elementsCountAlongArmToHand
aroundCount = elementsCountAroundArm
Expand All @@ -805,6 +810,7 @@ def generateBaseMesh(cls, region, options):
aroundCount = elementsCountAroundLeg
annotationAlongCounts.append(alongCount)
annotationAroundCounts.append(aroundCount)
annotationCoreBoundaryScalingMode.append(coreBoundaryScalingMode)

tubeNetworkMeshBuilder = BodyTubeNetworkMeshBuilder(
networkMesh,
Expand All @@ -818,6 +824,8 @@ def generateBaseMesh(cls, region, options):
elementsCountTransition=options['Number of elements across core transition'],
defaultElementsCountCoreBoxMinor=options["Number of elements across core box minor"],
annotationElementsCountsCoreBoxMinor=[],
defaultCoreBoundaryScalingMode=defaultCoreBoundaryScalingMode,
annotationCoreBoundaryScalingMode=annotationCoreBoundaryScalingMode,
useOuterTrimSurfaces=True)

meshDimension = 3
Expand Down
26 changes: 17 additions & 9 deletions src/scaffoldmaker/utils/eft_utils.py
Original file line number Diff line number Diff line change
Expand Up @@ -973,8 +973,8 @@ def getTricubicHermiteSerendipityElementNodeParameter(eft, scalefactors, nodePar
Currently assumes extra versions are not used, but this is not checked.
:param eft: Element field template.
:param scalefactors: List of scale factors.
:param nodeParameters: List over 8 (3-D)
local nodes in Zinc ordering of 4 parameter vectors x, d1, d2, d3 each with 3 components.
:param nodeParameters: List over 8 (3-D) local nodes in Zinc ordering of 4 parameter vectors
x, d1, d2, d3 each with 3 components.
:param localNodeIndex: 1-based local node index on original basis.
:param valueLabel: A value label mapped by CubicHermiteSerendipityFunctionOffset.
:return: Parameter.
Expand All @@ -1001,22 +1001,26 @@ def getTricubicHermiteSerendipityElementNodeParameter(eft, scalefactors, nodePar
return parameter


def addTricubicHermiteSerendipityEftParameterScaling(eft, scalefactors, nodeParameters, localNodeIndexes, valueLabel):
def addTricubicHermiteSerendipityEftParameterScaling(eft, scalefactors, nodeParameters, localNodeIndexes, valueLabel,
version=1):
"""
Scale tricubic hermite serendipity element field template value label at local nodes to fit actual parameters.
Parameter must currently be an unscaled single term expression at that local node.
Calculation assumes value label parameter is an unscaled single term expression at the local nodes.
:param eft: Element field template e.g. returned from determineCubicHermiteSerendipityEft.
Note: be careful to pass a separate eft to this function from any cached which do not needing scaling!
:param scalefactors: Existing scale factors for element
:param nodeParameters: As passed to determineCubicHermiteSerendipityEft: list over 8 (3-D)
local nodes in Zinc ordering of 4 parameter vectors x, d1, d2, d3 each with 3 components.
:param localNodeIndexes: Local node indexes to scale value label at. Currently must be in [5, 6, 7, 8].
:param valueLabel: Single value label to scale. Currently only implemened for D_DS3.
:return: Modified eft, scalefactors
:param version: Set to a number > 1 to use that derivative version instead of scale factors, and client is
responsible for assigning that derivative version in proportion to the returned scale factors.
:return: Modified eft, final scalefactors, add scalefactors (multiplying the affected derivatives; these are
included in final scalefactors if version == 1, otherwise client must use these to scale versioned derivatives).
"""
assert len(nodeParameters) == 8
assert valueLabel == Node.VALUE_LABEL_D_DS3
newScalefactors = copy.copy(scalefactors) if scalefactors else []
addScalefactors = []
for ln in localNodeIndexes:
assert ln in [5, 6, 7, 8]
na = ln - 5
Expand All @@ -1027,6 +1031,10 @@ def addTricubicHermiteSerendipityEftParameterScaling(eft, scalefactors, nodePara
db = nodeParameters[nb][3]
dbScaled = computeCubicHermiteEndDerivative(xa, da, xb, db)
scalefactor = magnitude(dbScaled) / magnitude(db)
newScalefactors.append(scalefactor)
addScaleEftNodesValueLabel(eft, localNodeIndexes, Node.VALUE_LABEL_D_DS3, 3)
return eft, newScalefactors
addScalefactors.append(scalefactor)
if version == 1:
newScalefactors = scalefactors + addScalefactors if scalefactors else addScalefactors
addScaleEftNodesValueLabel(eft, localNodeIndexes, Node.VALUE_LABEL_D_DS3, 3)
return eft, newScalefactors, addScalefactors
remapEftNodeValueLabelsVersion(eft, localNodeIndexes, [valueLabel], version)
return eft, scalefactors, addScalefactors
81 changes: 70 additions & 11 deletions src/scaffoldmaker/utils/tubenetworkmesh.py
Original file line number Diff line number Diff line change
Expand Up @@ -189,11 +189,50 @@ def getNewTrimAnnotationGroup(self):
self._trimAnnotationGroupCount += 1
return self.getOrCreateAnnotationGroup(("trim surface " + "{:03d}".format(self._trimAnnotationGroupCount), ""))

def resolveEftCoreBoundaryScaling(self, eft, scalefactors, nodeParameters, nodeIdentifiers, mode):
"""
Resolve differences in d3 scaling across core-shell boundary by one of 2 modes.
Works for 8-noded tricubic Hermite serendipity basis only.
:param eft: Element field template to modify.
:param scalefactors: Existing scalefactors for use with eft.
:param nodeParameters: List over 8 (3-D) local nodes in Zinc ordering of 4 parameter vectors
x, d1, d2, d3 each with 3 components.
:param nodeIdentifiers: List over 8 3-D local nodes giving global node identifiers.
:param mode: 1 to set scale factors, 2 to add version 2 to d3 for the boundary nodes and assigning
values to that version equal to the scale factors x version 1.
:return: New eft, new scalefactors.
"""
assert mode in (1, 2)
eft, scalefactors, addScalefactors = addTricubicHermiteSerendipityEftParameterScaling(
eft, scalefactors, nodeParameters, [5, 6, 7, 8], Node.VALUE_LABEL_D_DS3, version=mode)
if mode == 2:
nodetemplate = self._nodes.createNodetemplate()
n = 4
for nodeIdentifier, scalefactor in zip(nodeIdentifiers[4:], addScalefactors):
node = self._nodes.findNodeByIdentifier(nodeIdentifier)
nodetemplate.defineFieldFromNode(self._coordinates, node)
versionsCount = nodetemplate.getValueNumberOfVersions(self._coordinates, -1, Node.VALUE_LABEL_D_DS3)
if versionsCount == 1:
# make version 2 of d3 at the node and assign to it
nodetemplate.setValueNumberOfVersions(self._coordinates, -1, Node.VALUE_LABEL_D_DS3, 2)
# merge clears the current parameters so need to re-set
node.merge(nodetemplate)
self._fieldcache.setNode(node)
for valueLabel, value in zip(
(Node.VALUE_LABEL_VALUE, Node.VALUE_LABEL_D_DS1,
Node.VALUE_LABEL_D_DS2, Node.VALUE_LABEL_D_DS3), nodeParameters[n]):
self._coordinates.setNodeParameters(self._fieldcache, -1, valueLabel, 1, value)
d3_v2 = mult(nodeParameters[n][3], scalefactor)
self._coordinates.setNodeParameters(self._fieldcache, -1, Node.VALUE_LABEL_D_DS3, 2, d3_v2)
n += 1
return eft, scalefactors


class TubeNetworkMeshSegment(NetworkMeshSegment):

def __init__(self, networkSegment, pathParametersList, elementsCountAround, elementsCountThroughShell,
isCore=False, elementsCountCoreBoxMinor: int = 2, elementsCountTransition: int = 1):
isCore=False, elementsCountCoreBoxMinor: int=2, elementsCountTransition: int=1,
coreBoundaryScalingMode: int=1):
"""
:param networkSegment: NetworkSegment this is built from.
:param pathParametersList: [pathParameters] if 2-D or [outerPathParameters, innerPathParameters] if 3-D
Expand All @@ -210,6 +249,7 @@ def __init__(self, networkSegment, pathParametersList, elementsCountAround, elem
self._elementsCountCoreBoxMajor = (elementsCountAround // 2) - elementsCountCoreBoxMinor
self._elementsCountCoreBoxMinor = elementsCountCoreBoxMinor
self._elementsCountTransition = elementsCountTransition
self._coreBoundaryScalingMode = coreBoundaryScalingMode

assert elementsCountThroughShell > 0
self._elementsCountThroughShell = elementsCountThroughShell
Expand Down Expand Up @@ -241,6 +281,9 @@ def __init__(self, networkSegment, pathParametersList, elementsCountAround, elem
# [nAlong][nAcrossMajor][nAcrossMinor] format.
self._boxElementIds = None # [along][major][minor]

def getCoreBoundaryScalingMode(self):
return self._coreBoundaryScalingMode

def getElementsCountAround(self):
return self._elementsCountAround

Expand Down Expand Up @@ -1586,8 +1629,8 @@ def generateMesh(self, generateData: TubeNetworkMeshGenerateData, n2Only=None):
nodeLayouts.append(None)
eft, scalefactors = determineCubicHermiteSerendipityEft(mesh, nodeParameters, nodeLayouts)
if self._elementsCountTransition == 1:
eft, scalefactors = addTricubicHermiteSerendipityEftParameterScaling(
eft, scalefactors, nodeParameters, [5, 6, 7, 8], Node.VALUE_LABEL_D_DS3)
eft, scalefactors = generateData.resolveEftCoreBoundaryScaling(
eft, scalefactors, nodeParameters, nids, self._coreBoundaryScalingMode)
elementtemplate = mesh.createElementtemplate()
elementtemplate.setElementShapeType(Element.SHAPE_TYPE_CUBE)
elementtemplate.defineField(coordinates, -1, eft)
Expand Down Expand Up @@ -1625,8 +1668,8 @@ def generateMesh(self, generateData: TubeNetworkMeshGenerateData, n2Only=None):
for n1 in (e1, n1p):
nodeParameters.append(self.getRimCoordinates(n1, n2, n3))
eft = generateData.createElementfieldtemplate()
eft, scalefactors = addTricubicHermiteSerendipityEftParameterScaling(
eft, scalefactors, nodeParameters, [5, 6, 7, 8], Node.VALUE_LABEL_D_DS3)
eft, scalefactors = generateData.resolveEftCoreBoundaryScaling(
eft, scalefactors, nodeParameters, nids, self._coreBoundaryScalingMode)
elementtemplateTransition = mesh.createElementtemplate()
elementtemplateTransition.setElementShapeType(Element.SHAPE_TYPE_CUBE)
elementtemplateTransition.defineField(coordinates, -1, eft)
Expand Down Expand Up @@ -2851,8 +2894,8 @@ def _generateTransitionElements(self, s, n2, mesh, elementtemplate, coordinates,

eft, scalefactors = determineCubicHermiteSerendipityEft(mesh, nodeParameters, nodeLayouts)
if elementsCountTransition == 1:
eft, scalefactors = addTricubicHermiteSerendipityEftParameterScaling(
eft, scalefactors, nodeParameters, [5, 6, 7, 8], Node.VALUE_LABEL_D_DS3)
eft, scalefactors = generateData.resolveEftCoreBoundaryScaling(
eft, scalefactors, nodeParameters, nids, segment.getCoreBoundaryScalingMode())
elementtemplate.defineField(coordinates, -1, eft)
element = mesh.createElement(elementIdentifier, elementtemplate)
element.setNodesByIdentifier(eft, nids)
Expand Down Expand Up @@ -3037,8 +3080,8 @@ def generateMesh(self, generateData: TubeNetworkMeshGenerateData):
if lastTransition:
# need to generate eft again otherwise modifying object in eftList, mucking up outer layers
eft, scalefactors = determineCubicHermiteSerendipityEft(mesh, nodeParameters, nodeLayouts)
eft, scalefactors = addTricubicHermiteSerendipityEftParameterScaling(
eft, scalefactors, nodeParameters, [5, 6, 7, 8], Node.VALUE_LABEL_D_DS3)
eft, scalefactors = generateData.resolveEftCoreBoundaryScaling(
eft, scalefactors, nodeParameters, nids, segment.getCoreBoundaryScalingMode())

elementtemplate.defineField(coordinates, -1, eft)
element = mesh.createElement(elementIdentifier, elementtemplate)
Expand All @@ -3060,6 +3103,7 @@ def __init__(self, networkMesh: NetworkMesh, targetElementDensityAlongLongestSeg
defaultElementsCountAround: int=8, annotationElementsCountsAround: list=[],
elementsCountThroughShell: int=1, isCore=False, elementsCountTransition: int=1,
defaultElementsCountCoreBoxMinor: int=2, annotationElementsCountsCoreBoxMinor: list=[],
defaultCoreBoundaryScalingMode=1, annotationCoreBoundaryScalingMode=[],
useOuterTrimSurfaces=True):
"""
Builds contiguous tube network meshes with smooth element size transitions at junctions, optionally with solid
Expand All @@ -3085,6 +3129,9 @@ def __init__(self, networkMesh: NetworkMesh, targetElementDensityAlongLongestSeg
elements across core box minor/d3 direction for segments with any elements in the annotation group. Client must
ensure exclusive map from segments. Groups with zero value or past end of this list use the
defaultElementsCountCoreBoxMinor.
:param defaultCoreBoundaryScalingMode: Mode 1=scalefactor, 2=version.
:param annotationCoreBoundaryScalingMode: Boundary mode 1 or 2 in order of layoutAnnotationGroups,
or 0 to use default.
:param useOuterTrimSurfaces: Set to False to use separate trim surfaces on inner and outer tubes. Ignored if
no inner path.
"""
Expand All @@ -3097,6 +3144,8 @@ def __init__(self, networkMesh: NetworkMesh, targetElementDensityAlongLongestSeg
self._elementsCountTransition = elementsCountTransition
self._defaultElementsCountCoreBoxMinor = defaultElementsCountCoreBoxMinor
self._annotationElementsCountsCoreBoxMinor = annotationElementsCountsCoreBoxMinor
self._defaultCoreBoundaryScalingMode = defaultCoreBoundaryScalingMode
self._annotationCoreBoundaryScalingMode = annotationCoreBoundaryScalingMode
layoutFieldmodule = self._layoutRegion.getFieldmodule()
self._layoutInnerCoordinates = layoutFieldmodule.findFieldByName("inner coordinates").castFiniteElement()
if not self._layoutInnerCoordinates.isValid():
Expand All @@ -3114,6 +3163,7 @@ def createSegment(self, networkSegment):
elementsCountAround = self._defaultElementsCountAround
elementsCountCoreBoxMinor = self._defaultElementsCountCoreBoxMinor

coreBoundaryScalingMode = self._defaultCoreBoundaryScalingMode
i = 0
for layoutAnnotationGroup in self._layoutAnnotationGroups:
if i >= len(self._annotationElementsCountsAround):
Expand All @@ -3135,10 +3185,19 @@ def createSegment(self, networkSegment):
elementsCountCoreBoxMinor = self._annotationElementsCountsCoreBoxMinor[i]
break
i += 1

i = 0
for layoutAnnotationGroup in self._layoutAnnotationGroups:
if i >= len(self._annotationCoreBoundaryScalingMode):
break
if self._annotationCoreBoundaryScalingMode[i] > 0:
if networkSegment.hasLayoutElementsInMeshGroup(
layoutAnnotationGroup.getMeshGroup(self._layoutMesh)):
coreBoundaryScalingMode = self._annotationCoreBoundaryScalingMode[i]
break
i += 1
return TubeNetworkMeshSegment(networkSegment, pathParametersList, elementsCountAround,
self._elementsCountThroughShell, self._isCore, elementsCountCoreBoxMinor,
self._elementsCountTransition)
self._elementsCountTransition, coreBoundaryScalingMode)

def createJunction(self, inSegments, outSegments):
"""
Expand Down

0 comments on commit 47b22e5

Please sign in to comment.