ElementMesh (rendering?) issue

Question

The appearance of the following ElementMesh is bad as you can see for example on the inside surface of the shell. The same problem happens converting to a MeshRegion (conversion does not warn about problems).

Needs["NDSolve`FEM`"]
em = Import[
  "https://www.dropbox.com/s/wc9nd0kyf55t9bg/ElementMeshRenderingIssue.mx?dl=1"]
em["Wireframe"["MeshElement" -> "MeshElements", 
  "ElementMeshDirective" -> 
   Directive[Specularity[.2], EdgeForm@[email protected], 
    FaceForm@[email protected]], 
  Lighting -> {{"Ambient", GrayLevel[0.45]}, {"Directional", 
     GrayLevel[0.3], ImageScaled[{2, 0, 2}]}, {"Directional", 
     GrayLevel[0.33], ImageScaled[{2, 2, 2}]}, {"Directional", 
     GrayLevel[0.3], ImageScaled[{0, 2, 2}]}}]]

I don't understand if this is just a rendering problem or there is also a problem in the underlying mesh.

I have built this mesh by myself. As I stated in another question the vertices of the QuadElement are not exactly coplanar, and I don't know how to adjust the coordinates so that they eventually become coplanar.

I noticed some strange behavior using this mesh to solve a PDE with FEM: this can be another clue of some problem with the mesh.

To summarize I want to try to understand

1. if the visualization problem is caused by some problem with the mesh construction
2. if the problem is related with the non-exactly-planar faces, if there is some way to adjust

Any idea on how to proceed? Thanks

---

After more investigation, and thanks to @user21, I discovered that the problem of the failed ToBoundaryMesh/ToElementMesh roundtrip is probably related to the fact that the mesh uploaded and shown here is no a "complete" spherical shell (this is what I want for a figure to illustrate the process of mesh creation).

The mesh I use as a domain for my PDE is a complete spherical shell, so this should not be the reason why I get strange result when solving over this mesh.

I made another simpler mesh with the same algorith, starting with just 1 face of the cubed sphere (i.e. to get 1/6 of the spherical shell), and with a coarse grid.

mesh = Import[
  "https://www.dropbox.com/s/46e93mbacvgn2i6/\
ElementMeshRenderingIssue2.mx?dl=1"]

I also made a function to extract only part of this mesh.

ElementMeshExtract[mesh_, mkform : {memkform_, bemkform_, pemkform_}, 
  idform : {meidform_, beidform_, peidform_}] :=
 Module[{elems},
  elems = MapThread[
    Cases[#1, type_[el_, mk_] :> With[{sel = MatchQ[#2] /@ mk},
        If[! Or @@ sel, Nothing, type[Pick[el, sel], Pick[mk, sel]]]]
      , {1}] &, {mesh /@ {"MeshElements", "BoundaryElements", 
       "PointElements"}, mkform}];
  elems = MapThread[
    Cases[#1, 
      type_[el_, mk_] :> With[{sel = MatchQ[#2] /@ Range@Length@el},
        If[! Or @@ sel, Nothing, type[Pick[el, sel], Pick[mk, sel]]]]
      , {1}] &, {elems, idform}];
  incidents = Union@Flatten@DeleteCases[elems, {}][[All, 1, 1]];
  assoc = AssociationThread[incidents, Range@Length@incidents];
  elems[[All, All, 1]] = 
   Map[Lookup[assoc, #] &, elems[[All, All, 1]], {-2}];
  ToElementMesh[
   "Coordinates" -> mesh["Coordinates"][[Keys@assoc]],
   "MeshElements" -> elems[[1]],
   "BoundaryElements" -> elems[[2]],
   "PointElements" -> elems[[3]]
   ]
  ]

ElementMeshExtract[mesh_, mkform : {_, _, _}, idform_] := 
 ElementMeshExtract[mesh, mkform, {idform, idform, idform}]
ElementMeshExtract[mesh_, mkform_, idform : {_, _, _}] := 
 ElementMeshExtract[mesh, {mkform, mkform, mkform}, idform]
ElementMeshExtract[mesh_, mkform_, idform_] := 
 ElementMeshExtract[
  mesh, {mkform, mkform, mkform}, {idform, idform, idform}]
ElementMeshExtract[mesh_, mkform_] := 
 ElementMeshExtract[mesh, mkform, _]

Putting all together in a Manipulate I still cannot find anything wrong with that mesh that could explain the bad rendering or the strange results od NDSolve.

Manipulate[(
  lmesh = 
   ElementMeshExtract[
    mesh, {Alternatives @@ Range@layers, Except@_, Except@_}];
  lc = Total[Length@*First /@ lmesh["MeshElements"]];
  pmesh = 
   ElementMeshExtract[lmesh, _, {i_ /; 1 <= i <= progress, _, _}];
  Column[{
    Show[
     MeshRegion[pmesh, Method -> {"CoplanarityTolerance" -> -100}],
     ImageSize -> Medium, PlotRange -> lmesh["Bounds"], 
     Method -> {"ShrinkWrap" -> False}],
    BoxWhiskerChart[ElementMeshQuality[pmesh],
     BarOrigin -> Left, BarSpacing -> $MachineEpsilon,
     PlotLabel -> "Quality", AspectRatio -> 1/4, 
     PlotRange -> {{0.5, 1}, Automatic},
     Frame -> {{False, False}, {True, True}}, 
     GridLines -> {Automatic, None},
     ImageSize -> Small]
    }, Center]
  ),
 {{lmesh, lmesh}, None}, {{pmesh, pmesh}, None}, {{lc, lc}, None},
 {{layers, 2}, 1, 5, 1, 
  Appearance -> "Labeled"}, {{progress, Floor[lc/2]}, 1, lc, 1, 
  Appearance -> "Labeled"},
 ControlPlacement -> Bottom
 ]