Initial implementation of chtMultiRegionSimpleFoam, conjugate heat transfer support

CfdOF/CfdTools.py

@@ -166,15 +166,24 @@ def getInitialisationZoneObjects(analysis_object):
 
 def getZoneObjects(analysis_object):
     from CfdOF.Solve.CfdZone import CfdZone
-    return [i for i in analysis_object.Group if isinstance(i.Proxy, CfdZone)]
+    return [i for i in analysis_object.Group if isinstance(getattr(i, 'Proxy', None), CfdZone)]
 
 
 def getInitialConditions(analysis_object):
     return getModelOfType(analysis_object, 'CfdInitialVariables')
 
 
 def getMaterials(analysis_object):
-    return [i for i in analysis_object.Group if i.isDerivedFrom('App::MaterialObjectPython')]
+    from CfdOF.Solve.CfdSolidMaterial import CfdSolidMaterial
+    return [i for i in analysis_object.Group
+            if i.isDerivedFrom('App::MaterialObjectPython')
+            and not isinstance(getattr(i, 'Proxy', None), CfdSolidMaterial)]
+
+
+def getSolidMaterials(analysis_object):
+    from CfdOF.Solve.CfdSolidMaterial import CfdSolidMaterial
+    return [i for i in analysis_object.Group
+            if isinstance(getattr(i, 'Proxy', None), CfdSolidMaterial)]
 
 
 def getSolver(analysis_object):

CfdOF/Mesh/CfdMeshTools.py

@@ -35,6 +35,99 @@
 import Part
 
 
+
+def _get_interior_point(solid):
+    """Return a point guaranteed to be inside the solid material.
+
+    CenterOfMass can lie outside the material for non-convex solids (e.g. an air
+    shell surrounding an object, or a heatsink whose centroid falls between fins).
+    Falls back to face centroids (always on the surface, which passes isInside
+    with allowFace=True), then to a 3-D bounding-box grid.
+    """
+    tol = 1e-3
+    centroid = solid.CenterOfMass
+    if solid.isInside(centroid, tol, True):
+        return centroid
+    # Face centroids are guaranteed to lie on the solid surface and always pass
+    # isInside(allowFace=True).  This handles thin-feature solids (e.g. heatsink
+    # fins) where the centroid and a coarse grid both land outside the material.
+    for face in solid.Faces:
+        pt = face.CenterOfMass
+        if solid.isInside(pt, tol, True):
+            return pt
+    bb = solid.BoundBox
+    for xi in (0.1, 0.3, 0.5, 0.7, 0.9):
+        for yi in (0.1, 0.3, 0.5, 0.7, 0.9):
+            for zi in (0.1, 0.3, 0.5, 0.7, 0.9):
+                pt = FreeCAD.Vector(
+                    bb.XMin + xi * bb.XLength,
+                    bb.YMin + yi * bb.YLength,
+                    bb.ZMin + zi * bb.ZLength,
+                )
+                if solid.isInside(pt, tol, True):
+                    return pt
+    return centroid
+
+
+def _getCompoundLinks(part_obj):
+    """Return the list of child body objects from a compound or BooleanFragments shape object."""
+    if hasattr(part_obj, 'Links'):      # Part::Compound
+        return list(part_obj.Links)
+    if hasattr(part_obj, 'Shapes'):     # Part::BooleanFragments (older FreeCAD)
+        return list(part_obj.Shapes)
+    if hasattr(part_obj, 'Objects'):    # Part::BooleanFragments (FreeCAD 1.x)
+        return list(part_obj.Objects)
+    return [part_obj]
+
+
+def _get_shape_ref_solids(shape_ref):
+    """Return solid shapes from a FreeCAD LinkSub reference.
+
+    Solid materials can reference selected sub-solids in a compound/body.  Use
+    those sub-elements when present instead of broadening the selection to the
+    whole source object, otherwise unrelated volumes can be assigned to the
+    solid region after BooleanFragments.
+    """
+    obj = shape_ref[0]
+    sub_names = shape_ref[1] if len(shape_ref) > 1 else []
+    solids = []
+    for sub_name in sub_names:
+        try:
+            sub_shape = obj.Shape.getElement(sub_name)
+        except (Part.OCCError, ValueError, IndexError, AttributeError):
+            continue
+        solids.extend(sub_shape.Solids if sub_shape.Solids else [sub_shape])
+    if solids:
+        return solids
+    return list(obj.Shape.Solids if obj.Shape.Solids else [obj.Shape])
+
+
+def _solid_belongs_to_reference(result_solid, reference_solids):
+    """Return True if a BooleanFragments result solid belongs to a source solid.
+
+    The normal test is containment of an interior point, but OCC can occasionally
+    report false for solids that came from BooleanFragments.  In that case, fall
+    back to the boolean common volume, which is slower but much more reliable for
+    assigning CHT physical volumes.
+    """
+    interior = _get_interior_point(result_solid)
+    for ref_solid in reference_solids:
+        if ref_solid.isInside(interior, 1e-3, False):
+            return True
+
+    result_volume = abs(result_solid.Volume)
+    if result_volume <= 0:
+        return False
+    for ref_solid in reference_solids:
+        try:
+            common_volume = abs(result_solid.common(ref_solid).Volume)
+        except Part.OCCError:
+            continue
+        if common_volume / result_volume > 0.5:
+            return True
+    return False
+
+
 class CfdMeshTools:
     def __init__(self, cart_mesh_obj):
         self.mesh_obj = cart_mesh_obj
@@ -687,6 +780,52 @@ def writeMeshCase(self):
             self.gmsh_settings['ClMin'] = self.clmin
             sols = (''.join((str(n+1) + ', ') for n in range(len(self.mesh_obj.Part.Shape.Solids)))).rstrip(', ')
             self.gmsh_settings['Solids'] = sols
+
+            # Build per-region volume map for CHT (chtMultiRegionSimpleFoam/chtMultiRegionFoam)
+            solid_material_objs = CfdTools.getSolidMaterials(self.analysis)
+            if solid_material_objs:
+                region_volume_map = {}
+                mesh_part = self.mesh_obj.Part
+                links = _getCompoundLinks(mesh_part)
+                result_solids = mesh_part.Shape.Solids
+                solid_vol_indices_all = set()
+                # For each CfdSolidMaterial, find the result solid volumes that belong
+                # to it by geometric containment/overlap in its referenced shapes. Each
+                # material is handled independently so that multiple solid regions
+                # receive distinct volume index sets rather than the same combined set.
+                for solid_obj in solid_material_objs:
+                    rname = getattr(solid_obj, 'RegionName', '') or solid_obj.Label
+                    body_names = set()
+                    body_solids = {}
+                    for ref in solid_obj.ShapeRefs:
+                        body_names.add(ref[0].Name)
+                        body_solids.setdefault(ref[0].Name, []).extend(_get_shape_ref_solids(ref))
+                    for lnk in links:
+                        if lnk.Name in body_names and lnk.Name not in body_solids:
+                            body_solids[lnk.Name] = list(lnk.Shape.Solids if lnk.Shape.Solids else [lnk.Shape])
+                    flat_solids = [solid for solids in body_solids.values() for solid in solids]
+                    vol_indices = [i for i, s in enumerate(result_solids, start=1)
+                                   if _solid_belongs_to_reference(s, flat_solids)]
+                    solid_vol_indices_all.update(vol_indices)
+                    if vol_indices:
+                        region_volume_map[rname] = ', '.join(str(v) for v in vol_indices)
+                # Fluid volumes are those not claimed by any solid material
+                fluid_vol_indices = [i for i in range(1, len(result_solids) + 1)
+                                     if i not in solid_vol_indices_all]
+                mr_objs = CfdTools.getMeshRefinementObjs(self.mesh_obj)
+                internal_zones = [o for o in mr_objs if o.Internal]
+                if internal_zones:
+                    fluid_rname = internal_zones[0].Label
+                else:
+                    fluid_mats = CfdTools.getMaterials(self.analysis)
+                    fluid_rname = fluid_mats[0].Label if fluid_mats else 'fluid'
+                if fluid_vol_indices:
+                    region_volume_map[fluid_rname] = ', '.join(str(v) for v in fluid_vol_indices)
+                self.gmsh_settings['IsMultiRegion'] = True
+                self.gmsh_settings['RegionVolumeMap'] = region_volume_map
+            else:
+                self.gmsh_settings['IsMultiRegion'] = False
+
             self.gmsh_settings['BoundaryFaceMap'] = {}
             for k in range(len(self.patch_faces)):
                 for l in range(len(self.patch_faces[k])):