🚀 go-pugleaf

Howdy:

When compiling the file included below with g++-2.2.2 on a Sun
sparcstation running SunOS 4.1.1, we get the following output:

======================================================================
% gcc -v
Reading specs from /usr/local/gnu/gcc-2.2.2/lib/gcc-lib/sparc-sun-sunos4.1/2.2.2
/specs
gcc version 2.2.2
% g++ high.c
high.c: In method `void  topbox::split_all_crowded (auto int (class wholebox::*)
())':
high.c:688: internal compiler error 116, please report
======================================================================

As it turns out, there's a bug in the sample program (look for the
line marked with "//###"); when fixed, g++ doesn't complain.  We
thought we should report the problem anyway.

						--Cliff Krumvieda

high.c:
======================================================================
// high.c

//refine the octree this is the highest level call

#include <math.h>  //just for sqrt
#include <stream.h>

enum neighbor_direction {no_dir = -1, west = 0, north, up, east, south, down};
// The first three directions are perpendicular to one another,
// and their value corresponds to the bit position in corner_direction
// determining that direction + or -

enum corner_direction {  no_corner = -1,
	downsoutheast = 0, downsouthwest, downnortheast, downnorthwest, 
	upsoutheast, upsouthwest, upnortheast, upnorthwest  
};
#define DSE downsoutheast
#define DSW downsouthwest
#define DNE downnortheast
#define DNW downnorthwest
#define USE upsoutheast
#define USW upsouthwest
#define UNE upnortheast
#define UNW upnorthwest

static inline neighbor_direction opp_neighbor(int n) {
    return (neighbor_direction) (n > 2) ? n - 3 : n + 3;
}

static inline corner_direction opp_corner(corner_direction d) {
    return (corner_direction) ( ((int) d) ^ 7 );
}

static inline int is_sibling(corner_direction wc, neighbor_direction n) {
// n assumed in [0..2] = [west, north, up]
// returns 1 if samed sized neighbor in direction n must exist because
// of this->which_child() = wc
//
        return (wc >> n) ^ 1;  //return opposite of bit in position n
}

static inline int sibling_index(corner_direction wc, neighbor_direction n) {
// wc is this_octbox->which_child()
// asume n in [0..2]
// returns sibling in direction n, or direction n+3, wherever sibling exists
	return wc ^ (1 << n);
}

static inline int sibling_dir(corner_direction wc, neighbor_direction n) {
// wc is this_octbox->which_child()
// asume n in [0..2]
	return (is_sibling(wc, n)) ? n : n+3;
}

static inline int no_sibling_dir(corner_direction wc, neighbor_direction n) {
// wc is this_octbox->which_child()
// asume n in [0..2]
	return (is_sibling(wc, n)) ? n+3 : n;
}

#define INC(a) (a = a->rest())
#define num_children 8
#define num_neighbors 6
#define num_parent_neighbor 3

//forward declarations
class wholebox; //see wholebox.h
class octbox; class headbox : public octbox; class dupbox : public octbox; 

typedef octbox * childarray[num_children];
typedef headbox * headchildarray[num_children];

class octbox {
  private:
    int leaf;

  protected:

    dupbox * next;

    //CONSTRUCTORS
    // no one should construct a pure octbox, so constructors are only to be
    // called by the derived classes
    octbox() {leaf = 1; w = 0; next = 0;}
    dupbox * duplicate(headbox * h);  //verb
    void oct_split() {leaf = 0;}

  public:
    virtual int size() {return 0;} //size = level of octree

    // is that smaller that this ?
    int is_smaller( octbox * that ) {return that->size() < size();}

    wholebox * w; //wholebox and octbox point to each other

    //for duplicate list
    dupbox * rest() {return next;}
    octbox * last_box() { 
    for (octbox * n = this; n->rest(); INC(n)); 
    return n;
    }

    virtual headbox * child(int n) {return 0;}
    virtual headbox * neighbor(int n) {return 0;}
    virtual headbox * parent() {return 0;}
    virtual dupbox * duplicate() {return 0;}
 
    int is_top() {return (!parent());}
    int is_leaf() {return (leaf);}
    virtual int which_child() {return -1;}
    virtual int is_children() {return 0;} 

    // return parent()->child(c)
    headbox * sibling(int c);

    virtual void split(childarray * c = 0) {;}
    // if c != 0 then we fill up what it points to with the 
    // new children we make
    // if c == 0, then we asume the caller is uninterested in looking
    //  at the children, and we leave c alone, but still create the children
    // Same for split() in wholebox.h .c
};


class headbox : public octbox {
  private:
    friend dupbox;
    headbox * neighbors[num_neighbors];
    headbox * father;
    headchildarray children;
    int child_index;
    void makechildren( childarray * c );
    int level; //level of tree this resides in

  public:
    virtual int size() { return level; }
    virtual headbox * child(int n) {return children[n];}
    virtual headbox * neighbor(int n) {return neighbors[n];}
    virtual headbox * parent() {return father;}
    virtual dupbox * duplicate() {return octbox::duplicate(this);}

    virtual void split(childarray * c = 0) { 
        octbox::oct_split(); makechildren( c );
    }

    virtual int is_children() {return (!child(1));}
    virtual int which_child() {return child_index;}

    //CONSTRUCTOR
    headbox(headbox * father_in = 0, int which_child_in = -1); //see octbox.c
};

class dupbox : public octbox {
  private:
    friend headbox;
    headbox * head;
  public:
    virtual int size() { return head->level; }
    virtual headbox * child(int n) {return head->child(n);}
    virtual headbox * neighbor(int n) {return head->neighbor(n);}
    virtual headbox * parent() {return head->parent();}
    virtual dupbox * duplicate() {return octbox::duplicate(head);}

    virtual void split(childarray * c = 0) { 
        octbox::oct_split(); head->makechildren(c);
    }

    virtual int which_child() {return head->which_child();}
    virtual int is_children() {return head->is_children();}

    //CONSTRUCTOR
    dupbox(headbox * headin) {head = headin;}
};


//===================================================
//=====================  search =====================
//===================================================
//For traveling the octree in various ways

//should only be accessed by search and its derivatives
class stackitem {
 public:
  stackitem * next;
  octbox * item;
  stackitem(octbox * putbox) {item = putbox; next = 0;}
};

//can only be used as part of another class
class stack {
 protected:  
  //pop always off top
  octbox * pop();
  stackitem * top;
  stackitem * bottom;
  stack () {top = bottom = 0;}
  ~stack () { for ( ; pop(); ); }  //pop everything
  int emptystack() {return (top == 0);} 
};

class search : public stack {
 public:
  int not_done() {return (current != 0);}
 protected:
  octbox * current;
  search( octbox * top ) {current = top;}
  virtual void do_stack(octbox * putbox) {;}
 public:
  virtual void inc();  //main work done here, dependant on do_stack
  octbox * o() {return current;}
  wholebox * w() {return current->w;}
};

class depthfs : public search {
 protected:
  virtual void do_stack(octbox * putbox);
 public:
  depthfs( octbox * top ) : search(top) {;}
  virtual void inc() { search::inc(); }
};

class depthfs_leaf : public depthfs {
 public:
  depthfs_leaf( octbox * top) : depthfs(top) {;}
  virtual void inc();
};

class breadthfs : public search {
 protected:
  virtual void do_stack(octbox * putbox);
 public:
  breadthfs( octbox * top ) : search(top) {;}
  virtual void inc() {search::inc();}
};

class breadthfs_leaf : public breadthfs {
 public:
  breadthfs_leaf( octbox * top) : breadthfs(top) {;}
  virtual void inc();
};

typedef long idtype;

class out {
private:
    idtype id;
    static idtype maxid;   

public:
    // CONSTRUCTORS
    out(idtype idin = 0) { id = idin; }

    friend void start_label() {  /* needs no this pointer passed */
	out::maxid = 0;
    }
    void label_it() { if (id==0) id = ++maxid; }
    idtype label() { return id; } 
    void dup(out* that) {
	id = that->label();
    }
};

idtype out::maxid;  //can't set to zero here, as included several times

//coordinates for points
typedef double coor;
#define MAX(a,b) (a>b ? a:b)
#define MIN(a,b) (a<b ? a:b)
#define MAX3(a,b,c) (a>b ? (a>c ? a:c ) : (b>c ? b:c))
#define MIN3(a,b,c) (a<b ? (a<c ? a:c ) : (b<c ? b:c))
#define SQR(a) ((a)*(a))

class point {
public:
    coor x,y,z;
    point(coor ix = -99, coor iy = -99, coor iz = -99) {x=ix; y=iy; z=iz;}
    point(point & inpt) {x=inpt.x; y=inpt.y; z=inpt.z;}

 
    //-- VECTOR ARITHMATIC --
  
    // comparison
    int operator==(point& that) {
        return (this->x == that.x && this->y == that.y && this->z == that.z);
    }
    int operator!=(point& that) {
        return (this->x != that.x || this->y != that.y || this->z != that.z);
    }

    // vector sum
    point operator+(point& that) {
        return point(this->x + that.x, this->y + that.y, this->z + that.z);
    }
    void operator+=(point& that) {
        this->x += that.x; this->y += that.y; this->z += that.z;
    }
    point operator-() {
        return point(- this->x, - this->y, -this->z);
    }
    point operator-(point& that) {
        return point(this->x - that.x, this->y - that.y, this->z - that.z);
    }
    void operator-=(point& that) {
        this->x -= that.x; this->y -= that.y; this->z -= that.z;
    }

    // dot product
    coor operator*(point& that) {
        return (coor) this->x * (coor) that.x +
               (coor) this->y * (coor) that.y +
               (coor) this->z * (coor) that.z;
    }  

    // scalar product
    point operator*(coor a) {
        return point( a * x, a * y, a * z );
    }
    point operator*(int a) {
        return point( a * x, a * y, a * z );
    }
    point operator/(coor a) {
        return point( x / a, y /a, z /a );
    }

    coor normsqr() {return SQR(*this);}  //dot product
    coor norm() {return sqrt( normsqr() );}
    point unit() {return  *this / norm();}  //unit vector with same direction

};

inline point unitvector(coor scale = 1) {return point(scale,scale,scale);};

inline point midpoint(point & a, point & b) {
    return point((a.x + b.x) / 2, (a.y + b.y)/2); };

//linked list of faces, formerly points
class geo; //forward
class lattice; //forward
class geolist {  
  private:
    geolist * next;  

  public:

    #define INC(a) (a = a->rest())

    // static public members - like global vars - for readdata 
    // shared by all geolist class objects
    static point minpoint, MAXpoint;

    // the geometric face
    geo * g;
    geo * first() { return g; }

    //CONSTRUCTOR
    geolist(geo * ip = 0, geolist * in = 0) { next=in; g=ip; }
    // create a new geolist, put on list after this but before this->next
    geolist* addafter( geo * ip ) { 
        this->next = new geolist(ip,this->next);
        return this->next;
    }
    //add geolist with g1, g2 to this if they aren't already there
    //returns head of list
    geolist * addunique( geo * g1, geo * g2 = 0 );

    //DESTRUCTOR
    //base class list destructor automatically called
    ~geolist() { if (next != 0) delete next; }

    geolist * rest() { return next; }
    geolist * pop() {
        geolist * retval = next;
        next = 0;
        delete this;
        return retval;
    }


    // for setting and searching for g->mark in this and this->rest() etc
    void markgeos();
    int find_unmarked(); // returns 1 if there is at least one unmarked geo on this list
    void unmarkgeos();

    int contains_face(geo * F);

    // return lattice of geolist p and edges e
    friend void readdata(char * fname, lattice & P);
};


class geo {

  friend geolist;
  protected:    
    int dimension;
    int mark;

  private: 

  public:
    geolist *superfaces, *subfaces;

    //CONSTRUCTOR
    geo(int dim=0) {dimension=dim; subfaces = 0; superfaces = 0; mark = 0;}
    int what_dim() {return dimension;}

    // defined differently for each face dimension
    // at present, returns zero if same dimensional function as this is called
    // I will change this if the need arises, but it will take more storage
    virtual geolist * facets() {return 0;}
    virtual geolist * edges() {return 0;}
    virtual geolist * vertices() {return 0;}

    virtual void describe(ostream os);  //for drawing

    // set this to contain argument, and vice versa
    void contains( geo * sub ) {
        subfaces = new geolist(sub,subfaces);  //add sub to subfaces
        // add this to sub's superfaces
        sub->superfaces = new geolist(this,sub->superfaces);
     }
    
    // set this to be contained in argument
    void contained( geo * sup ) {
        superfaces = new geolist(sup,superfaces); //add sup to superfaces
        // add this to sup's subfaces
        sup->subfaces = new geolist(this,sup->subfaces);
    }
};


class vertex: public geo, public point, public out {
  private:
  public:
    vertex() :geo(0) {;}
    vertex(coor x, coor y, coor z) : geo(0) , point(x,y,z) {;}
    vertex(point &p) : geo(0) , point(p) {;}

    geolist * facetlist;
    virtual geolist * facets() {return facetlist;}
    virtual geolist * edges() {return superfaces;}
    virtual geolist * vertices() {return 0;}

    virtual void describe(ostream os);
};

class edge:public geo {
  public:
    edge(geo * vertex1 = 0, geo * vertex2 = 0) : geo(1) {
        if (vertex1 != 0) { contains(vertex1); }
        if (vertex2 != 0) { contains(vertex2); }
    }

    vertex * endpt1() { return (vertex*) subfaces->first(); }
    vertex * endpt2() { return (vertex*) subfaces->rest()->first(); }

    virtual geolist * facets() {return superfaces;}
    virtual geolist * edges() {return 0;}
    virtual geolist * vertices() {return subfaces;}
    virtual void describe(ostream os);
};

class facet:public geo {
  private:
    point perpsav;
  protected:
    geolist * verts;
    point normal; //perpendicular vector
    friend void readdata(char * fname, lattice & P);
  public:
    facet() : geo(2) {verts=0; normal=unitvector();} 
    // also call geo constuctor

    point perp() {return normal;} 
    //returns perpendicular vector, pointing to interior

    point base() {return * (vertex *) (verts->g);} 

    //returns a point contained in the plane

    virtual geolist * facets() {return 0;}
    virtual geolist * edges() {return subfaces;}
    virtual geolist * vertices() {return verts;}

    virtual void describe(ostream os);
};

//this is only used for the output lattice, not the input lattice
class tetra:public geo {
  protected:
    geolist * verts;
    friend void readdata(char * fname, lattice & P);
  public:
    tetra() : geo(3) {verts = 0;}
    virtual void describe(ostream os);

    //edges(), vertices(), facets() not implemented
};

class lattice:public geo {
  private: 
     geolist * edgs, * verts, * facts;
  public:
    lattice() : geo(4) {verts=0; edgs=0; facts=0;}
    geolist * vertices() {return verts;}
    geolist * edges()    {return edgs;}
    geolist * facets() {return facts;}
    geolist * tetra() {return subfaces;}

    friend void readdata(char * fname, lattice& P);
    virtual void describe(ostream os);

}; //lattice

class pbox {
  private:
    geolist * Pvertlist;
    geolist * Pedgelist;
    geolist * Pfacetlist;

  public:

    geolist * Pvertices() {return Pvertlist;}
    geolist * Pedges() {return Pedgelist;}
    geolist * Pfacets() {return Pfacetlist;}


    //CONSTRUCTOR
    pbox() { Pvertlist = 0; Pedgelist = 0; Pfacetlist = 0; }

    void distribute_contains( pbox * c ); 
    //this is "parent", c is "child"

//  int adjacent( pbox * b ); 
    // this is neighboring box. direction available if needed in implement.
    // b is box to test adjacency with

    int has_Pvertex() {return (Pvertlist != 0);}
    int has_Pedge() {return (Pedgelist != 0);}
    int has_Pfacet() {return (Pfacetlist != 0);}

    // Does this contain a  dim dimensional face?
    int has_Pface(int dim = 5);

    int has_particular_Pface(geo * F);

    int contains_non_superface( geo * F );
};

// pointer to member function taking no args and returning int:
// has_Pvertex(), for example
//See prog in C++ book by Dewhurst page 72 bottom
//binfn returns 1 if want premature exit from to_all_balance_adj
class wholebox; //forward
typedef int (wholebox::*boolfn ) ();  
typedef void (wholebox::*voidfn)();
typedef int (wholebox::*binfn)( wholebox * ); 

class wholebox : public pbox {
   private:
    int frozen;
    pbox * extended;
 
   public:
    octbox * o;
    int is_frozen() {return frozen;}

    //CONSTRUCTORS
    //there is *ALWAYS* an octbox associated with a wholebox
    wholebox(octbox * ob) {
        frozen = 0; 
    extended = 0;
        o = ob;
        ob->w = this;
    }

    wholebox() {
    // for bounding box construction
    frozen = 0; 
    extended = 0;
        o = new headbox();
        o->w = this;
    }


    //split balance adjacent boxes to this size
    void balance();
    private:
      int balance_split ( wholebox * box_to_split );
    public:

    void split(childarray * c = 0); 
    // "return" children created as c, 
    // but each c is merely the first in the list
    // of children created, as children are duplicated here as well

    int merge( wholebox * that ) {return 0;} //binfn
    // Called by topbox::merge_neighbors
    // That is a box balance adjacent to this
    // Use a static pointer to keep track of when merge is being
    //  called for a new box.
    // Move out the faces of this box.
    // There is no need to change contains information - all balance_adj
    //  boxes just contain superfaces of a vertex in this
    // Set some flag in the merged boxes, so that when we share faces,
    //  we will know to look in this box instead.

    void duplicate_for_components_and_extended() {;}
    //called by wholebox::split
    //duplicates children for each connected component, and also
    //distributes faces of P, and makes extended box
        //find intersection of Pfaces in parent->contains and boundary of this box
        //identify components
        //find intersection of Pfaces in parent->extended->contains and . . .
        //boundary of this->extended
        //
        //for each component of P cap bdy this
        //  make a new box (if second or more component)
        //  put faces in the list of P faces ("distribute contains")
        //  make a new extended box
        //  identify corresponding component of extended box and put P faces
        //  into extended box
        //

    void balance_same(geo * F); //make every box balance_adj same size, 
                                // may split this!
    int is_smaller_balance_adj();
    int is_smaller( wholebox * that ) { return o->is_smaller( that->o ); }

    void freeze_all( boolfn should_freeze ); //recursively protect boxes
    //example call,  protect( & pbox::has_Pvertex )
    // still to do, change the extended box information 
    //  or we could leave it as is and have worse constant size bound
    private:
      int freeze( wholebox * dummy = 0 ) { frozen = 1; return 0; }
    public:

    int adjacent( pbox * b ) {return 0;}  //heavy duty goemetry needed

    int to_all_balance_adj( binfn fn, int bigger_only = 0);

    int is_vertex_crowded();
    int is_edge_crowded();
    int is_facet_crowded();

};

class topbox : public wholebox {
 private:
  void split_all_crowded( boolfn is_crowded );
  void balance_same_all( boolfn has_inphase_Pface );

  //cleanup for merging vertex boxes with same sized neighbors
  void merge_neighbors( boolfn has_inphase_Pface );
  void share_all_faces();

 public:
  void refine();  //refine the octree, this is the only public call 
};

void topbox::refine() {

  // ============  Vertex Phase  ============
  split_all_crowded( &is_vertex_crowded );
  balance_same_all( &wholebox::has_Pvertex );
  freeze_all( &wholebox::has_Pvertex );

  // ============  Edge Phase  ============
  split_all_crowded(  &wholebox::is_edge_crowded );
  freeze_all( &wholebox::has_Pedge );
  
  // ============  Facet Phase  ============
  split_all_crowded(  &wholebox::is_facet_crowded );
  freeze_all( &wholebox::has_Pfacet );  //needed??

  //cleanup, merge vertex frozen boxes with same sized neighbors
  merge_neighbors( &wholebox::has_Pvertex );

  //make sure all balance adjacent boxes share their common faces,
  // so that there are not copies of any face.
  share_all_faces();
}

//=======================================================================
void topbox::split_all_crowded( boolfn is_crowded ) {
     for (breadthfs_leaf bfs(o); bfs.not_done(); bfs.inc() ) {

//### If the line below is changed to 
//###   if ((bfs.w()->*is_crowded)()) {
//### g++ doesn't get an internal compiler error.

        if (bfs.w()->*is_crowded()) {  
            //creates children, put on stack eventually by inc
            bfs.w()->balance();
            bfs.w()->split();    //w calls o.split
	      //children will be duplicated for connected components
        }
    }
}

//following only works for vertices at present because of balance_same call
// using vertex of this box
void topbox::balance_same_all( boolfn has_inphase_Pface ) {
    for ( depthfs_leaf bfs(o) ; bfs.not_done(); bfs.inc() ) {
        if ((bfs.w()->*has_inphase_Pface)()) {
            bfs.w()->balance_same( Pvertices()->g ); 
            //Note any box created by this will not be one we would need to
            // consider later in the search, so we don't care if the boxes
            // split by this are stacked correctly or not!
        }
    }
}

//cleanup for merging vertex boxes with same sized neighbors
void topbox::merge_neighbors( boolfn has_inphase_Pface ) {
    //search order irrelevant
    for ( depthfs_leaf bfs(o); bfs.not_done(); bfs.inc() ) {
        if ((bfs.w()->*has_inphase_Pface)()) {
            bfs.w()->to_all_balance_adj( & wholebox::merge );
        }
    }
}

void topbox::share_all_faces() {
    //search order irrelevant
//    for ( depthfs_leaf bfs(o); bfs.not_done(); bfs.inc() ) {
//          bfs.w()->to_all_balance_adj( & wholebox::share_faces );
//    }
;
}
======================================================================