finite element method – 3D FEM Meshing Internal Regions


I am trying to set up a 3D element mesh with intersecting regions having different mesh densities. I am having difficulty setting up the defining boundary meshes from which I will then apply ToElementMesh. I understand how to do it in 2D but I do not know the best way to do it for 3D. My code below has been cut down to try and show the basic problem I have. I need to set up the boundary mesh on the green problem volume so the intersections with the “e-core” region on the x=z=0 axis can be meshed consistent with the finer mesh to be used in the e-core region volume. Although I have shown the full core, due to symmetry in the problem I will only use 1/4 of it, i.e, that which intersects with the green volume.

Please note I only have MM 10.4, so I do not have access to FEMAddons. However, I would also be interested to see how it could be done if I upgraded in the future.

Clear("Global`*");
Needs("NDSolve`FEM`");

eCore(cw_, ch_, cd_, ww_, wh_) := 
  Module((*cw = core width, ch = core height, cd = core depth, www = 
   window width, w = window height*){vertices, topFace, reg},
   vertices = {{-cw/2, 0}, {-cw/4 - ww/2, 0}, {-cw/4 - ww/2, 
      wh}, {-cw/4 + ww/2, wh}, {-cw/4 + ww/2, 0}, {cw/4 - ww/2, 
      0}, {cw/4 - ww/2, wh}, {cw/4 + ww/2, wh}, {cw/4 + ww/2, 
      0}, {cw/2, 0}, {cw/2, ch}, {-cw/2, ch}};
   topFace = 
    BoundaryMeshRegion(vertices, 
     Line({1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 1}));
   reg = RegionProduct(topFace, 
     MeshRegion({{-ch/2}, {ch/2}}, Line({1, 2}))); reg);

(*Create an e-core using above function and rotate/translate position 
as required*)
regCore1 = 
  TransformedRegion(
   TransformedRegion(eCore(0.065, 0.033, .027, .013, .022), 
    RotationTransform(0, {0, 0, 1})), 
   TranslationTransform({0, 0.002, 0})) ;
bmeshCore1 = 
  BoundaryDiscretizeRegion(regCore1, 
   MaxCellMeasure -> {"Length" -> 0.005}, Axes -> True, 
   AxesLabel -> {x, y, z});
(*get coordinates of 1/4 core1 mesh in problem volume*)
core1Coord = 
  Cases(DeleteDuplicates(MeshCoordinates(bmeshCore1)), {x_, y_, z_} /;
     x (GreaterSlantEqual) 0 && z (LessSlantEqual) 0);

(*Create air region that defines the problem boundaries allowing for 
symmetry in the problem*)
radiusAir = 0.15;
regAir1 = 
  RegionIntersection(
   Cuboid({0, 0, -radiusAir}, {radiusAir, radiusAir, 0}), 
   Ball({0, 0, 0}, radiusAir));
bmeshAir1 = 
  BoundaryDiscretizeRegion(regAir1, 
   MaxCellMeasure -> {"Length" -> 0.01}, Axes -> True, 
   AxesLabel -> {x, y, z});
RegionPlot3D({regCore1, regAir1}, Axes -> True, 
 AxesLabel -> {x, y, z}, PlotStyle -> {Blue, Green})

Region Plot

I guess I want the 3D equivalent of the Wolfram 2D example given under Element Mesh Generation. Here I have modified it to have a higher mesh density on the internal line boundary.

(*2D Example of open line boundary within a closed rectangular 
boundary - modified from Wolfram FEM Meshing example*)
n = 20; 
lineCoord = 
 DeleteDuplicates(
  Join(Table({1/6. + (i - 1)*4/(6.*(n - 1)), 1/6.}, {i, 1, n}), 
   Table({5/6., 1/6. + (i - 1)*4/(6.*(n - 1))}, {i, 1, n})));
bmesh = ToBoundaryMesh(
   "Coordinates" -> Join({{0, 0}, {1, 0}, {1, 1}, {0, 1}}, lineCoord),
    "BoundaryElements" -> {LineElement({{1, 2}, {2, 3}, {3, 4}, {4, 
        1}}), LineElement(
      Partition(Delete(Last(FindShortestTour(lineCoord)), 1), 2, 1) + 
       4)});
mesh = ToElementMesh(bmesh, MaxCellMeasure -> {"Length" -> 0.5});
mesh("Wireframe")

2D Mesh Example

Any help would be much appreciated.