c++ - Generating outside supporters into mesh for 3D printing

Question

Prologue

This is my attempt to re-ask the closed Generating supporters for 3D printing as it is interesting question but lacking important details ... This is intended as Q&A and currently I am working on the code for the answer but feel free to answer (I accept the best answer).

Problem description

OK here some basic info about the problem:

• Supports in 3D Printing: A technology overview

As this is a huge problem I will focus on the generic mesh/support-pattern merging geometry problem.

In a nutshell If we want to print any mesh we can do it only if it is connected to the starting plane up to angle of ~45 degrees (+/- for different printing technologies). So if we got parts that are not connected to this plane we need to create a bridge that will hold/connect it to it. Something like this (image taken from the page linked above):

Of coarse we need to add as small amount of material possible and still has it strong enough to hold our mesh in place without bending. On top of all this we need to weaken the support near the mesh so it can be break off after printing easily.

Do not forget that shape and placement is dependent on many things like material and technology used, heat flow.

Question:

To narrow this huge topic to answerable question let us focus solely on this problem:

How to merge 3D triangulated mesh (boundary representation like STL) with predefined support pattern (like 3 side prism) connecting it from defined plane perpendicularly ?

Using simple C++.

Answer 1

OK lets start with the absolute basics.

1. support shape

   You can use any shape to meet the specifics of used printing technology. The easiest to generate within STL is 3 side prism like shape which contains 2 triangular bases (top and bottom) and 3 sides all of which have 2 triangles. So 8 triangles total.

   

   This shape will start at some base plane (Z=0) and will go up until it hits the mesh. However to make this work the support must have a small gap between mesh and itself where we will add our weakened joint structure with mesh latter on.

2. support pattern

   there are a lot of options here so I chose the simplest (not foul proof however) and that is to place the supports in a uniform grid with constant distance grid between the supports.

   so simply cast a ray from each grid position on the base plane in up direction and check for intersection with mesh. If found place the support at that position with height just gap below the intersection point.

3. Joints

   The idea is to join fan of very thin supports in cone like shape connecting & covering the supported surface above main support prism with less than 45 deg angle (so the gap should be big enough to cover grid distance in such manner).

   The main problem here is that we must subdivide the triangles we are connecting to so we meet the STL mesh properties. To solve the connection problem (avoid holes or breaking connection requirements of the STL) we can use different solid for supports and different for our mesh. That will also allow us to touch surfaces without re-triangulating them making this a lot easier task.

   For simplicity I chose tetrahedron shape which is simple to construct from triangles and also present weakness at the mesh/support joint.

So let us take some test STL mesh and place it above our base plane:

and place our main supports:

and also the joints:

Here VCL/C++ code for this STL3D.h:

//---------------------------------------------------------------------------
//--- simple STL 3D mesh ----------------------------------------------------
//---------------------------------------------------------------------------
#ifndef _STL3D_h
#define _STL3D_h
//---------------------------------------------------------------------------
#ifdef ComctrlsHPP
TProgressBar *progress=NULL;        // loading progress bar for realy big STL files
#endif
void _progress_init(int n);
void _progress     (int ix);
void _progress_done();
//---------------------------------------------------------------------------
class STL3D                         // STL 3D mesh
    {                                                                      
public:
    double center[3],size[3],rmax;  // bbox center,half sizes, max(size[])
    struct _fac
        {
        float p[3][3];              // triangle vertexes CCW order
        float n[3];                 // triangle unit normal pointing out
        WORD attr;
        _fac()  {}
        _fac(_fac& a)   { *this=a; }
        ~_fac() {}
        _fac* operator = (const _fac *a) { *this=*a; return this; }
        //_fac* operator = (const _fac &a) { ...copy... return this; }
        void compute()                                  // compute normal
            {
            float a[3],b[3];
            vectorf_sub(a,p[1],p[0]);
            vectorf_sub(b,p[2],p[1]);
            vectorf_mul(n,a,b);
            vectorf_one(n,n);
            }
        double intersect_ray(double *pos,double *dir)   // return -1 or distance to triangle and unit ray intersection
            {
            double p0[3],p1[3],p2[3];                   // input triangle vertexes
            double e1[3],e2[3],pp[3],qq[3],rr[3];       // dir must be unit vector !!!
            double t,u,v,det,idet;
            // get points
            vector_ld(p0,p[0][0],p[0][1],p[0][2]);
            vector_ld(p1,p[1][0],p[1][1],p[1][2]);
            vector_ld(p2,p[2][0],p[2][1],p[2][2]);
            //compute ray triangle intersection
            vector_sub(e1,p1,p0);
            vector_sub(e2,p2,p0);
            // Calculate planes normal vector
            vector_mul(pp,dir,e2);
            det=vector_mul(e1,pp);
            // Ray is parallel to plane
            if (fabs(det)<1e-8) return -1.0;
            idet=1.0/det;
            vector_sub(rr,pos,p0);
            u=vector_mul(rr,pp)*idet;
            if ((u<0.0)||(u>1.0)) return -1.0;
            vector_mul(qq,rr,e1);
            v=vector_mul(dir,qq)*idet;
            if ((v<0.0)||(u+v>1.0)) return -1.0;
            // distance
            t=vector_mul(e2,qq)*idet;
            if (t<0.0) t=-1.0;
            return t;
            }
        };
    List<_fac> fac;                         // faces

    STL3D() { reset(); }
    STL3D(STL3D& a) { *this=a; }
    ~STL3D() {}
    STL3D* operator = (const STL3D *a) { *this=*a; return this; }
    //STL3D* operator = (const STL3D &a) { ...copy... return this; }

    void reset(){ fac.num=0; compute(); }   // clear STL
    void draw();                            // render STL mesh (OpenGL)
    void draw_normals(float size);          // render STL normals (OpenGL)
    void compute();                         // compute bbox
    void compute_normals();                 // recompute normals from points
    void supports(reper &obj);              // compute supports with obj placement above base plane z=0
    void load(AnsiString name);
    void save(AnsiString name);
    };
//---------------------------------------------------------------------------
void STL3D::draw()
    {
    _fac *f; int i,j; BYTE r,g,b;
    glBegin(GL_TRIANGLES);
    for (f=fac.dat,i=0;i<fac.num;i++,f++)
        {
        glNormal3fv(f->n);
        if (f->attr<32768)
            {
            r= f->attr     &31; r<<=3;
            g=(f->attr>> 5)&31; g<<=3;
            b=(f->attr>>10)&31; b<<=3;
            glColor3ub(r,g,b);
            }
        for (j=0;j<3;j++) glVertex3fv(f->p[j]);
        }
    glEnd();
    }
//---------------------------------------------------------------------------
void STL3D::draw_normals(float size)
    {
    _fac *f;
    int i; float a[3],b[3];
    glBegin(GL_LINES);
    for (f=fac.dat,i=0;i<fac.num;i++,f++)
        {
        vectorf_add(a,f->p[0],f->p[1]);
        vectorf_add(a,a      ,f->p[2]);
        vectorf_mul(a,a,1.0/3.0);
        vectorf_mul(b,f->n,size); glVertex3fv(a);
        vectorf_add(b,b,a);       glVertex3fv(b);
        }
    glEnd();
    }
//---------------------------------------------------------------------------
void STL3D::compute()
    {
    _fac *f;
    int i,j,k;
    double p0[3],p1[3];
    vector_ld(center,0.0,0.0,0.0);
    vector_ld(size,0.0,0.0,0.0);
    rmax=0.0;
    if (fac.num==0) return;
    // bbox
    for (k=0;k<3;k++) p0[k]=fac.dat[0].p[0][k];
    for (k=0;k<3;k++) p1[k]=fac.dat[0].p[0][k];
    for (f=fac.dat,i=0;i<fac.num;i++,f++)
     for (j=0;j<3;j++)
      for (k=0;k<3;k++)
        {
        if (p0[k]>f->p[j][k]) p0[k]=f->p[j][k];
        if (p1[k]<f->p[j][k]) p1[k]=f->p[j][k];
        }
    vector_add(center,p0,p1); vector_mul(center,center,0.5);
    vector_sub(size  ,p1,p0); vector_mul(size  ,size  ,0.5);
                      rmax=size[0];
    if (rmax<size[1]) rmax=size[1];
    if (rmax<size[2]) rmax=size[2];
    // attr repair
    for (f=fac.dat,i=0;i<fac.num;i++,f++)
     if (f->attr==0) f->attr=32768;
    }
//---------------------------------------------------------------------------
void STL3D::compute_normals()
    {
    _fac *f; int i;
    for (f=fac.dat,i=0;i<fac.num;i++,f++) f->compute();
    }
//---------------------------------------------------------------------------
void STL3D::supports(reper &obj)
    {
    _fac *f,ff;
    int i,j,k;
    double p[3],dp[3],x0,y0,h0,x1,y1,x2,y2,h1,t;
    // some config values first
    const WORD   attr0=31<<10;              // support attr should be different than joint
    const WORD   attr1=31<<5;               // joint attr should be different than mesh,support
    const double grid0=8.0;                 // distance between supports
    const double grid1=2.0;                 // distance between joints
    const double gap=grid0/tan(45.0*deg);// distance between main support and mesh (joint size)
    const double ha=1.0;                    // main support side size
    // do not mess with these
    const double hx=    ha*cos(60.0*deg);   // half size of main support in x
    const double hy=0.5*ha*sin(60.0*deg);   // half size of main support in y
    const double grid2=0.4*hy;              // distance between joints bases
    const double ga=2.0*grid2*grid1/grid0;  // main support side size
    const double gx=hx*grid2/grid0;         // half size of joint support in x
    const double gy=hy*grid2/grid0;         // half size of joint support in y

    // apply placement obj (may lose some accuracy) not needed if matrices are not used
    for (f=fac.dat,i=0;i<fac.num;i++,f++)
        {
        for (j=0;j<3;j++)
            {
            for (k=0;k<3;k++) p[k]=f->p[j][k];  // float->double
            obj.l2g(p,p);
            for (k=0;k<3;k++) f->p[j][k]=p[k];  // double->float
            }
        for (k=0;k<3;k++) p[k]=f->n[k]; // float->double
        obj.l2g_dir(p,p);
        for (k=0;k<3;k++) f->n[k]=p[k]; // double->float
        } compute();

    // create supports
    for (x0=center[0]-size[0]+(0.5*grid0);x0<=center[0]+size[0]-(0.5*grid0);x0+=grid0)
     for (y0=center[1]-size[1]+(0.5*grid0);y0<=center[1]+size[1]-(0.5*grid0);y0+=grid0)
        {
        // cast ray x0,y0,0 in Z+ direction to check for mesh intersection to compute the support height h0
        h0=center[2]+size[2]+1e6;
        vector_ld(p,x0,y0,0.0);
        vector_ld(dp,0.0,0.0,+1.0);
        for (f=fac.dat,i=0;i<fac.num;i++,f++)
            {
            t=f->intersect_ray(p,dp);
            if ((t>=0.0)&&(t<h0)) h0=t;
            }
        if (h0>center[2]+size[2]+1e5) continue; // skip non intersected rays
        h0-=gap; if (h0<0.0) h0=0.0;
        // main suport prism
        ff.attr=attr0;
        // sides
        ff.attr=attr0;
        vectorf_ld(ff.p[0],x0-hx,y0-hy,0.0);
        vectorf_ld(ff.p[1],x0+hx,y0-hy,0.0);
        vectorf_ld(ff.p[2],x0-hx,y0-hy, h0); ff.compute(); fac.add(ff);
        vectorf_ld(ff.p[0],x0+hx,y0-hy,0.0);
        vectorf_ld(ff.p[1],x0+hx,y0-hy, h0);
        vectorf_ld(ff