Unexpected behaviour when splitting a rotated cylinder

[jannisEspelage]

Hi all,

i want to split a rotate cylinder (Solid) at a certain height. I do so by adding the cylinder to a Workplane. Then i use the transformed function with the offset argument to create a plane where i want to split the cylinder. However, i experienced an unexpected behaviour. For a certain cylinder and rotations about y and x axis (by ry and rx angles) the split function seems to return the whole cylinder as Top AND Bottom instead of properly splitting it (the splitting plane should intersect the cylinder for this case).

This behaviour can be recreated by the below minimum example which results in the following output on my machine:

Python version: 3.11.9 | packaged by Anaconda, Inc. | (main, Apr 19 2024, 16:40:41) [MSC v.1916 64 bit (AMD64)]
Cadquery version: 2.4.0
Numpy version: 1.26.4
Volumes are NOT equal!
[3500. 3500. 3500. 7000.]

import sys
import cadquery as cq
import numpy as np

print("Python version:", sys.version)
print("Cadquery version:", cq.__version__)
print("Numpy version:", np.__version__)

# cylinder properties
r = 9.8
h = 11.6
h_split = 3.265
# base vectors for roation
e0 = (0, 0, 0)
ex = (1, 0, 0)
ey = (0, 1, 0)
# rotations
rx1 = np.degrees(1.0e-6)                # this rotation results in expected behaviour
rx2 = np.degrees(1.490e-8)              # this rotation results in expected behaviour
rx3 = np.degrees(1.490116e-8)           # this rotation results in unexpected behaviour
ry = np.degrees(0.14009)

t0 = cq.Solid.makeCylinder(radius=r, height=h)
v0 = t0.Volume()
# rotate around y axis by ry > rotate about x axis by rx1
t1 = t0.rotate(e0, ey, ry).rotate(e0, ex, rx1)
# rotate around y axis by ry > rotate about x axis by rx2
t2 = t0.rotate(e0, ey, ry).rotate(e0, ex, rx2)
# rotate around y axis by ry > rotate about x axis by rx3
t3 = t0.rotate(e0, ey, ry).rotate(e0, ex, rx3)

# split t1 at h_split and return both upper and lower part
# convert into list of solids and create compound
# this should result in t0
t1_split = cq.Compound.makeCompound(
    cq.Workplane()
    .add(t1)
    .transformed(offset=(0, 0, h_split))
    .split(keepTop=True, keepBottom=True)
    .solids().vals()
)
v1 = t1_split.Volume()

# same for t2
t2_split = cq.Compound.makeCompound(
    cq.Workplane()
    .add(t2).transformed(offset=(0, 0, h_split))
    .split(keepTop=True, keepBottom=True)
    .solids().vals()
)
v2 = t2_split.Volume()

# same for t3
t3_split = cq.Compound.makeCompound(
    cq.Workplane()
    .add(t3).transformed(offset=(0, 0, h_split))
    .split(keepTop=True, keepBottom=True)
    .solids().vals()
)
v3 = t3_split.Volume()

# compare volumes of t0, t1, t2 and t3
if not all(np.isclose([v1, v2, v3], v0)):
    print("Volumes are NOT equal!")
    print(np.round([v0, v1, v2, v3], 0))
else:
    print("All fine")

[lorenzncode]

I reproduced the issue. It appears to be an OCCT kernel issue.

import cadquery as cq

r = 9.8
h = 11.6
h_split = 3.265

cyl1 = cq.Solid.makeCylinder(r, h)
cyl1 = cyl1.rotate((0, 0, 0), (0, 1, 0), 8.03)

# invalid shape unless commented out:
cyl1 = cyl1.rotate((0, 0, 0), (1, 0, 0), 8.54e-7)

result = (
    cq.Workplane()
    .add(cyl1)
    .transformed(offset=(0, 0, h_split))
    .split(keepTop=True, keepBottom=True)
)

result_1 = result.item(0)
result_2 = result.item(1)


print(result_1.val().isValid())
print(result_2.val().isValid())

or similar result with Free function API:

from cadquery.occ_impl.shapes import *

r = 9.8
h = 11.6
h_split = 3.265

cyl1 = cylinder(r * 2, h)
cyl1 = cyl1.rotate((0, 0, 0), (0, 1, 0), 8.03)

# invalid shape unless commented out:
cyl1 = cyl1.rotate((0, 0, 0), (1, 0, 0), 8.54e-7)


tool1 = face(rect(30, 30))
tool1 = tool1.translate((0, 0, h_split))
result = split(cyl1, tool1)

print(result.isValid())
result_1 = result.Solids()[0]
result_2 = result.Solids()[1]

Reproduced with OCCT Draw (tested with 7.8.0, 7.8.1):

pload MODELING VISUALIZATION

dset r 9.8
dset h 11.6 
dset h_split 3.265

pcylinder cyl1 r h
trotate cyl1 0 0 0 0 1 0 8.03

# invalid shape unless commented out:
trotate cyl1 0 0 0 1 0 0 8.54e-7

polyline plane1 0 0 0 30 0 0 30 30 0 0 30 0 0 0 0 
mkplane plane1 plane1 
ttranslate plane1 -15 -15 h_split

baddobjects cyl1 plane1
bapisplit result
explode result So

puts [vprops result_1]
puts [vprops result_2]
puts [checkshape result_1 -short]

vinit
vdisplay -dispMode 1 result
vfit

[jannisEspelage]

Hi @lorenzncode,

thank you for confirming this behaviour! Should this be issued in OCCT so they can fix this?