Union of touching polyhedra

edit

Let me first define your polyhedron over the same base ring:

sage: cube_a = Polyhedron(vertices=[(0,0,0), (0,0,1), (0,1,0), (0,1,1), (1,0,0), (1,0,1), (1,1,0), (1,1,1)], base_ring=QQ) 
sage: cube_b = Polyhedron(vertices=[(1,1/2,0), (1,1/2,1), (1,3/2,0), (1,3/2,1), (2,1/2,0), (2,1/2,1), (2,3/2,0), (2,3/2,1)], base_ring=QQ)                                                          
sage: cube_a                                                                                         
A 3-dimensional polyhedron in QQ^3 defined as the convex hull of 8 vertices
sage: cube_b                                                                                         
A 3-dimensional polyhedron in QQ^3 defined as the convex hull of 8 vertices

You may define this function which returns the convex hull of a list of polyhedron:

def convex_hull(list_of_polyhedron): 
    V = sum((p.vertices() for p in list_of_polyhedron), tuple()) 
    return Polyhedron(V)

In your case, you obtain:

sage: convex_hull([cube_a, cube_b])
A 3-dimensional polyhedron in QQ^3 defined as the convex hull of 12 vertices

If the union is convex, the convex hull is the union. In your case, the convex hull does not seem to be equal to the union of the two cubes. Comparing their volume, we obtain:

sage: cube_a.volume()                                                                                
1
sage: cube_b.volume()                                                                                
1
sage: convex_hull([cube_a, cube_b]).volume()                                                         
5/2

Recently, I needed to check whether the union of a list of polyhedron is convex before creating their union using convex hull. I wrote that in my package which can be installed with sage -pip install slabbe. In your case, it returns False, so the convex hull is not equal to the union of the two polyhedron:

sage: from slabbe.polyhedron_partition import is_union_convex   
sage: is_union_convex([cube_a, cube_b])                                                              
False

Polyhedral complex may be a way to go - e.g., try:

P = PolyhedralComplex( [cube_a, cube_b] )
P.plot()