FreeCAD compatible Fenics mshr geometry building

Fenics geometry building API is based on CGAL, it can import CAD, but can not export CAD file like step. Here is a solution, not fully tested. This file may be further developed in github.

#FreeCAD path
import sys
sys.path.append("/usr/lib/freecad/lib")
import FreeCAD
App = FreeCAD
import Part
FreeCAD.ParamGet("User parameter:BaseApp/Preferences/Units/").SetInt("UserSchema", 1) # MKS ==0
# standard length should be in metre

# to be more compatible with Fenics mshr package
from FreeCAD import Vector as Point 
def Box(p1, p2):  # corners points
    box = App.ActiveDocument.addObject("Part::Box", 'Box')
    #new_center = (p2+p1).scale(0.5, 0.5, 0.5)
    new_pos = p1  # assume each component of p2 > p1
    box.Placement.move( new_pos )
    dim = p2-p1
    box.Length = math.fabs(dim[0])
    box.Width = math.fabs(dim[1])
    box.Height = math.fabs(dim[2])
    return box


def Cylinder(p1, p2, r1, r2):
    c = App.ActiveDocument.addObject("Part::Cylinder", 'Cylinder')
    #new_center = (p2+p1).scale(0.5, 0.5, 0.5)
    new_pos = p1  # assume each component of p2 > p1
    scaler = 1.0/(p2-p1).Length
    axis_vector = (p2-p1).scale(scaler, scaler, scaler)
    c.Placement = App.Placement(new_pos, App.Rotation(axis_vector, 0))  # axis setting does not work
    #c.Placement.move( new_pos )
    #c.Placement.move( new_pos )
    c.Height = (p2-p1).Length
    c.Radius = r1
    return c

def CSGUnion(base, tool):
    out = App.activeDocument().addObject("Part::Fuse","Union")
    out.Base = base
    out.Tool = tool
    return out

def CSGIntersection(base, tool):
    out = App.activeDocument().addObject("Part::Common","Intersecton")
    out.Base = base
    out.Tool = tool
    return out

def CSGDifference(base, tool):
    out = App.activeDocument().addObject("Part::Cut","Cut")
    out.Base = base
    out.Tool = tool
    return out

##################################################
def test():
    box1 = Box(Point(0, 0, 0), Point(10, 20, 30))
    cylinder = Cylinder(Point(10, 0, 25), Point(20, 0, 25), 5, 5)
    box2 = Box(Point(05, 15, 25), Point(15, 25, 35))
    box = CSGUnion(box1, box2)
    return CSGUnion(box, cylinder)
    
def export_geom(docname):
    App.newDocument(docname)
    App.setActiveDocument(docname)
    App.ActiveDocument=App.getDocument(docname)

    goem = test()
    App.ActiveDocument.recompute()
    __objs__ = []

    __objs__.append( goem )
    Part.export(__objs__,"{}.step".format(docname))

    App.ActiveDocument.saveAs("/tmp/{}.fcstd".format(docname))

def main():
    export_geom("myGeom")

if __name__ == '__main__':
     main()