sandbox/Antoonvh/tree/tc.c

    Example of the space colonization algorithm

    Growing branches kill attraction points

    howto

    First the relevant header files are included

    #include "grid/octree.h"
    #include "view.h"
    #define BVIEW 1
    int n_part;
    #include "../scatter.h"
    #include "colonization.h"

    A function is defined that plots the growing root and the remaining attraction points.

    trace
    void movie_maker (Tnode * tnodes, int tn, coord * atr, int na) {
      Branch * trees = calloc (499, sizeof(Branch)); //?
      nb = nodestotree (tnodes, tn, &trees);
      scalar cs[], R[];
      cs.prolongation = cs.refine = fraction_refine;
      face vector fs[];
      tree_interface (trees, cs, fs, R); //This maybe expensive
      boundary({cs});
      
      view (fov = 31, psi = -pi/4, bg = {0.6, 0.45, 0.33});
    #if dimension == 2
      draw_vof ("cs", "fs", filled = -1, fc = {0.8,0.8,0.8});
    #else
      view (quat = {0.433278,0.701835,0.0938225,0.557577}, ty = -.3);
      draw_vof ("cs", "fs", fc = {0.8,0.8,0.8});
    #endif
      n_part = na;
      scatter (atr);
      save ("mov.mp4");
      adapt_wavelet ({cs}, (double[]){0.01}, 9, 5);
      free (trees);
      dump();
    }

    We setup a scenario with a single seed node at the edge of the domain that is attracted to 500 attraction points, randomly placed in a piramid.

    int main() {
      init_grid (128);
      L0 = 2.1;
      X0 = Y0 = Z0 = -L0/2;
      int tn = 1;
      Tnode * nodes = calloc (tn, sizeof(Tnode));
      foreach_dimension()
        nodes[0].x = 1;
      int na = 250;
      coord atr[na];
      for (int k = 0; k < na; k++) {
        coord ap;
        do {
          ap = (coord){noise(), noise(), noise()};
        } while (ap.x + ap.y + ap.z < 0);
        atr[k] = ap;
      }
          
    #if _MPI
      MPI_Bcast (&atr[0], na*sizeof(coord), MPI_BYTE, 0, MPI_COMM_WORLD); 
    #endif

    The nodes are computed and the movie_maker function is called each iteration.

      tn =  colonize (atr, na, &nodes, tn, 0.5, .025, 0.1,
    	    myfun = movie_maker);
      free (nodes);
    }