*Not an answer, but too long for a comment.*

I don't know exactly about the details of the meshing algorithm that _Mathematica_ uses. But I expect that it generates a region whose boundary vertices all lie on the true boundary surface (up to a very small error due to solving the boundary equations iteratively). As the body at hand is convex (it is the intersection of convex cylinders), and as all triangles/tetrahedra have their vertices inside the body (or on its boundary),
this would imply that the discretized body lies completely in the _inside_ of the true body. Hence _Mathematica_'s result for the volume should be _lower_ than the true volume (plus some rounding errors). 
Because `4.40045` is greater than `22/5` (and since I deem `0.00045` as too much for a rounding error), I would say that the 22/5 conjecture is not true.

**Edit**

When you use `DiscretizeRegion` and `BoundaryDiscretizeRegion` one should prescibe a maximal edge length of the generated tets and triangles. For example like this:

    M = BoundaryDiscretizeRegion[reg, MaxCellMeasure -> (1 -> 0.01)];
    Volume[M]

> 4.40038

If plot `M` and zoom-in, then you will see that the sharp edges of the body are cut off. Thus, the actual volume should be larger than `4.40038`.


**Edit 2**

Here is a function that can be used to check whether a vertex lies within the body. It returns `1` if the 3-vector `X` lies inside and `0` otherwise.

    cf = Compile[{{X, _Real, 1}},
       Block[{x, y, z},
        x = Compile`GetElement[X, 1];
        y = Compile`GetElement[X, 2];
        z = Compile`GetElement[X, 3];
        
        Boole[
         (0.3333333333333333` x + 0.6666666666666666` y - 
              0.3333333333333333` z)^2 + (0.6666666666666666` x + 
              0.3333333333333333` y + 
              0.3333333333333333` z)^2 + (0.3333333333333333` x - 
              0.3333333333333333` y + 0.6666666666666666` z)^2 <= 
           1.` && (0.6666666666666666` x - 0.3333333333333333` y - 
              0.3333333333333333` z)^2 + (-0.3333333333333333` x + 
              0.6666666666666666` y - 
              0.3333333333333333` z)^2 + (-0.3333333333333333` x - 
              0.3333333333333333` y + 0.6666666666666666` z)^2 <= 
           1.` && (0.6666666666666666` x + 0.3333333333333333` y - 
              0.3333333333333333` z)^2 + (0.3333333333333333` x + 
              0.6666666666666666` y + 
              0.3333333333333333` z)^2 + (-0.3333333333333333` x + 
              0.3333333333333333` y + 0.6666666666666666` z)^2 <= 
           1.` && (0.6666666666666666` x - 0.3333333333333333` y + 
              0.3333333333333333` z)^2 + (-0.3333333333333333` x + 
              0.6666666666666666` y + 
              0.3333333333333333` z)^2 + (0.3333333333333333` x + 
              0.3333333333333333` y + 0.6666666666666666` z)^2 <= 1.` && 
          x^2 + y^2 <= 1.`
         ]
        ],
       CompilationTarget -> "C",
       RuntimeAttributes -> {Listable},
       Parallelization -> True
       ];

If I apply this to the mesh `M` from above then I get

    b = cf[MeshCoordinates[M]];
    Count[b, 0]
    Count[b, 1]

> 109

> 399393

So the vast majority of points lies inside. But few of them loe outside. I _guess_ they deviate so little that this would not matter. Maybe one can move them just a little inside? For simplicity, I just move _all_ the points of the discrete region a little bit towards the origin:

    Mscaled = TransformedRegion[M, {x, y, z} |-> 0.99999999 {x, y, z}];
    Count[cf[MeshCoordinates[Mscaled]], 0]
    Volume[Mscaled]

> 0

> 4.40038

So now all of the discrete regions vertices lies within the convex body, but the volume enclosed by the mesh is still greater than 22/5.


**Edit 3**

By the way, in order to have at least on higher resolution picture of the body:

[![enter image description here][1]][1]


  [1]: https://i.sstatic.net/3Udl9l.png