Not an answer, just a illustration of what I mean in my comment above :

Unfortunately, the OpenCascadeLink is a relatively new addition and I do not believe there is any way now to pre-assign element markers to boundaries. Currently, Mathematica will do a good job preserving continuity across internal surfaces, but there is limited support at assign boundary conditions aside for DirichletConditions. Additionally, the default boundary condition is a zero flux NeumannValue. For a wide variety of physics problems and symmetry conditions, zero flux applies. Therefore, if one excludes zero flux and symmetry conditions, often that leaves only a few surfaces left to assign.

I too, found no great way to display boundary ElementMarkers in 3D, so I wrote the following manipulate code to display the ElementMarkers one-at-a-time. Perhaps, you will find it useful.

Manipulate[Show[{bmesh["Edgeframe"],
     bmesh[
      "Wireframe"[ElementMarker == #, 
       "MeshElementStyle" -> FaceForm[colors[[#]]]]]}, 
    PlotLabel -> 
     Style[StringTemplate["ElementMarker = ``"][#], 18, 
      colors[[#, 2]]]] &[group], {group, First@groups, Last@groups, 
  1}, ControlPlacement -> Top]

Update 1: Using CheckBox control

To view multiple surface at once, it is convenient to use a CheckboxBar control like so:

surfaces = 
  AssociationThread[groups, 
   bmesh["Wireframe"[ElementMarker == #, 
       "MeshElementStyle" -> FaceForm[colors[[#]]]]] & /@ groups];
Manipulate[Show[{bmesh["Edgeframe"],
   choices /. surfaces}], {{choices, groups}, groups, CheckboxBar}, 
 ControlPlacement -> Top]

The following displays all the internal surfaces:

Update 2: Surface Identification Using Annotations

Here is a technique to identify a surface group with a PlotLabel by hovering over a GraphicsComplex that represents the group ID.

First, we will define some shorthands to extract mesh information:

(* Shorthand functions to extract mesh info *)
ebif = ElementIncidents[#["BoundaryElements"]][[1]] &;
ebmf = ElementMarkers[#["BoundaryElements"]][[1]] &;
eif = ElementIncidents[#["MeshElements"]][[1]] &;
emf = ElementMarkers[#["MeshElements"]][[1]] &;
epif = Flatten@ElementIncidents[#["PointElements"]] &;
epmf = Flatten@ElementMarkers[#["PointElements"]] &;
UF = Union@Flatten[#, Infinity] &;

Next, we define a function, meshByBoundaryID, that generates an independent BoundaryMesh by group ID:

(* Function to Extract Boundary Mesh by Marker *)
Clear[meshByBoundaryID]
meshByBoundaryID[m_][marker_] := Module[
  {inc = ebif[m], mrk = ebmf[m], crd = m["Coordinates"], tinc, uniq, 
   newcrd, bm},
  tinc = Extract[inc, Position[mrk, marker]];
  uniq = UF@tinc;
  newcrd = crd[[UF@tinc]];
  tinc = tinc /. AssociationThread[uniq -> Range[Length@uniq]];
  bm = ToBoundaryMesh["Coordinates" -> newcrd, 
    "BoundaryElements" -> {TriangleElement[tinc, 
       ConstantArray[marker, Length@tinc]]}, "MeshOrder" -> 1];
  bm
  ]

(* Create instance of function based on grouped BoundaryMesh *)
bm = meshByBoundaryID[bmesh];

Next, we create the visualization:

(* BarLegend of groups *)
legendBar = 
  BarLegend[{"BrightBands", MinMax[groups]}, (Length[groups] - 2), 
   LegendLabel -> Style["group", Opacity[0.6`]]];
(* Convert Individual bmeshes into GraphicsComplex for Annotation *)
gcs = Graphics3D[{Directive[EdgeForm[Black], FaceForm[colors[[#]]]], 
      Annotation[ElementMeshToGraphicsComplex[bm[#]], 
       Style[StringTemplate["ElementMarker = ``"][#], 18, 
        colors[[#]]], "Mouse"]}, 
     PlotLabel -> Dynamic[MouseAnnotation[""]], Boxed -> False] & /@ 
   groups;
(* Surface Association for CheckBoxBar *)
surfaces = AssociationThread[groups, gcs[[#]] & /@ groups];
(* Visualization *)
Manipulate[Legended[Show[{
    {choices /. surfaces}}], legendBar], {{choices, groups}, groups, 
  CheckboxBar}, ControlPlacement -> Top]

The animation augmented with Camtasia annotations at the bottom shows the following workflow:

• Orient Model to view symmetry planes
• Sweep mouse over symmetry surfaces to identify group IDs (1, 2, 3, 4, 6)
• Toggle visibility of (1, 2, 3, 4, 6)
• Identify Curved Boundary (5)
• Toggle visibility of (5)
• Remaining surfaces are internal

One could envision making a more dynamic application that could toggle visibility with MouseDown events, but this simple approach took less than 20 seconds to identify key surfaces.

Unless it is explicitly documented somewhere, I don't think you can count on the ordering of internal surfaces last to be stable version-to-version.