`volume_factory.loft` stl bug

adamvas1 commented 3 years ago

Running the example I found in the documentation, it seems that the volume_factory.loft function has some kind of problem. The produced stl has a weird thing protruding. I assume it is some kind of intersection happening. This makes problems when trying to process the output with cgal.

example:

using: Splipy 1.4.0 Scipy1.5.4 Numpy 1.19.4

for reference the code I ran

import matplotlib.pyplot as plt
import numpy as np
from splipy import surface_factory, volume_factory
from splipy.io.stl import STL
import os 

srf1 = surface_factory.disc(r=1)
srf2 = 1.5 * srf1 + [0,0,1]
srf3 = 1.0 * srf1 + [0,0,2]
srf4 = 1.3 * srf1 + [0,0,4]
vol = volume_factory.loft(srf1, srf2, srf3, srf4)

# alternatively you can provide all input curves as a list
#all_my_surfaces = [srf1, srf2, srf3, srf4]
#vol = volume_factory.loft(all_my_surfaces)

SAVE = True
if SAVE:
    with STL('test.stl') as myfile:
        myfile.write(vol.swap()) # swap() is to make sure normals are pointing out of the object

os.system('fstl ./test.stl')

TheBB commented 3 years ago

Thanks for the report. I think the problem is in the STL writer, not the loft function. I don't see the same problem in other formats.

adamvas1 commented 3 years ago

Don't know if this helps but here is what freecad analysis reports

VikingScientist commented 3 years ago

Good catch. It would seem the bug is in the extraction of edges from the volumetric data: You model this as a solid trivariate (volumetric) spline, but STL does not support this representation, so it extracts boundary-information automatically for you and this is what is stored in the stl file. It is this particular logic which fails when dealing with periodic (i.e. the circluar disc) representations

This is for sure a bug and it will be fixed within the next release on PyPi. For now I propose three workarounds, and you can take your pick:

1 Don't use periodic (polar) representations

srf1 = surface_factory.disc(r=1, type='square')

2 Break the periodicity

vol = vol.split(0, 'v')

Note that this is still going to create the 'sem' as two overlapping surfaces on the inside of your object. Mesh analysis tools or 3D printing software is not going to be very happy with this.

3 Manually extract the relevant faces yourself

faces = list(vol.faces()) # convert from tuple to list so we can make some changes
del faces[2:4] # drop the faces corresponding to the sem
del faces[0] # drop the faces corresponding the center pole (line going up along the z-axis with area=0) 
for f in faces:
    myfile.write(f)

This also allows you to manually reorient individual faces to make sure all normals point outward. With this particular setup you will need to swap face[0] and face[2] (indexing refer to after deleting things from the list of faces).

SINTEF / Splipy