Tutorials (back to the list of tutorials)

## Cell Division and Growth Algorithm 2 (requires iGeo version 0.9.1.5)

### Triangular Link for Surface Formation and Peripheral Edge Division

The tutorial codes in this page extend the cell division and growth algorithm in the previous page by adding data structure of a triangular face which defines relationship of three cell agents in addition to a cell link which defines relationship of two, to explore surface formation through cell division and growth.

The first code in this page adds a new class CellFace which contains references to three cells and three links. Cell class also stores what CellFace instances contain the cell in the variable array field faces. CellLink also stores what CellFace instances contain the cell link in faces field.

The code below adds an algorithm to divide a face in Cell class's divide method. Conceptually the division of a face is done by dividing one of edge of the face triangle and insert a cell there and divide the original triangle into two with a new link in the middle as depicted in the diagram below. In the code, it's done by deleting the edge link to divide and the original face, and then add three new links and two new faces. When a cell contains only one link (only one line connecting two cells), it creates two links and one face with the new child cell to have a first face to be divided later.

```import igeo.*;
import processing.opengl.*;

void setup(){
size(480,360,IG.GL);
IConfig.syncDrawAndDynamics=true; //not to crash when some geometry is deleted while drawing
Cell cell = new Cell(new IVec(0,0,0), 10);
}

class Cell extends IParticle{
int growthDuration = 1400; //duration of growth and division
int growthInterval = 10;
int divisionInterval = 50;
double growthSpeed=0.1;

ArrayList< CellFace > faces; //store faces
boolean active = false;

super(pos, new IVec(0,0,0));
faces = new ArrayList< CellFace >();
fric(0.2);
}

void interact(ArrayList< IDynamics > agents){
IVec neighborCenter = new IVec(0,0,0);
int neighborCount=0;
for(int i=0; i < agents.size(); i++){
if(agents.get(i) instanceof Cell){
Cell c = (Cell)agents.get(i);
if(c != this){
// push if closer than two radii
IVec dif = c.pos().dif(pos());
c.push(dif);
}
// count neighbors and calculate their center
neighborCount++;
}
}
}
}
if(neighborCount >= 1){ // push from center of neighbors
neighborCenter.div(neighborCount);
IVec dif = pos().dif(neighborCenter).len(50); // constant force
this.push(dif);
}
}

void update(){
if(IG.time() < growthDuration){
if(time() > 0 && time()%divisionInterval==0){
if(active){ // divide when active flag is on
divide();
}
if(IRand.pct(50)){ // random division
active=true;
}
}
if(time()%growthInterval==0){
grow();
}
if(active){ clr(1.0,0.5,0); } // update color by being active
else{ clr(1.0,0,0); }
}
}

void grow(){ // growing cell size
}
}

void divide(){ // cell division
Cell child = createChild(IRand.dir());
}
Cell child = createChild(IRand.dir());
}
else{ // strip state
active = false;
return; // if link has two faces, skip division
}
IVec dir = dividingLink.oppositeDir(this); // vector to the other cell
Cell child = createChild(dir); // make a child toward the link direction

Cell c1 = f1.oppositeCell(dividingLink); // opposite cell from the link on the face
// delete 1 link + 1 face, create 3 links + 2 faces
new CellFace(this, c1, child);
new CellFace(child, c1, c0);
f1.del();
}
}

Cell createChild(IVec dir){
radius *= 0.5; //make both cell size half
Cell child = new Cell(pos().cp(dir), radius);
pos().sub(dir); // move to the other way
active=false;  //reset activation
return child;
}
}

class CellLink extends IAgent{ // a link to connect 2 cells with spring force
double maxForce = 100;
Cell cell1, cell2;
ArrayList< CellFace > faces;
ICurve line;

cell1 = c1; cell2 = c2;
faces = new ArrayList< CellFace >();
line = new ICurve(c1.pos(), c2.pos()).clr(1.,0,0);
}

void interact(ArrayList< IDynamics > agents){
// spring force
IVec dif = cell1.pos().dif(cell2.pos());
if(force > maxForce){ force=maxForce; } // avoid putting too much force
else if(force < -maxForce){ force=-maxForce; }
dif.len(force);
cell1.pull(dif);
cell2.push(dif);
}

boolean contains(Cell c){ //check if link contains the cell
if(c==cell1 || c==cell2) return true;
return false;
}

void del(){
line.del(); // delete line geometry
super.del(); // stop agent
}

Cell oppositeCell(Cell c){ // find other cell on the link
if(cell1==c) return cell2;
if(cell2==c) return cell1;
IG.err("Link does not contain the input cell");
return null;
}

IVec oppositeDir(Cell c){ //calculate a vector to the other cell
return oppositeCell(c).pos().dif(c.pos());
}
}

class CellFace{ // a triangle grouping 3 cells and 3 links
Cell cell1, cell2, cell3;

CellFace(Cell c1, Cell c2, Cell c3){
cell1 = c1;
cell2 = c2;
cell3 = c3;
}

IVec center(){ // calc center of triangle
return IVec.center(cell1.pos(), cell2.pos(), cell3.pos());
}

}

for(int i=0; i < c1.links.size() ; i++){
}
return null;
}

Cell oppositeCell(CellLink l){ // find a cell opposite of a link edge in a triangle
if(cell1!=l.cell1 && cell1!=l.cell2){ return cell1; }
if(cell2!=l.cell1 && cell2!=l.cell2){ return cell2; }
return cell3;
}

void del(){
cell1.faces.remove(this); // unregister this from cells and links
cell2.faces.remove(this);
cell3.faces.remove(this);
}
}
```

### Visualization of Cell with Spherical Body

The code below has a sphere as the geometric representation of a cell instead of a circle. The code uses IG.meshSphere producing IMesh instance instead of using ISphere (NURBS surface sphere) to keep the code's memory usage lighter.

```import igeo.*;
import processing.opengl.*;

void setup(){
size(480,360,IG.GL);
IConfig.syncDrawAndDynamics=true; //not to crash when some geometry is deleted while drawing
Cell cell = new Cell(new IVec(0,0,0), 10);
cell.clr(0,0,1.0);
}

class Cell extends IParticle{
int growthDuration = 1400; //duration of growth and division
int growthInterval = 10;
int divisionInterval = 50;
double growthSpeed=0.1;

ArrayList< CellFace > faces; //store faces
boolean active = false;
IMesh sphere;

super(pos, new IVec(0,0,0));
faces = new ArrayList< CellFace >();
fric(0.2);
}

void interact(ArrayList< IDynamics > agents){
IVec neighborCenter = new IVec(0,0,0);
int neighborCount=0;
for(int i=0; i < agents.size(); i++){
if(agents.get(i) instanceof Cell){
Cell c = (Cell)agents.get(i);
if(c != this){
// push if closer than two radii
IVec dif = c.pos().dif(pos());
c.push(dif);
}
// count neighbors and calculate their center
neighborCount++;
}
}
}
}
if(neighborCount >= 1){ // push from center of neighbors
neighborCenter.div(neighborCount);
IVec dif = pos().dif(neighborCenter).len(50); // constant force
this.push(dif);
}
}

void update(){
if(IG.time() < growthDuration){
if(time() > 0 && time()%divisionInterval==0){
if(active){ // divide when active flag is on
divide();
}
if(IRand.pct(50)){ // random division
active=true;
}
}
if(time()%growthInterval==0){
grow();
}
}
// update geometry
if(sphere!=null) sphere.del();
}

void grow(){ // growing cell size
}
}

void divide(){ // cell division
Cell child = createChild(IRand.dir());
}
Cell child = createChild(IRand.dir());
}
else{ // strip state
active = false;
return;  // if link has two faces, skip division
}
IVec dir = dividingLink.oppositeDir(this); // vector to the other cell
Cell child = createChild(dir); // make a child toward the link direction

Cell c1 = f1.oppositeCell(dividingLink); // opposite cell from the link on the face
// delete 1 link + 1 face, create 3 links + 2 faces
new CellFace(this, c1, child);
new CellFace(child, c1, c0);
f1.del();
}
}

Cell createChild(IVec dir){
radius *= 0.5; //make both cell size half
Cell child = new Cell(pos().cp(dir), radius);
child.hsb(hue()-0.05,1,0.8);
pos().sub(dir); // move to the other way
active=false;  //reset activation
return child;
}
}

class CellLink extends IAgent{ // a link to connect 2 cells with spring force
double maxForce = 100;
Cell cell1, cell2;
ArrayList< CellFace > faces;
ICurve line;

cell1 = c1; cell2 = c2;
faces = new ArrayList< CellFace >();
line = new ICurve(c1.pos(), c2.pos()).clr(1.,0,0);
}

void interact(ArrayList< IDynamics > agents){
// spring force
IVec dif = cell1.pos().dif(cell2.pos());
if(force > maxForce){ force=maxForce; } // avoid putting too much force
else if(force < -maxForce){ force=-maxForce; }
dif.len(force);
cell1.pull(dif);
cell2.push(dif);
}

boolean contains(Cell c){ //check if link contains the cell
if(c==cell1 || c==cell2) return true;
return false;
}

void del(){
line.del(); // delete line geometry
super.del(); // stop agent
}

Cell oppositeCell(Cell c){ // find other cell on the link
if(cell1==c) return cell2;
if(cell2==c) return cell1;
IG.err("Link does not contain the input cell");
return null;
}

IVec oppositeDir(Cell c){ //calculate a vector to the other cell
return oppositeCell(c).pos().dif(c.pos());
}
}

class CellFace{ // a triangle grouping 3 cells and 3 links
Cell cell1, cell2, cell3;

CellFace(Cell c1, Cell c2, Cell c3){
cell1 = c1;
cell2 = c2;
cell3 = c3;
}

IVec center(){ // calc center of triangle
return IVec.center(cell1.pos(), cell2.pos(), cell3.pos());
}

}

for(int i=0; i < c1.links.size() ; i++){
}
return null;
}

Cell oppositeCell(CellLink l){ // find a cell opposite of a link edge in a triangle
if(cell1!=l.cell1 && cell1!=l.cell2){ return cell1; }
if(cell2!=l.cell1 && cell2!=l.cell2){ return cell2; }
return cell3;
}

void del(){
cell1.faces.remove(this); // unregister this from cells and links
cell2.faces.remove(this);
cell3.faces.remove(this);
}
}
```

### Building Polygon Mesh by Cells and Links

The next code generates a polygon mesh out of cells, links and faces. A global variable mesh is defined and initialized in setup method, and then mesh faces are added at CellFace's constructor and also sometimes mesh faces are deleted at del method. Cell class has a mesh vertex instance vertex for CellFace class to instantiate a mesh face (IFace) instance. Cell class and CellFace get nml method to calculate a normal vector of a face or of a vertex. To keep the normal direction consistent in adjacent faces, CellFace class adjust the order of cells using cellOrder method.

```import igeo.*;
import processing.opengl.*;

IMesh mesh; // global variable of mesh geometry

void setup(){
size(480,360,IG.GL);
IConfig.syncDrawAndDynamics=true; //not to crash when some geometry is deleted while drawing
mesh = new IMesh().clr(0.7,0,0);
Cell cell = new Cell(new IVec(0,0,0), 10);
}

class Cell extends IParticle{
int growthDuration = 1400; //duration of growth and division
int growthInterval = 10;
int divisionInterval = 50;
double growthSpeed=0.1;

ArrayList< CellFace > faces; //store faces
boolean active = false;
IVertex vertex;

super(pos, new IVec(0,0,0));
faces = new ArrayList< CellFace >();
vertex = new IVertex();
vertex.pos = pos; // set reference of particle pos directly
fric(0.2);
}

void interact(ArrayList< IDynamics > agents){
IVec neighborCenter = new IVec(0,0,0);
int neighborCount=0;
for(int i=0; i < agents.size(); i++){
if(agents.get(i) instanceof Cell){
Cell c = (Cell)agents.get(i);
if(c != this){
// push if closer than two radii
IVec dif = c.pos().dif(pos());
c.push(dif);
}
// count neighbors and calculate their center
neighborCount++;
}
}
}
}
if(neighborCount >= 1){ // push from center of neighbors
neighborCenter.div(neighborCount);
IVec dif = pos().dif(neighborCenter).len(50); // constant force
this.push(dif);
}
}

void update(){
if(IG.time() < growthDuration){
if(time() > 0 && time()%divisionInterval==0){
if(active){ // divide when active flag is on
divide();
}
if(IRand.pct(50)){ // random division
active=true;
}
}
if(time()%growthInterval==0){
grow();
}
if(active){ clr(1.0,0.5,0); } // update color by being active
else{ clr(mesh.clr()); }
}
vertex.nml(nml()); // update mesh vertex normal
}

void grow(){ // growing cell size
}
}

void divide(){ // cell division
Cell child = createChild(IRand.dir());
}
Cell child = createChild(IRand.dir());
}
else{ // strip state
active = false;
return;  // if link has two faces, skip division
}
IVec dir = dividingLink.oppositeDir(this); // vector to the other cell
dir.projectToPlane(nml()); // division dir is projected on normal plane
Cell child = createChild(dir); // make a child toward the link direction

Cell c1 = f1.oppositeCell(dividingLink); // opposite cell from the link on the face
// delete 1 link + 1 face, create 3 links + 2 faces
new CellFace(this, c1, child);
new CellFace(child, c1, c0);
f1.del();
}
}

IVec nml(){ // calc vertex normal from face normal
if(faces.size()==0){ return new IVec(0,0,1); }
IVec n = new IVec(0,0,0);
for(int i=0; i < faces.size(); i++){ n.add(faces.get(i).nml()); }
return n.unit();
}

Cell createChild(IVec dir){
radius *= 0.5; //make both cell size half
Cell child = new Cell(pos().cp(dir), radius);
pos().sub(dir); // move to the other way
active=false;  //reset activation
return child;
}
}

class CellLink extends IAgent{ // a link to connect 2 cells with spring force
double maxForce = 100;
Cell cell1, cell2;
ArrayList< CellFace > faces;

cell1 = c1; cell2 = c2;
faces = new ArrayList< CellFace >();
}

void interact(ArrayList< IDynamics > agents){
// spring force
IVec dif = cell1.pos().dif(cell2.pos());
if(force > maxForce){ force=maxForce; } // avoid putting too much force
else if(force < -maxForce){ force=-maxForce; }
dif.len(force);
cell1.pull(dif);
cell2.push(dif);
}

boolean contains(Cell c){ //check if link contains the cell
if(c==cell1 || c==cell2) return true;
return false;
}

void del(){
super.del(); // stop agent
}

Cell oppositeCell(Cell c){ // find other cell on the link
if(cell1==c) return cell2;
if(cell2==c) return cell1;
IG.err("Link does not contain the input cell");
return null;
}

IVec oppositeDir(Cell c){ //calculate a vector to the other cell
return oppositeCell(c).pos().dif(c.pos());
}
}

class CellFace{ // a triangle grouping 3 cells and 3 links
Cell cell1, cell2, cell3;
IFace face;

CellFace(Cell c1, Cell c2, Cell c3){
// keep order of cells in right hand screw for consistent normal
cell1 = c2; cell2 = c1; cell3 = c3;
}
else{
cell1 = c1; cell2 = c2; cell3 = c3;
}
face = new IFace(cell1.vertex,cell2.vertex,cell3.vertex);
}

IVec center(){ // calc center of triangle
return IVec.center(cell1.pos(), cell2.pos(), cell3.pos());
}

IVec nml(){ // normal vector
return cell1.pos().nml(cell2.pos(), cell3.pos()).unit();
}

// true if ordr of c1-c2 follows the order of cell1-cell2-cell3, otherwise false
boolean cellOrder(Cell c1, Cell c2){
if(cell1==c1 && cell2==c2 || cell2==c1 && cell3==c2 || cell3==c1 && cell1==c2 ){ return true; }
return false;
}

}

for(int i=0; i < c1.links.size() ; i++){
}
return null;
}

Cell oppositeCell(CellLink l){ // find a cell opposite of a link edge in a triangle
if(cell1!=l.cell1 && cell1!=l.cell2){ return cell1; }
if(cell2!=l.cell1 && cell2!=l.cell2){ return cell2; }
return cell3;
}

void del(){
cell1.faces.remove(this); // unregister this from cells and links
cell2.faces.remove(this);
cell3.faces.remove(this);
mesh.deleteFace(face);
}
}
```

### Interior Edge Division on Open Surface

The previous face division algorithm always divides a peripheral edge which has only one adjacent face and never divides an interior edge which has two adjacent faces. The next code divides an interior edge in addition to peripheral edge. This is done in divide method and the use of the peripheral edge division or the interior edge division is determined by the number of faces the edge has (dividingLink.faces.size()).

```import igeo.*;
import processing.opengl.*;

IMesh mesh; // global variable of mesh geometry

void setup(){
size(480,360,IG.GL);
IConfig.syncDrawAndDynamics=true; //not to crash when some geometry is deleted while drawing
mesh = new IMesh().clr(0.7,0,0);
Cell cell = new Cell(new IVec(0,0,0), 10);
}

class Cell extends IParticle{
int growthDuration = 1000; //duration of growth and division
int growthInterval = 10;
int divisionInterval = 50;
double growthSpeed=0.1;

ArrayList< CellFace > faces; //store faces
boolean active = false;
IVertex vertex;

super(pos, new IVec(0,0,0));
faces = new ArrayList< CellFace >();
vertex = new IVertex();
vertex.pos = pos; // set reference of particle pos directly
fric(0.2);
}

void interact(ArrayList< IDynamics > agents){
IVec neighborCenter = new IVec(0,0,0);
int neighborCount=0;
for(int i=0; i < agents.size(); i++){
if(agents.get(i) instanceof Cell){
Cell c = (Cell)agents.get(i);
if(c != this){
// push if closer than two radii
IVec dif = c.pos().dif(pos());
c.push(dif);
}
// count neighbors and calculate their center
neighborCount++;
}
}
}
}
if(neighborCount >= 1){ // push from center of neighbors
neighborCenter.div(neighborCount);
IVec dif = pos().dif(neighborCenter).len(50); // constant force
this.push(dif);
}
}

void update(){
if(IG.time() < growthDuration){
if(time() > 0 && time()%divisionInterval==0){
if(active){ // divide when active flag is on
divide();
}
if(IRand.pct(50)){ // random division
active=true;
}
}
if(time()%growthInterval==0){
grow();
}
if(active){ clr(1.0,0.5,0); } // update color by being active
else{ clr(mesh.clr()); }
}
vertex.nml(nml()); // update mesh vertex normal
}

void grow(){ // growing cell size
}
}

void divide(){ // cell division
Cell child = createChild(IRand.dir());
}
Cell child = createChild(IRand.dir());
}
else{ // strip state

IVec dir = dividingLink.oppositeDir(this); // vector to the other cell
dir.projectToPlane(nml()); // division dir is projected on normal plane
Cell child = createChild(dir); // make a child toward the link direction

Cell c1 = f1.oppositeCell(dividingLink); // opposite cell from the link on the face
// delete 1 link + 1 face, create 3 links + 2 faces
new CellFace(this, c1, child);
new CellFace(child, c1, c0);
f1.del();
}
// delete 1 link + 1 face, create 4 links + 4 faces
new CellFace(this, c1, child);
new CellFace(this, child, c2);
new CellFace(child, c1, c0);
new CellFace(child, c0, c2);
f1.del();
f2.del();
}
}
}

IVec nml(){ // calc vertex normal from face normal
if(faces.size()==0){ return new IVec(0,0,1); }
IVec n = new IVec(0,0,0);
for(int i=0; i < faces.size(); i++){ n.add(faces.get(i).nml()); }
return n.unit();
}

Cell createChild(IVec dir){
radius *= 0.5; //make both cell size half
Cell child = new Cell(pos().cp(dir), radius);
pos().sub(dir); // move to the other way
active=false;  //reset activation
return child;
}
}

class CellLink extends IAgent{ // a link to connect 2 cells with spring force
double maxForce = 100;
Cell cell1, cell2;
ArrayList< CellFace > faces;

cell1 = c1; cell2 = c2;
faces = new ArrayList< CellFace >();
}

void interact(ArrayList< IDynamics > agents){
// spring force
IVec dif = cell1.pos().dif(cell2.pos());
if(force > maxForce){ force=maxForce; } // avoid putting too much force
else if(force < -maxForce){ force=-maxForce; }
dif.len(force);
cell1.pull(dif);
cell2.push(dif);
}

boolean contains(Cell c){ //check if link contains the cell
if(c==cell1 || c==cell2) return true;
return false;
}

void del(){
super.del(); // stop agent
}

Cell oppositeCell(Cell c){ // find other cell on the link
if(cell1==c) return cell2;
if(cell2==c) return cell1;
IG.err("Link does not contain the input cell");
return null;
}

IVec oppositeDir(Cell c){ //calculate a vector to the other cell
return oppositeCell(c).pos().dif(c.pos());
}
}

class CellFace{ // a triangle grouping 3 cells and 3 links
Cell cell1, cell2, cell3;
IFace face;

CellFace(Cell c1, Cell c2, Cell c3){
// keep order of cells in right hand screw for consistent normal
cell1 = c2; cell2 = c1; cell3 = c3;
}
else{
cell1 = c1; cell2 = c2; cell3 = c3;
}
face = new IFace(cell1.vertex,cell2.vertex,cell3.vertex);
}

IVec center(){ // calc center of triangle
return IVec.center(cell1.pos(), cell2.pos(), cell3.pos());
}

IVec nml(){ // normal vector
return cell1.pos().nml(cell2.pos(), cell3.pos()).unit();
}

// true if ordr of c1-c2 follows the order of cell1-cell2-cell3, otherwise false
boolean cellOrder(Cell c1, Cell c2){
if(cell1==c1 && cell2==c2 || cell2==c1 && cell3==c2 || cell3==c1 && cell1==c2 ){ return true; }
return false;
}

}

for(int i=0; i < c1.links.size() ; i++){
}
return null;
}

Cell oppositeCell(CellLink l){ // find a cell opposite of a link edge in a triangle
if(cell1!=l.cell1 && cell1!=l.cell2){ return cell1; }
if(cell2!=l.cell1 && cell2!=l.cell2){ return cell2; }
return cell3;
}

void del(){
cell1.faces.remove(this); // unregister this from cells and links
cell2.faces.remove(this);
cell3.faces.remove(this);
mesh.deleteFace(face);
}
}
```

### Closed Surface Formation and Face Division

There is another way to divide a triangular face by insertion of one new node of a cell other than dividing one edge and dividing a triangle into two. If you insert a new cell node inside of a triangle and add three new edges, you can divide one face into three. The code below implements this face division algorithm.

The code generates closed polygon mesh by forming a tetrahedron after cells form a triangle. Once the network of cells, links and faces form a tetrahedron, the division algorithm never makes any holes and keep the whole mesh surface closed.

```import igeo.*;
import processing.opengl.*;

IMesh mesh; // global variable of mesh geometry

void setup(){
size(480,360,IG.GL);
IConfig.syncDrawAndDynamics=true; //not to crash when some geometry is deleted while drawing
mesh = new IMesh().clr(0.7,0,0);
Cell cell = new Cell(new IVec(0,0,0), 10);
}

class Cell extends IParticle{
int growthDuration = 1000; //duration of growth and division
int growthInterval = 10;
int divisionInterval = 50;
double growthSpeed=0.1;

ArrayList< CellFace > faces; //store faces
boolean active = false;
IVertex vertex;

super(pos, new IVec(0,0,0));
faces = new ArrayList< CellFace >();
vertex = new IVertex();
vertex.pos = pos; // set reference of particle pos directly
fric(0.2);
}

void interact(ArrayList< IDynamics > agents){
IVec neighborCenter = new IVec(0,0,0);
int neighborCount=0;
for(int i=0; i < agents.size(); i++){
if(agents.get(i) instanceof Cell){
Cell c = (Cell)agents.get(i);
if(c != this){
// push if closer than two radii
IVec dif = c.pos().dif(pos());
c.push(dif);
}
// count neighbors and calculate their center
neighborCount++;
}
}
}
}
if(neighborCount >= 1){ // push from center of neighbors
neighborCenter.div(neighborCount);
IVec dif = pos().dif(neighborCenter).len(50); // constant force
this.push(dif);
}
}

void update(){
if(IG.time() < growthDuration){
if(time() > 0 && time()%divisionInterval==0){
if(active){ // divide when active flag is on
divide();
}
if(IRand.pct(50)){ // random division
active=true;
}
}
if(time()%growthInterval==0){
grow();
}
if(active){ clr(1.0,0.5,0); } // update color by being active
else{ clr(mesh.clr()); }
}
vertex.nml(nml()); // update mesh vertex normal
}

void grow(){ // growing cell size
}
}

void divide(){ // cell division
Cell child = createChild(IRand.dir());
}
// making a triangle loop
Cell child = createChild(IRand.dir());
}
// making a tetrahedron enclosure
Cell child = createChild(IRand.dir());
CellFace f = faces.get(0);
IVec center = IVec.center(child.pos(),f.cell1.pos(),f.cell2.pos(),f.cell3.pos());
if(f.center().dif(center).dot(f.nml())<0){ f.flipNml(); } // adjust normal to be outward
CellFace f1 = new CellFace(f.cell1, f.cell2, child);
CellFace f2 = new CellFace(f.cell2, f.cell3, child);
CellFace f3 = new CellFace(f.cell3, f.cell1, child);
if(f1.center().dif(center).dot(f1.nml())<0){ f1.flipNml(); }
if(f2.center().dif(center).dot(f2.nml())<0){ f2.flipNml(); }
if(f3.center().dif(center).dot(f3.nml())<0){ f3.flipNml(); }
}
// divide one face
CellFace dividingFace = faces.get(IRand.getInt(0,faces.size()-1));
Cell[] c = dividingFace.oppositeCell(this);
IVec dir = dividingFace.center().dif(pos());
dir.projectToPlane(nml());
Cell child = createChild(dir);
// delete 1 face, create 3 links + 3 faces
new CellFace(this, c[0], child);
new CellFace(this, child, c[1]);
new CellFace(child, c[0], c[1]);
dividingFace.del();
}
}

IVec nml(){ // calc vertex normal from face normal
if(faces.size()==0){ return new IVec(0,0,1); }
IVec n = new IVec(0,0,0);
for(int i=0; i < faces.size(); i++){ n.add(faces.get(i).nml()); }
return n.unit();
}

Cell createChild(IVec dir){
radius *= 0.5; //make both cell size half
Cell child = new Cell(pos().cp(dir), radius);
pos().sub(dir); // move to the other way
active=false;  //reset activation
return child;
}
}

class CellLink extends IAgent{ // a link to connect 2 cells with spring force
double maxForce = 100;
Cell cell1, cell2;
ArrayList< CellFace > faces;

cell1 = c1; cell2 = c2;
faces = new ArrayList< CellFace >();
}

void interact(ArrayList< IDynamics > agents){
// spring force
IVec dif = cell1.pos().dif(cell2.pos());
if(force > maxForce){ force=maxForce; } // avoid putting too much force
else if(force < -maxForce){ force=-maxForce; }
dif.len(force);
cell1.pull(dif);
cell2.push(dif);
}

boolean contains(Cell c){ //check if link contains the cell
if(c==cell1 || c==cell2) return true;
return false;
}

void del(){
super.del(); // stop agent
}

Cell oppositeCell(Cell c){ // find other cell on the link
if(cell1==c) return cell2;
if(cell2==c) return cell1;
IG.err("Link does not contain the input cell");
return null;
}

IVec oppositeDir(Cell c){ //calculate a vector to the other cell
return oppositeCell(c).pos().dif(c.pos());
}
}

class CellFace{ // a triangle grouping 3 cells and 3 links
Cell cell1, cell2, cell3;
IFace face;

CellFace(Cell c1, Cell c2, Cell c3){
// keep order of cells in right hand screw for consistent normal
cell1 = c2; cell2 = c1; cell3 = c3;
}
else{
cell1 = c1; cell2 = c2; cell3 = c3;
}
face = new IFace(cell1.vertex,cell2.vertex,cell3.vertex);
}

IVec center(){ // calc center of triangle
return IVec.center(cell1.pos(), cell2.pos(), cell3.pos());
}

IVec nml(){ // normal vector
return cell1.pos().nml(cell2.pos(), cell3.pos()).unit();
}

void flipNml(){ // flip normal
Cell tmp = cell2;
cell2 = cell3;
cell3 = tmp;
}

// true if ordr of c1-c2 follows the order of cell1-cell2-cell3, otherwise false
boolean cellOrder(Cell c1, Cell c2){
if(cell1==c1 && cell2==c2 || cell2==c1 && cell3==c2 || cell3==c1 && cell1==c2 ){ return true; }
return false;
}

}

for(int i=0; i < c1.links.size() ; i++){
}
return null;
}

Cell[] oppositeCell(Cell c){ // find 2 cells opposite of a cell in a triangle
if(cell1==c){ return  new Cell[]{ cell2, cell3 }; }
if(cell2==c){ return  new Cell[]{ cell3, cell1 };  }
if(cell3==c){ return  new Cell[]{ cell1, cell2 };  }
return null;
}

Cell oppositeCell(CellLink l){ // find a cell opposite of a link edge in a triangle
if(cell1!=l.cell1 && cell1!=l.cell2){ return cell1; }
if(cell2!=l.cell1 && cell2!=l.cell2){ return cell2; }
return cell3;
}

void del(){
cell1.faces.remove(this); // unregister this from cells and links
cell2.faces.remove(this);
cell3.faces.remove(this);
mesh.deleteFace(face);
}
}
```

### Interior Edge Division on Closed Surface

The previous face division algorithm tends to generate sharp corners. The code below uses the interior edge division instead of the face division, which is the same algorithm with the one in the previous section.

```import igeo.*;
import processing.opengl.*;

IMesh mesh; // global variable of mesh geometry

void setup(){
size(480,360,IG.GL);
IConfig.syncDrawAndDynamics=true; //not to crash when some geometry is deleted while drawing
mesh = new IMesh().clr(0.7,0,0);
Cell cell = new Cell(new IVec(0,0,0), 10);
}

class Cell extends IParticle{
int growthDuration = 1000; //duration of growth and division
int growthInterval = 10;
int divisionInterval = 50;
double growthSpeed=0.1;

ArrayList< CellFace > faces; //store faces
boolean active = false;
IVertex vertex;

super(pos, new IVec(0,0,0));
faces = new ArrayList< CellFace >();
vertex = new IVertex();
vertex.pos = pos; // set reference of particle pos directly
fric(0.2);
}

void interact(ArrayList< IDynamics > agents){
IVec neighborCenter = new IVec(0,0,0);
int neighborCount=0;
for(int i=0; i < agents.size(); i++){
if(agents.get(i) instanceof Cell){
Cell c = (Cell)agents.get(i);
if(c != this){
// push if closer than two radii
IVec dif = c.pos().dif(pos());
c.push(dif);
}
// count neighbors and calculate their center
neighborCount++;
}
}
}
}
if(neighborCount >= 1){ // push from center of neighbors
neighborCenter.div(neighborCount);
IVec dif = pos().dif(neighborCenter).len(50); // constant force
this.push(dif);
}
}

void update(){
if(IG.time() < growthDuration){
if(time() > 0 && time()%divisionInterval==0){
if(active){ // divide when active flag is on
divide();
}
if(IRand.pct(50)){ // random division
active=true;
}
}
if(time()%growthInterval==0){
grow();
}
if(active){ clr(1.0,0.5,0); } // update color by being active
else{ clr(mesh.clr()); }
}
vertex.nml(nml()); // update mesh vertex normal
}

void grow(){ // growing cell size
}
}

void divide(){ // cell division
Cell child = createChild(IRand.dir());
}
// making a triangle loop
Cell child = createChild(IRand.dir());
}
// making a tetrahedron enclosure
Cell child = createChild(IRand.dir());
CellFace f = faces.get(0);
IVec center = IVec.center(child.pos(),f.cell1.pos(),f.cell2.pos(),f.cell3.pos());
if(f.center().dif(center).dot(f.nml())<0){ f.flipNml(); } // adjust normal to be outward
CellFace f1 = new CellFace(f.cell1, f.cell2, child);
CellFace f2 = new CellFace(f.cell2, f.cell3, child);
CellFace f3 = new CellFace(f.cell3, f.cell1, child);
if(f1.center().dif(center).dot(f1.nml())<0){ f1.flipNml(); }
if(f2.center().dif(center).dot(f2.nml())<0){ f2.flipNml(); }
if(f3.center().dif(center).dot(f3.nml())<0){ f3.flipNml(); }
}
dir.projectToPlane(nml());
Cell child = createChild(dir);
// delete 1 link + 2 faces, create 4 links + 4 faces
new CellFace(this, c1, child);
new CellFace(this, child, c2);
new CellFace(child, c1, c0);
new CellFace(child, c0, c2);
f1.del();
f2.del();
}
}

IVec nml(){ // calc vertex normal from face normal
if(faces.size()==0){ return new IVec(0,0,1); }
IVec n = new IVec(0,0,0);
for(int i=0; i < faces.size(); i++){ n.add(faces.get(i).nml()); }
return n.unit();
}

Cell createChild(IVec dir){
radius *= 0.5; //make both cell size half
Cell child = new Cell(pos().cp(dir), radius);
pos().sub(dir); // move to the other way
active=false;  //reset activation
return child;
}
}

class CellLink extends IAgent{ // a link to connect 2 cells with spring force
double maxForce = 100;
Cell cell1, cell2;
ArrayList< CellFace > faces;

cell1 = c1; cell2 = c2;
faces = new ArrayList< CellFace >();
}

void interact(ArrayList< IDynamics > agents){
// spring force
IVec dif = cell1.pos().dif(cell2.pos());
if(force > maxForce){ force=maxForce; } // avoid putting too much force
else if(force < -maxForce){ force=-maxForce; }
dif.len(force);
cell1.pull(dif);
cell2.push(dif);
}

boolean contains(Cell c){ //check if link contains the cell
if(c==cell1 || c==cell2) return true;
return false;
}

void del(){
super.del(); // stop agent
}

Cell oppositeCell(Cell c){ // find other cell on the link
if(cell1==c) return cell2;
if(cell2==c) return cell1;
IG.err("Link does not contain the input cell");
return null;
}

IVec oppositeDir(Cell c){ //calculate a vector to the other cell
return oppositeCell(c).pos().dif(c.pos());
}
}

class CellFace{ // a triangle grouping 3 cells and 3 links
Cell cell1, cell2, cell3;
IFace face;

CellFace(Cell c1, Cell c2, Cell c3){
// keep order of cells in right hand screw for consistent normal
cell1 = c2; cell2 = c1; cell3 = c3;
}
else{
cell1 = c1; cell2 = c2; cell3 = c3;
}
face = new IFace(cell1.vertex,cell2.vertex,cell3.vertex);
}

IVec center(){ // calc center of triangle
return IVec.center(cell1.pos(), cell2.pos(), cell3.pos());
}

IVec nml(){ // normal vector
return cell1.pos().nml(cell2.pos(), cell3.pos()).unit();
}

void flipNml(){ // flip normal
Cell tmp = cell2;
cell2 = cell3;
cell3 = tmp;
}

// true if ordr of c1-c2 follows the order of cell1-cell2-cell3, otherwise false
boolean cellOrder(Cell c1, Cell c2){
if(cell1==c1 && cell2==c2 || cell2==c1 && cell3==c2 || cell3==c1 && cell1==c2 ){ return true; }
return false;
}

}

for(int i=0; i < c1.links.size() ; i++){
}
return null;
}

Cell[] oppositeCell(Cell c){ // find 2 cells opposite of a cell in a triangle
if(cell1==c){ return  new Cell[]{ cell2, cell3 }; }
if(cell2==c){ return  new Cell[]{ cell3, cell1 };  }
if(cell3==c){ return  new Cell[]{ cell1, cell2 };  }
return null;
}

Cell oppositeCell(CellLink l){ // find a cell opposite of a link edge in a triangle
if(cell1!=l.cell1 && cell1!=l.cell2){ return cell1; }
if(cell2!=l.cell1 && cell2!=l.cell2){ return cell2; }
return cell3;
}

void del(){
cell1.faces.remove(this); // unregister this from cells and links
cell2.faces.remove(this);
cell3.faces.remove(this);
mesh.deleteFace(face);
}
}
```

### Finding Edges to Divide on Closed Surface

The face division and the edge division are not only way to divide a cell and recreate links and connection. To integrate different ways, the next codes implements an algorithm to select two edges on a dividing cell and insert a new child cell, three new edges and two new faces. The way of division is changed by the interval of two selected edges. If the interval is one, it's equivalent to the face division. If the interval is two, it's equivalent to the edge division.

This division algorithm is organized in several methods in addition to divide method. The two edges for division is selected at findDividingLinks method. To check the interval of two edges, the code needs to know which edge is next to each other and order of edges. For this purpose, sortLinks methods sorts the stored CellLink instances. insertChild method reconnect links and faces with a new child cell by deleting existing links and faces and creating new ones.

```import igeo.*;
import processing.opengl.*;

IMesh mesh; // global variable of mesh geometry

void setup(){
size(480,360,IG.GL);
IConfig.syncDrawAndDynamics=true; //not to crash when some geometry is deleted while drawing
mesh = new IMesh().clr(0.7,0,0);
Cell cell = new Cell(new IVec(0,0,0), 10);
}

class Cell extends IParticle{
int growthDuration = 1200; //duration of growth and division
int growthInterval = 10;
int divisionInterval = 50;
double growthSpeed=0.1;

ArrayList< CellFace > faces; //store faces
boolean active = false;
IVertex vertex;

super(pos, new IVec(0,0,0));
faces = new ArrayList< CellFace >();
vertex = new IVertex();
vertex.pos = pos; // set reference of particle pos directly
fric(0.2);
}

void interact(ArrayList< IDynamics > agents){
IVec neighborCenter = new IVec(0,0,0);
int neighborCount=0;
for(int i=0; i < agents.size(); i++){
if(agents.get(i) instanceof Cell){
Cell c = (Cell)agents.get(i);
if(c != this){
// push if closer than two radii
IVec dif = c.pos().dif(pos());
c.push(dif);
}
// count neighbors and calculate their center
neighborCount++;
}
}
}
}
if(neighborCount >= 1){ // push from center of neighbors
neighborCenter.div(neighborCount);
IVec dif = pos().dif(neighborCenter).len(50); // constant force
this.push(dif);
}
this.push(nml().len(20)); // pressure toward normal
}

void update(){
if(IG.time() < growthDuration){
if(time() > 0 && time()%divisionInterval==0){
if(active){ // divide when active flag is on
divide();
}
if(IRand.pct(50)){ // random division
active=true;
}
}
if(time()%growthInterval==0){
grow();
}
if(active){ clr(1.0,0.5,0); } // update color by being active
else{ clr(mesh.clr()); }
}
vertex.nml(nml()); // update mesh vertex normal
}

void grow(){ // growing cell size
}
}

void divide(){ // cell division
Cell child = createChild(IRand.dir());
}
// making a triangle loop
Cell child = createChild(IRand.dir());
}
// making a tetrahedron enclosure
Cell child = createChild(IRand.dir());
CellFace f = faces.get(0);
IVec center = IVec.center(child.pos(),f.cell1.pos(),f.cell2.pos(),f.cell3.pos());
if(f.center().dif(center).dot(f.nml())<0){ f.flipNml(); } // adjust normal to be outward
CellFace f1 = new CellFace(f.cell1, f.cell2, child);
CellFace f2 = new CellFace(f.cell2, f.cell3, child);
CellFace f3 = new CellFace(f.cell3, f.cell1, child);
if(f1.center().dif(center).dot(f1.nml())<0){ f1.flipNml(); }
if(f2.center().dif(center).dot(f2.nml())<0){ f2.flipNml(); }
if(f3.center().dif(center).dot(f3.nml())<0){ f3.flipNml(); }
}
// pick two links and insert two faces between
}
IVec dir = new IVec();
}
dir.projectToPlane(nml());
Cell child = createChild(dir);
}
}

IVec nml(){ // calc vertex normal from face normal
if(faces.size()==0){ return new IVec(0,0,1); }
IVec n = new IVec(0,0,0);
for(int i=0; i < faces.size(); i++){ n.add(faces.get(i).nml()); }
return n.unit();
}

Cell createChild(IVec dir){
radius *= 0.5; //make both cell size half
Cell child = new Cell(pos().cp(dir), radius);
pos().sub(dir); // move to the other way
active=false;  //reset activation
return child;
}

if(links.size() <= 3) return; // no need to sort
for(int i=0; i < links.size(); i++){
}
}

int[] findDividingLinks(){ // find splittable two links and return the index numbers

return new int[]{ idx, idx+linkInterval }; // index can be larger than links.size()
}

ArrayList< CellFace > removingFaces = new ArrayList< CellFace >();
for(int i=0; i < faces.size(); i++){
for(int j=0; j < linksToChild.size()-1; j++){
break;
}
}
}
for(int i=1; i < linksToChild.size()-1; i++){
}
for(int i=0; i < removingFaces.size(); i++){
removingFaces.get(i).del(); // replace face by deleting and recreating
Cell[] cells = removingFaces.get(i).oppositeCell(this);
new CellFace(child, cells[0], cells[1]);
}
new CellFace(this, cell1, child);
new CellFace(this, child, cell2);
}
}

class CellLink extends IAgent{ // a link to connect 2 cells with spring force
double maxForce = 100;
Cell cell1, cell2;
ArrayList< CellFace > faces;

cell1 = c1; cell2 = c2;
faces = new ArrayList< CellFace >();
}

void interact(ArrayList< IDynamics > agents){
// spring force
IVec dif = cell1.pos().dif(cell2.pos());
if(force > maxForce){ force=maxForce; } // avoid putting too much force
else if(force < -maxForce){ force=-maxForce; }
dif.len(force);
cell1.pull(dif);
cell2.push(dif);
}

boolean contains(Cell c){ //check if link contains the cell
if(c==cell1 || c==cell2) return true;
return false;
}

void del(){
super.del(); // stop agent
}

CellFace oppositeFace(CellFace f){ // find other face on the link
if(faces.size()!=2) return null;
if(faces.get(0)==f) return faces.get(1);
return faces.get(0);
}

Cell oppositeCell(Cell c){ // find other cell on the link
if(cell1==c) return cell2;
if(cell2==c) return cell1;
IG.err("Link does not contain the input cell");
return null;
}

IVec oppositeDir(Cell c){ //calculate a vector to the other cell
return oppositeCell(c).pos().dif(c.pos());
}
}

class CellFace{ // a triangle grouping 3 cells and 3 links
Cell cell1, cell2, cell3;
IFace face;

CellFace(Cell c1, Cell c2, Cell c3){
// keep order of cells in right hand screw for consistent normal
cell1 = c2; cell2 = c1; cell3 = c3;
}
else{
cell1 = c1; cell2 = c2; cell3 = c3;
}
face = new IFace(cell1.vertex,cell2.vertex,cell3.vertex);
}

IVec center(){ // calc center of triangle
return IVec.center(cell1.pos(), cell2.pos(), cell3.pos());
}

IVec nml(){ // normal vector
return cell1.pos().nml(cell2.pos(), cell3.pos()).unit();
}

void flipNml(){ // flip normal
Cell tmp = cell2;
cell2 = cell3;
cell3 = tmp;
}

// true if ordr of c1-c2 follows the order of cell1-cell2-cell3, otherwise false
boolean cellOrder(Cell c1, Cell c2){
if(cell1==c1 && cell2==c2 || cell2==c1 && cell3==c2 || cell3==c1 && cell1==c2 ){ return true; }
return false;
}

}

for(int i=0; i < c1.links.size() ; i++){
}
return null;
}

Cell[] oppositeCell(Cell c){ // find 2 cells opposite of a cell in a triangle
if(cell1==c){ return  new Cell[]{ cell2, cell3 }; }
if(cell2==c){ return  new Cell[]{ cell3, cell1 };  }
if(cell3==c){ return  new Cell[]{ cell1, cell2 };  }
return null;
}

return null;
}

void del(){
cell1.faces.remove(this); // unregister this from cells and links
cell2.faces.remove(this);
cell3.faces.remove(this);
mesh.deleteFace(face);
}
}
```

### Division Control by Maximum Link Count

The next code brings a rule to limit number of links in a cell. The number of links is checked in findDividingLinks method not only with the dividing cell but also with adjacent cells which would increase links if new links and faces are inserted. findDividingLinks also searches a pair of edges which satisfy the maximum link number condition by checking all possible pairs with the specified link interval until it finds one.

```import igeo.*;
import processing.opengl.*;

IMesh mesh; // global variable of mesh geometry

void setup(){
size(480,360,IG.GL);
IConfig.syncDrawAndDynamics=true; //not to crash when some geometry is deleted while drawing
mesh = new IMesh().clr(0.7,0,0);
Cell cell = new Cell(new IVec(0,0,0), 10);
}

class Cell extends IParticle{
int growthDuration = 20000; //duration of growth and division
int growthInterval = 10;
int divisionInterval = 50;
double growthSpeed=0.1;

ArrayList< CellFace > faces; //store faces
boolean active = false;
IVertex vertex;

super(pos, new IVec(0,0,0));
faces = new ArrayList< CellFace >();
vertex = new IVertex();
vertex.pos = pos; // set reference of particle pos directly
fric(0.2);
}

void interact(ArrayList< IDynamics > agents){
IVec neighborCenter = new IVec(0,0,0);
int neighborCount=0;
for(int i=0; i < agents.size(); i++){
if(agents.get(i) instanceof Cell){
Cell c = (Cell)agents.get(i);
if(c != this){
// push if closer than two radii
IVec dif = c.pos().dif(pos());
c.push(dif);
}
// count neighbors and calculate their center
neighborCount++;
}
}
}
}
if(neighborCount >= 1){ // push from center of neighbors
neighborCenter.div(neighborCount);
IVec dif = pos().dif(neighborCenter).len(50); // constant force
this.push(dif);
}
this.push(nml().len(20)); // pressure toward normal
}

void update(){
if(IG.time() < growthDuration){
if(time() > 0 && time()%divisionInterval==0){
if(active){ // divide when active flag is on
divide();
}
if(IRand.pct(50)){ // random division
active=true;
}
}
if(time()%growthInterval==0){
grow();
}
if(active){ clr(1.0,0.5,0); } // update color by being active
else{ clr(mesh.clr()); }
}
vertex.nml(nml()); // update mesh vertex normal
}

void grow(){ // growing cell size
}
}

void divide(){ // cell division
Cell child = createChild(IRand.dir());
}
// making a triangle loop
Cell child = createChild(IRand.dir());
}
// making a tetrahedron enclosure
Cell child = createChild(IRand.dir());
CellFace f = faces.get(0);
IVec center = IVec.center(child.pos(),f.cell1.pos(),f.cell2.pos(),f.cell3.pos());
if(f.center().dif(center).dot(f.nml())<0){ f.flipNml(); } // adjust normal to be outward
CellFace f1 = new CellFace(f.cell1, f.cell2, child);
CellFace f2 = new CellFace(f.cell2, f.cell3, child);
CellFace f3 = new CellFace(f.cell3, f.cell1, child);
if(f1.center().dif(center).dot(f1.nml())<0){ f1.flipNml(); }
if(f2.center().dif(center).dot(f2.nml())<0){ f2.flipNml(); }
if(f3.center().dif(center).dot(f3.nml())<0){ f3.flipNml(); }
}
// pick two links and insert two faces between
if(linkIdx==null){ // no valid edge to split. skip division.
active = false;
return;
}
}
IVec dir = new IVec();
}
dir.projectToPlane(nml());
Cell child = createChild(dir);
}
}

IVec nml(){ // calc vertex normal from face normal
if(faces.size()==0){ return new IVec(0,0,1); }
IVec n = new IVec(0,0,0);
for(int i=0; i < faces.size(); i++){ n.add(faces.get(i).nml()); }
return n.unit();
}

Cell createChild(IVec dir){
radius *= 0.5; //make both cell size half
Cell child = new Cell(pos().cp(dir), radius);
pos().sub(dir); // move to the other way
active=false;  //reset activation
return child;
}

if(links.size() <= 3) return; // no need to sort
for(int i=0; i < links.size(); i++){
}
}

int[] findDividingLinks(){ // find splittable two links and return the index numbers

return null;
}
for(int i=0; i < links.size(); i++){ // check all pairs
return new int[]{ i+idx, i+idx+linkInterval }; // index can be larger than links.size()
}
}
return null;
}

ArrayList< CellFace > removingFaces = new ArrayList< CellFace >();
for(int i=0; i < faces.size(); i++){
for(int j=0; j < linksToChild.size()-1; j++){
break;
}
}
}
for(int i=1; i < linksToChild.size()-1; i++){
}
for(int i=0; i < removingFaces.size(); i++){
removingFaces.get(i).del(); // replace face by deleting and recreating
Cell[] cells = removingFaces.get(i).oppositeCell(this);
new CellFace(child, cells[0], cells[1]);
}
new CellFace(this, cell1, child);
new CellFace(this, child, cell2);
}
}

class CellLink extends IAgent{ // a link to connect 2 cells with spring force
double maxForce = 100;
Cell cell1, cell2;
ArrayList< CellFace > faces;

cell1 = c1; cell2 = c2;
faces = new ArrayList< CellFace >();
}

void interact(ArrayList< IDynamics > agents){
// spring force
IVec dif = cell1.pos().dif(cell2.pos());
if(force > maxForce){ force=maxForce; } // avoid putting too much force
else if(force < -maxForce){ force=-maxForce; }
dif.len(force);
cell1.pull(dif);
cell2.push(dif);
}

boolean contains(Cell c){ //check if link contains the cell
if(c==cell1 || c==cell2) return true;
return false;
}

void del(){
super.del(); // stop agent
}

CellFace oppositeFace(CellFace f){ // find other face on the link
if(faces.size()!=2) return null;
if(faces.get(0)==f) return faces.get(1);
return faces.get(0);
}

Cell oppositeCell(Cell c){ // find other cell on the link
if(cell1==c) return cell2;
if(cell2==c) return cell1;
IG.err("Link does not contain the input cell");
return null;
}

IVec oppositeDir(Cell c){ //calculate a vector to the other cell
return oppositeCell(c).pos().dif(c.pos());
}
}

class CellFace{ // a triangle grouping 3 cells and 3 links
Cell cell1, cell2, cell3;
IFace face;

CellFace(Cell c1, Cell c2, Cell c3){
// keep order of cells in right hand screw for consistent normal
cell1 = c2; cell2 = c1; cell3 = c3;
}
else{
cell1 = c1; cell2 = c2; cell3 = c3;
}
face = new IFace(cell1.vertex,cell2.vertex,cell3.vertex);
}

IVec center(){ // calc center of triangle
return IVec.center(cell1.pos(), cell2.pos(), cell3.pos());
}

IVec nml(){ // normal vector
return cell1.pos().nml(cell2.pos(), cell3.pos()).unit();
}

void flipNml(){ // flip normal
Cell tmp = cell2;
cell2 = cell3;
cell3 = tmp;
}

// true if ordr of c1-c2 follows the order of cell1-cell2-cell3, otherwise false
boolean cellOrder(Cell c1, Cell c2){
if(cell1==c1 && cell2==c2 || cell2==c1 && cell3==c2 || cell3==c1 && cell1==c2 ){ return true; }
return false;
}

}

for(int i=0; i < c1.links.size() ; i++){
}
return null;
}

Cell[] oppositeCell(Cell c){ // find 2 cells opposite of a cell in a triangle
if(cell1==c){ return  new Cell[]{ cell2, cell3 }; }
if(cell2==c){ return  new Cell[]{ cell3, cell1 };  }
if(cell3==c){ return  new Cell[]{ cell1, cell2 };  }
return null;
}

return null;
}

void del(){
cell1.faces.remove(this); // unregister this from cells and links
cell2.faces.remove(this);
cell3.faces.remove(this);
mesh.deleteFace(face);
}
}
```

Although the limitation of 6 maximum links tends to bring a certain regularity in the surface formation, it also tends to limit the growth because number of links can be easily maxed out and a dividable pair of edges cannot be found. To ease the condition, the next code accepts one more link randomly with a certain probability (extraLinkAllowance field). This change happens in findDividingLinks method.

```import igeo.*;
import processing.opengl.*;

IMesh mesh; // global variable of mesh geometry

void setup(){
size(480,360,IG.GL);
IConfig.syncDrawAndDynamics=true; //not to crash when some geometry is deleted while drawing
mesh = new IMesh().clr(0.7,0,0);
Cell cell = new Cell(new IVec(0,0,0), 10);
}

class Cell extends IParticle{
int growthDuration = 2000; //duration of growth and division
int growthInterval = 10;
int divisionInterval = 50;
double growthSpeed=0.1;

ArrayList< CellFace > faces; //store faces
boolean active = false;
IVertex vertex;

super(pos, new IVec(0,0,0));
faces = new ArrayList< CellFace >();
vertex = new IVertex();
vertex.pos = pos; // set reference of particle pos directly
fric(0.2);
}

void interact(ArrayList< IDynamics > agents){
IVec neighborCenter = new IVec(0,0,0);
int neighborCount=0;
for(int i=0; i < agents.size(); i++){
if(agents.get(i) instanceof Cell){
Cell c = (Cell)agents.get(i);
if(c != this){
// push if closer than two radii
IVec dif = c.pos().dif(pos());
c.push(dif);
}
// count neighbors and calculate their center
neighborCount++;
}
}
}
}
if(neighborCount >= 1){ // push from center of neighbors
neighborCenter.div(neighborCount);
IVec dif = pos().dif(neighborCenter).len(50); // constant force
this.push(dif);
}
this.push(nml().len(20)); // pressure toward normal
}

void update(){
if(IG.time() < growthDuration){
if(time() > 0 && time()%divisionInterval==0){
if(active){ // divide when active flag is on
divide();
}
if(IRand.pct(50)){ // random division
active=true;
}
}
if(time()%growthInterval==0){
grow();
}
if(active){ clr(1.0,0.5,0); } // update color by being active
else{ clr(mesh.clr()); }
}
vertex.nml(nml()); // update mesh vertex normal
}

void grow(){ // growing cell size
}
}

void divide(){ // cell division
Cell child = createChild(IRand.dir());
}
// making a triangle loop
Cell child = createChild(IRand.dir());
}
// making a tetrahedron enclosure
Cell child = createChild(IRand.dir());
CellFace f = faces.get(0);
IVec center = IVec.center(child.pos(),f.cell1.pos(),f.cell2.pos(),f.cell3.pos());
if(f.center().dif(center).dot(f.nml())<0){ f.flipNml(); } // adjust normal to be outward
CellFace f1 = new CellFace(f.cell1, f.cell2, child);
CellFace f2 = new CellFace(f.cell2, f.cell3, child);
CellFace f3 = new CellFace(f.cell3, f.cell1, child);
if(f1.center().dif(center).dot(f1.nml())<0){ f1.flipNml(); }
if(f2.center().dif(center).dot(f2.nml())<0){ f2.flipNml(); }
if(f3.center().dif(center).dot(f3.nml())<0){ f3.flipNml(); }
}
// pick two links and insert two faces between
if(linkIdx==null){ // no valid edge to split. skip division.
active = false;
return;
}
}
IVec dir = new IVec();
}
dir.projectToPlane(nml());
Cell child = createChild(dir);
}
}

IVec nml(){ // calc vertex normal from face normal
if(faces.size()==0){ return new IVec(0,0,1); }
IVec n = new IVec(0,0,0);
for(int i=0; i < faces.size(); i++){ n.add(faces.get(i).nml()); }
return n.unit();
}

Cell createChild(IVec dir){
radius *= 0.5; //make both cell size half
Cell child = new Cell(pos().cp(dir), radius);
pos().sub(dir); // move to the other way
active=false;  //reset activation
return child;
}

if(links.size() <= 3) return; // no need to sort
for(int i=0; i < links.size(); i++){
}
}

int[] findDividingLinks(){ // find splittable two links and return the index numbers

return null;
}
for(int i=0; i < links.size(); i++){ // check all pairs
return new int[]{ i+idx, i+idx+linkInterval }; // index can be larger than links.size()
}
}
return null;
}

ArrayList< CellFace > removingFaces = new ArrayList< CellFace >();
for(int i=0; i < faces.size(); i++){
for(int j=0; j < linksToChild.size()-1; j++){
break;
}
}
}
for(int i=1; i < linksToChild.size()-1; i++){
}
for(int i=0; i < removingFaces.size(); i++){
removingFaces.get(i).del(); // replace face by deleting and recreating
Cell[] cells = removingFaces.get(i).oppositeCell(this);
new CellFace(child, cells[0], cells[1]);
}
new CellFace(this, cell1, child);
new CellFace(this, child, cell2);
}
}

class CellLink extends IAgent{ // a link to connect 2 cells with spring force
double maxForce = 100;
Cell cell1, cell2;
ArrayList< CellFace > faces;

cell1 = c1; cell2 = c2;
faces = new ArrayList< CellFace >();
}

void interact(ArrayList< IDynamics > agents){
// spring force
IVec dif = cell1.pos().dif(cell2.pos());
if(force > maxForce){ force=maxForce; } // avoid putting too much force
else if(force < -maxForce){ force=-maxForce; }
dif.len(force);
cell1.pull(dif);
cell2.push(dif);
}

boolean contains(Cell c){ //check if link contains the cell
if(c==cell1 || c==cell2) return true;
return false;
}

void del(){
super.del(); // stop agent
}

CellFace oppositeFace(CellFace f){ // find other face on the link
if(faces.size()!=2) return null;
if(faces.get(0)==f) return faces.get(1);
return faces.get(0);
}

Cell oppositeCell(Cell c){ // find other cell on the link
if(cell1==c) return cell2;
if(cell2==c) return cell1;
IG.err("Link does not contain the input cell");
return null;
}

IVec oppositeDir(Cell c){ //calculate a vector to the other cell
return oppositeCell(c).pos().dif(c.pos());
}
}

class CellFace{ // a triangle grouping 3 cells and 3 links
Cell cell1, cell2, cell3;
IFace face;

CellFace(Cell c1, Cell c2, Cell c3){
// keep order of cells in right hand screw for consistent normal
cell1 = c2; cell2 = c1; cell3 = c3;
}
else{
cell1 = c1; cell2 = c2; cell3 = c3;
}
face = new IFace(cell1.vertex,cell2.vertex,cell3.vertex);
}

IVec center(){ // calc center of triangle
return IVec.center(cell1.pos(), cell2.pos(), cell3.pos());
}

IVec nml(){ // normal vector
return cell1.pos().nml(cell2.pos(), cell3.pos()).unit();
}

void flipNml(){ // flip normal
Cell tmp = cell2;
cell2 = cell3;
cell3 = tmp;
}

// true if ordr of c1-c2 follows the order of cell1-cell2-cell3, otherwise false
boolean cellOrder(Cell c1, Cell c2){
if(cell1==c1 && cell2==c2 || cell2==c1 && cell3==c2 || cell3==c1 && cell1==c2 ){ return true; }
return false;
}

}

for(int i=0; i < c1.links.size() ; i++){
}
return null;
}

Cell[] oppositeCell(Cell c){ // find 2 cells opposite of a cell in a triangle
if(cell1==c){ return  new Cell[]{ cell2, cell3 }; }
if(cell2==c){ return  new Cell[]{ cell3, cell1 };  }
if(cell3==c){ return  new Cell[]{ cell1, cell2 };  }
return null;
}

return null;
}

void del(){
cell1.faces.remove(this); // unregister this from cells and links
cell2.faces.remove(this);
cell3.faces.remove(this);
mesh.deleteFace(face);
}
}
```

### Cascading Activation Control on Division

The next three codes below bring a type of division control by activation parameter (active field in the Cell class) and controlling it by a parent cell when it creates a new child cell in division method. The next code use the cascading logic to activate a child cell and deactivate the parent cell after it finishes division, in the same way with the example in the previous page.

```import igeo.*;
import processing.opengl.*;

IMesh mesh; // global variable of mesh geometry

void setup(){
size(480,360,IG.GL);
IConfig.syncDrawAndDynamics=true; //not to crash when some geometry is deleted while drawing
mesh = new IMesh().clr(0.7,0,0);
Cell cell = new Cell(new IVec(0,0,0), 10);
cell.active = true;
}

class Cell extends IParticle{
int growthDuration = 10000; //duration of growth and division
int growthInterval = 10;
int divisionInterval = 20;
double growthSpeed=0.1;

ArrayList< CellFace > faces; //store faces
boolean active = false;
IVertex vertex;

super(pos, new IVec(0,0,0));
faces = new ArrayList< CellFace >();
vertex = new IVertex();
vertex.pos = pos; // set reference of particle pos directly
fric(0.2);
}

void interact(ArrayList< IDynamics > agents){
IVec neighborCenter = new IVec(0,0,0);
int neighborCount=0;
for(int i=0; i < agents.size(); i++){
if(agents.get(i) instanceof Cell){
Cell c = (Cell)agents.get(i);
if(c != this){
// push if closer than two radii
IVec dif = c.pos().dif(pos());
c.push(dif);
}
// count neighbors and calculate their center
neighborCount++;
}
}
}
}
if(neighborCount >= 1){ // push from center of neighbors
neighborCenter.div(neighborCount);
IVec dif = pos().dif(neighborCenter).len(100); // constant force
this.push(dif);
}
this.push(nml().len(20)); // pressure toward normal
}

void update(){
if(IG.time() < growthDuration){
if(time() > 0 && time()%divisionInterval==0){
if(active){ // divide when active flag is on
divide();
}
}
if(time()%growthInterval==0){
grow();
}
if(active){ clr(1.0,0.5,0); } // update color by being active
else{ clr(mesh.clr()); }
}
vertex.nml(nml()); // update mesh vertex normal
}

void grow(){ // growing cell size
}
}

void divide(){ // cell division
Cell child = createChild(IRand.dir());
}
// making a triangle loop
Cell child = createChild(IRand.dir());
}
// making a tetrahedron enclosure
Cell child = createChild(IRand.dir());
CellFace f = faces.get(0);
IVec center = IVec.center(child.pos(),f.cell1.pos(),f.cell2.pos(),f.cell3.pos());
if(f.center().dif(center).dot(f.nml())<0){ f.flipNml(); } // adjust normal to be outward
CellFace f1 = new CellFace(f.cell1, f.cell2, child);
CellFace f2 = new CellFace(f.cell2, f.cell3, child);
CellFace f3 = new CellFace(f.cell3, f.cell1, child);
if(f1.center().dif(center).dot(f1.nml())<0){ f1.flipNml(); }
if(f2.center().dif(center).dot(f2.nml())<0){ f2.flipNml(); }
if(f3.center().dif(center).dot(f3.nml())<0){ f3.flipNml(); }
}
// pick two links and insert two faces between
if(linkIdx==null){ // no valid edge to split. skip division.
active = false;
return;
}
}
IVec dir = new IVec();
}
dir.projectToPlane(nml());
Cell child = createChild(dir);
}
}

IVec nml(){ // calc vertex normal from face normal
if(faces.size()==0){ return new IVec(0,0,1); }
IVec n = new IVec(0,0,0);
for(int i=0; i < faces.size(); i++){ n.add(faces.get(i).nml()); }
return n.unit();
}

Cell createChild(IVec dir){
radius *= 0.5; //make both cell size half
Cell child = new Cell(pos().cp(dir), radius);
pos().sub(dir); // move to the other way
active=false; //reset activation
child.active = true; //activate child always
return child;
}

if(links.size() <= 3) return; // no need to sort
for(int i=0; i < links.size(); i++){
}
}

int[] findDividingLinks(){ // find splittable two links and return the index numbers
int linkInterval = 1; // face division only

return null;
}
for(int i=0; i < links.size(); i++){ // check all pairs
return new int[]{ i+idx, i+idx+linkInterval }; // index can be larger than links.size()
}
}
return null;
}

ArrayList< CellFace > removingFaces = new ArrayList< CellFace >();
for(int i=0; i < faces.size(); i++){
for(int j=0; j < linksToChild.size()-1; j++){
break;
}
}
}
for(int i=1; i < linksToChild.size()-1; i++){
}
for(int i=0; i < removingFaces.size(); i++){
removingFaces.get(i).del(); // replace face by deleting and recreating
Cell[] cells = removingFaces.get(i).oppositeCell(this);
new CellFace(child, cells[0], cells[1]);
}
new CellFace(this, cell1, child);
new CellFace(this, child, cell2);
}
}

class CellLink extends IAgent{ // a link to connect 2 cells with spring force
double maxForce = 100;
Cell cell1, cell2;
ArrayList< CellFace > faces;

cell1 = c1; cell2 = c2;
faces = new ArrayList< CellFace >();
}

void interact(ArrayList< IDynamics > agents){
// spring force
IVec dif = cell1.pos().dif(cell2.pos());
if(force > maxForce){ force=maxForce; } // avoid putting too much force
else if(force < -maxForce){ force=-maxForce; }
dif.len(force);
cell1.pull(dif);
cell2.push(dif);
}

boolean contains(Cell c){ //check if link contains the cell
if(c==cell1 || c==cell2) return true;
return false;
}

void del(){
super.del(); // stop agent
}

CellFace oppositeFace(CellFace f){ // find other face on the link
if(faces.size()!=2) return null;
if(faces.get(0)==f) return faces.get(1);
return faces.get(0);
}

Cell oppositeCell(Cell c){ // find other cell on the link
if(cell1==c) return cell2;
if(cell2==c) return cell1;
IG.err("Link does not contain the input cell");
return null;
}

IVec oppositeDir(Cell c){ //calculate a vector to the other cell
return oppositeCell(c).pos().dif(c.pos());
}
}

class CellFace{ // a triangle grouping 3 cells and 3 links
Cell cell1, cell2, cell3;
IFace face;

CellFace(Cell c1, Cell c2, Cell c3){
// keep order of cells in right hand screw for consistent normal
cell1 = c2; cell2 = c1; cell3 = c3;
}
else{
cell1 = c1; cell2 = c2; cell3 = c3;
}
face = new IFace(cell1.vertex,cell2.vertex,cell3.vertex);
}

IVec center(){ // calc center of triangle
return IVec.center(cell1.pos(), cell2.pos(), cell3.pos());
}

IVec nml(){ // normal vector
return cell1.pos().nml(cell2.pos(), cell3.pos()).unit();
}

void flipNml(){ // flip normal
Cell tmp = cell2;
cell2 = cell3;
cell3 = tmp;
}

// true if ordr of c1-c2 follows the order of cell1-cell2-cell3, otherwise false
boolean cellOrder(Cell c1, Cell c2){
if(cell1==c1 && cell2==c2 || cell2==c1 && cell3==c2 || cell3==c1 && cell1==c2 ){ return true; }
return false;
}

}

for(int i=0; i < c1.links.size() ; i++){
}
return null;
}

Cell[] oppositeCell(Cell c){ // find 2 cells opposite of a cell in a triangle
if(cell1==c){ return  new Cell[]{ cell2, cell3 }; }
if(cell2==c){ return  new Cell[]{ cell3, cell1 };  }
if(cell3==c){ return  new Cell[]{ cell1, cell2 };  }
return null;
}

return null;
}

void del(){
cell1.faces.remove(this); // unregister this from cells and links
cell2.faces.remove(this);
cell3.faces.remove(this);
mesh.deleteFace(face);
}
}
```

### Branching Activation Control on Division

The division algorithm in the next code activates a child cell and also keeps the parent cell activated to cause branching, like the example in the previous page. Parent cells are maintained to be active only in 10% probability (90% deactivated) in the code not to generate exponentially many branches.

```import igeo.*;
import processing.opengl.*;

IMesh mesh; // global variable of mesh geometry

void setup(){
size(480,360,IG.GL);
IConfig.syncDrawAndDynamics=true; //not to crash when some geometry is deleted while drawing
mesh = new IMesh().clr(0.7,0,0);
Cell cell = new Cell(new IVec(0,0,0), 10);
cell.active = true;
}

class Cell extends IParticle{
int growthDuration = 2000; //duration of growth and division
int growthInterval = 10;
int divisionInterval = 20;
double growthSpeed=0.1;

ArrayList< CellFace > faces; //store faces
boolean active = false;
IVertex vertex;

super(pos, new IVec(0,0,0));
faces = new ArrayList< CellFace >();
vertex = new IVertex();
vertex.pos = pos; // set reference of particle pos directly
fric(0.2);
}

void interact(ArrayList< IDynamics > agents){
IVec neighborCenter = new IVec(0,0,0);
int neighborCount=0;
for(int i=0; i < agents.size(); i++){
if(agents.get(i) instanceof Cell){
Cell c = (Cell)agents.get(i);
if(c != this){
// push if closer than two radii
IVec dif = c.pos().dif(pos());
c.push(dif);
}
// count neighbors and calculate their center
neighborCount++;
}
}
}
}
if(neighborCount >= 1){ // push from center of neighbors
neighborCenter.div(neighborCount);
IVec dif = pos().dif(neighborCenter).len(100); // constant force
this.push(dif);
}
this.push(nml().len(20)); // pressure toward normal
}

void update(){
if(IG.time() < growthDuration){
if(time() > 0 && time()%divisionInterval==0){
if(active){ // divide when active flag is on
divide();
}
}
if(time()%growthInterval==0){
grow();
}
if(active){ clr(1.0,0.5,0); } // update color by being active
else{ clr(mesh.clr()); }
}
vertex.nml(nml()); // update mesh vertex normal
}

void grow(){ // growing cell size
}
}

void divide(){ // cell division
Cell child = createChild(IRand.dir());
}
// making a triangle loop
Cell child = createChild(IRand.dir());
}
// making a tetrahedron enclosure
Cell child = createChild(IRand.dir());
CellFace f = faces.get(0);
IVec center = IVec.center(child.pos(),f.cell1.pos(),f.cell2.pos(),f.cell3.pos());
if(f.center().dif(center).dot(f.nml())<0){ f.flipNml(); } // adjust normal to be outward
CellFace f1 = new CellFace(f.cell1, f.cell2, child);
CellFace f2 = new CellFace(f.cell2, f.cell3, child);
CellFace f3 = new CellFace(f.cell3, f.cell1, child);
if(f1.center().dif(center).dot(f1.nml())<0){ f1.flipNml(); }
if(f2.center().dif(center).dot(f2.nml())<0){ f2.flipNml(); }
if(f3.center().dif(center).dot(f3.nml())<0){ f3.flipNml(); }
}
// pick two links and insert two faces between
if(linkIdx==null){ // no valid edge to split. skip division.
active = false;
return;
}
}
IVec dir = new IVec();
}
dir.projectToPlane(nml());
Cell child = createChild(dir);
}
}

IVec nml(){ // calc vertex normal from face normal
if(faces.size()==0){ return new IVec(0,0,1); }
IVec n = new IVec(0,0,0);
for(int i=0; i < faces.size(); i++){ n.add(faces.get(i).nml()); }
return n.unit();
}

Cell createChild(IVec dir){
radius *= 0.5; //make both cell size half
Cell child = new Cell(pos().cp(dir), radius);
pos().sub(dir); // move to the other way
if(IRand.pct(90)){ active=false; }//10% stay active
child.active = true; //activate child always
return child;
}

if(links.size() <= 3) return; // no need to sort
for(int i=0; i < links.size(); i++){
}
}

int[] findDividingLinks(){ // find splittable two links and return the index numbers

return null;
}
for(int i=0; i < links.size(); i++){ // check all pairs
return new int[]{ i+idx, i+idx+linkInterval }; // index can be larger than links.size()
}
}
return null;
}

ArrayList< CellFace > removingFaces = new ArrayList< CellFace >();
for(int i=0; i < faces.size(); i++){
for(int j=0; j < linksToChild.size()-1; j++){
break;
}
}
}
for(int i=1; i < linksToChild.size()-1; i++){
}
for(int i=0; i < removingFaces.size(); i++){
removingFaces.get(i).del(); // replace face by deleting and recreating
Cell[] cells = removingFaces.get(i).oppositeCell(this);
new CellFace(child, cells[0], cells[1]);
}
new CellFace(this, cell1, child);
new CellFace(this, child, cell2);
}
}

class CellLink extends IAgent{ // a link to connect 2 cells with spring force
double maxForce = 100;
Cell cell1, cell2;
ArrayList< CellFace > faces;

cell1 = c1; cell2 = c2;
faces = new ArrayList< CellFace >();
}

void interact(ArrayList< IDynamics > agents){
// spring force
IVec dif = cell1.pos().dif(cell2.pos());
if(force > maxForce){ force=maxForce; } // avoid putting too much force
else if(force < -maxForce){ force=-maxForce; }
dif.len(force);
cell1.pull(dif);
cell2.push(dif);
}

boolean contains(Cell c){ //check if link contains the cell
if(c==cell1 || c==cell2) return true;
return false;
}

void del(){
super.del(); // stop agent
}

CellFace oppositeFace(CellFace f){ // find other face on the link
if(faces.size()!=2) return null;
if(faces.get(0)==f) return faces.get(1);
return faces.get(0);
}

Cell oppositeCell(Cell c){ // find other cell on the link
if(cell1==c) return cell2;
if(cell2==c) return cell1;
IG.err("Link does not contain the input cell");
return null;
}

IVec oppositeDir(Cell c){ //calculate a vector to the other cell
return oppositeCell(c).pos().dif(c.pos());
}
}

class CellFace{ // a triangle grouping 3 cells and 3 links
Cell cell1, cell2, cell3;
IFace face;

CellFace(Cell c1, Cell c2, Cell c3){
// keep order of cells in right hand screw for consistent normal
cell1 = c2; cell2 = c1; cell3 = c3;
}
else{
cell1 = c1; cell2 = c2; cell3 = c3;
}
face = new IFace(cell1.vertex,cell2.vertex,cell3.vertex);
}

IVec center(){ // calc center of triangle
return IVec.center(cell1.pos(), cell2.pos(), cell3.pos());
}

IVec nml(){ // normal vector
return cell1.pos().nml(cell2.pos(), cell3.pos()).unit();
}

void flipNml(){ // flip normal
Cell tmp = cell2;
cell2 = cell3;
cell3 = tmp;
}

// true if ordr of c1-c2 follows the order of cell1-cell2-cell3, otherwise false
boolean cellOrder(Cell c1, Cell c2){
if(cell1==c1 && cell2==c2 || cell2==c1 && cell3==c2 || cell3==c1 && cell1==c2 ){ return true; }
return false;
}

}

for(int i=0; i < c1.links.size() ; i++){
}
return null;
}

Cell[] oppositeCell(Cell c){ // find 2 cells opposite of a cell in a triangle
if(cell1==c){ return  new Cell[]{ cell2, cell3 }; }
if(cell2==c){ return  new Cell[]{ cell3, cell1 };  }
if(cell3==c){ return  new Cell[]{ cell1, cell2 };  }
return null;
}

return null;
}

void del(){
cell1.faces.remove(this); // unregister this from cells and links
cell2.faces.remove(this);
cell3.faces.remove(this);
mesh.deleteFace(face);
}
}
```

### Parent/Child Activation Behavior Control

The code below introduces two parameters in Cell class; stayActive and activateChild. They control whether a parent keeps activated and whether a child cell is activated when it's created in cell division. Now it's possible to control the activation control type to be cascading, branching or others outside of createChild method by setting the parameters in a cell.

```import igeo.*;
import processing.opengl.*;

IMesh mesh; // global variable of mesh geometry

void setup(){
size(480,360,IG.GL);
IConfig.syncDrawAndDynamics=true; //not to crash when some geometry is deleted while drawing
mesh = new IMesh().clr(0.7,0,0);
Cell cell = new Cell(new IVec(0,0,0), 10);
cell.stayActive = true;
cell.activateChild = true;
cell.active = true;
}

class Cell extends IParticle{
int growthDuration = 1800; //duration of growth and division
int growthInterval = 10;
int divisionInterval = 20;
double growthSpeed=0.1;

ArrayList< CellFace > faces; //store faces
boolean active = false;
boolean stayActive = false; //activation control variables
boolean activateChild = false;
IVertex vertex;

super(pos, new IVec(0,0,0));
faces = new ArrayList< CellFace >();
vertex = new IVertex();
vertex.pos = pos; // set reference of particle pos directly
fric(0.2);
}

void interact(ArrayList< IDynamics > agents){
IVec neighborCenter = new IVec(0,0,0);
int neighborCount=0;
for(int i=0; i < agents.size(); i++){
if(agents.get(i) instanceof Cell){
Cell c = (Cell)agents.get(i);
if(c != this){
// push if closer than two radii
IVec dif = c.pos().dif(pos());
c.push(dif);
}
// count neighbors and calculate their center
neighborCount++;
}
}
}
}
if(neighborCount >= 1){ // push from center of neighbors
neighborCenter.div(neighborCount);
IVec dif = pos().dif(neighborCenter).len(50); // constant force
this.push(dif);
}
this.push(nml().len(100)); // pressure toward normal
}

void update(){
if(IG.time() < growthDuration){
if(time() > 0 && time()%divisionInterval==0){
if(active){ // divide when active flag is on
divide();
}
}
if(time()%growthInterval==0){
grow();
}
if(active){ clr(1.0,0.5,0); } // update color by being active
else{ clr(mesh.clr()); }
}
vertex.nml(nml()); // update mesh vertex normal
}

void grow(){ // growing cell size
}
}

void divide(){ // cell division
Cell child = createChild(IRand.dir());
}
// making a triangle loop
Cell child = createChild(IRand.dir());
}
// making a tetrahedron enclosure
Cell child = createChild(IRand.dir());
CellFace f = faces.get(0);
IVec center = IVec.center(child.pos(),f.cell1.pos(),f.cell2.pos(),f.cell3.pos());
if(f.center().dif(center).dot(f.nml())<0){ f.flipNml(); } // adjust normal to be outward
CellFace f1 = new CellFace(f.cell1, f.cell2, child);
CellFace f2 = new CellFace(f.cell2, f.cell3, child);
CellFace f3 = new CellFace(f.cell3, f.cell1, child);
if(f1.center().dif(center).dot(f1.nml())<0){ f1.flipNml(); }
if(f2.center().dif(center).dot(f2.nml())<0){ f2.flipNml(); }
if(f3.center().dif(center).dot(f3.nml())<0){ f3.flipNml(); }
}
// pick two links and insert two faces between
if(linkIdx==null){ // no valid edge to split. skip division.
active = false;
return;
}
}
IVec dir = new IVec();
}
dir.projectToPlane(nml());
Cell child = createChild(dir);
}
}

IVec nml(){ // calc vertex normal from face normal
if(faces.size()==0){ return new IVec(0,0,1); }
IVec n = new IVec(0,0,0);
for(int i=0; i < faces.size(); i++){ n.add(faces.get(i).nml()); }
return n.unit();
}

Cell createChild(IVec dir){
radius *= 0.5; //make both cell size half
Cell child = new Cell(pos().cp(dir), radius);
pos().sub(dir); // move to the other way
if(stayActive){ active=true; } //activation control
else{ active=false; }
if(activateChild){ child.active=true; }
else{ child.active=false; }
child.stayActive = stayActive; //activation behavior parameters are passed to child
child.activateChild = activateChild;
return child;
}

if(links.size() <= 3) return; // no need to sort
for(int i=0; i < links.size(); i++){
}
}

int[] findDividingLinks(){ // find splittable two links and return the index numbers

return null;
}
for(int i=0; i < links.size(); i++){ // check all pairs
return new int[]{ i+idx, i+idx+linkInterval }; // index can be larger than links.size()
}
}
return null;
}

ArrayList< CellFace > removingFaces = new ArrayList< CellFace >();
for(int i=0; i < faces.size(); i++){
for(int j=0; j < linksToChild.size()-1; j++){
break;
}
}
}
for(int i=1; i < linksToChild.size()-1; i++){
}
for(int i=0; i < removingFaces.size(); i++){
removingFaces.get(i).del(); // replace face by deleting and recreating
Cell[] cells = removingFaces.get(i).oppositeCell(this);
new CellFace(child, cells[0], cells[1]);
}
new CellFace(this, cell1, child);
new CellFace(this, child, cell2);
}
}

class CellLink extends IAgent{ // a link to connect 2 cells with spring force
double maxForce = 100;
Cell cell1, cell2;
ArrayList< CellFace > faces;

cell1 = c1; cell2 = c2;
faces = new ArrayList< CellFace >();
}

void interact(ArrayList< IDynamics > agents){
// spring force
IVec dif = cell1.pos().dif(cell2.pos());
if(force > maxForce){ force=maxForce; } // avoid putting too much force
else if(force < -maxForce){ force=-maxForce; }
dif.len(force);
cell1.pull(dif);
cell2.push(dif);
}

boolean contains(Cell c){ //check if link contains the cell
if(c==cell1 || c==cell2) return true;
return false;
}

void del(){
super.del(); // stop agent
}

CellFace oppositeFace(CellFace f){ // find other face on the link
if(faces.size()!=2) return null;
if(faces.get(0)==f) return faces.get(1);
return faces.get(0);
}

Cell oppositeCell(Cell c){ // find other cell on the link
if(cell1==c) return cell2;
if(cell2==c) return cell1;
IG.err("Link does not contain the input cell");
return null;
}

IVec oppositeDir(Cell c){ //calculate a vector to the other cell
return oppositeCell(c).pos().dif(c.pos());
}
}

class CellFace{ // a triangle grouping 3 cells and 3 links
Cell cell1, cell2, cell3;
IFace face;

CellFace(Cell c1, Cell c2, Cell c3){
// keep order of cells in right hand screw for consistent normal
cell1 = c2; cell2 = c1; cell3 = c3;
}
else{
cell1 = c1; cell2 = c2; cell3 = c3;
}
face = new IFace(cell1.vertex,cell2.vertex,cell3.vertex);
}

IVec center(){ // calc center of triangle
return IVec.center(cell1.pos(), cell2.pos(), cell3.pos());
}

IVec nml(){ // normal vector
return cell1.pos().nml(cell2.pos(), cell3.pos()).unit();
}

void flipNml(){ // flip normal
Cell tmp = cell2;
cell2 = cell3;
cell3 = tmp;
}

// true if ordr of c1-c2 follows the order of cell1-cell2-cell3, otherwise false
boolean cellOrder(Cell c1, Cell c2){
if(cell1==c1 && cell2==c2 || cell2==c1 && cell3==c2 || cell3==c1 && cell1==c2 ){ return true; }
return false;
}

}

for(int i=0; i < c1.links.size() ; i++){
}
return null;
}

Cell[] oppositeCell(Cell c){ // find 2 cells opposite of a cell in a triangle
if(cell1==c){ return  new Cell[]{ cell2, cell3 }; }
if(cell2==c){ return  new Cell[]{ cell3, cell1 };  }
if(cell3==c){ return  new Cell[]{ cell1, cell2 };  }
return null;
}

return null;
}

void del(){
cell1.faces.remove(this); // unregister this from cells and links
cell2.faces.remove(this);
cell3.faces.remove(this);
mesh.deleteFace(face);
}
}
```

### Cell Division and Growth from Imported Mesh

The code below import a file which contains several meshes and build a network of cells, cell links and cell faces for each mesh. The method cellsFromMesh generates cells. You can specify which cells to be fixed adn which cells to be activated by locating points in spefici layers. The codes of Cell, CellLink, CellFace class add some fields and methods to deal with exsisting mesh and both open and closed mesh.

The sample mesh file used in the example is below.

```import igeo.*;
import processing.opengl.*;

double sizeRatio; // to adjustment forces based on size

void setup(){
size(480,360,IG.GL);
IConfig.syncDrawAndDynamics=true; //not to crash when some geometry is deleted while drawing
IG.open("mesh2.3dm");

IPoint[] fixedPts = IG.layer("fix points").points(); // points to speficy fixed cells
IPoint[] activePts = IG.layer("active points").points(); // points to speficy active cells
for(int i=0; i < IG.meshNum(); i++){
IMesh m = IG.mesh(i);
m.triangulate(); // quad mesh not accepted

Cell[] cells = cellsFromMesh(m, fixedPts);

for(int j=0; j < cells.length; j++){
for(int k=0; activePts!=null && k < activePts.length; k++){
if(cells[j].pos().eq(activePts[k])){ // activate cell if active point is located
cells[j].stayActiveProbability = 80;
cells[j].activateChild=true;
cells[j].active=true;
}
}
}
}
IG.del(activePts);
}

Cell[] cellsFromMesh(IMesh mesh, IPoint[] fixedPoints){ // create cells, cell links, cell faces from a mesh
Cell[] cells = new Cell[mesh.vertexNum()];
averageLength=0;
for(int i=0; i < mesh.edgeNum(); i++){
averageLength += mesh.edge(i).vertex(0).dist(mesh.edge(i).vertex(1));
}
averageLength /= mesh.edgeNum();
sizeRatio = averageLength/10;

for(int i=0; i < mesh.vertexNum(); i++){
for(int j=0; j < mesh.vertex(i).edgeNum(); j++){
}
cells[i] = new Cell(mesh.vertex(i), radius, mesh);

boolean done=false;
for(int j=0; j < fixedPoints.length&&!done; j++){
if(cells[i].pos.eq(fixedPoints[j])){
cells[i].fix();
}
}
}

for(int i=0; i < mesh.edgeNum(); i++){
int idx1 = mesh.mesh.vertices.indexOf(mesh.edge(i).vertex(0));
int idx2 = mesh.mesh.vertices.indexOf(mesh.edge(i).vertex(1));
}

CellFace[] faces = new CellFace[mesh.faceNum()];
for(int i=0; i < mesh.faceNum(); i++){
int idx1 = mesh.mesh.vertices.indexOf(mesh.face(i).vertex(0));
int idx2 = mesh.mesh.vertices.indexOf(mesh.face(i).vertex(1));
int idx3 = mesh.mesh.vertices.indexOf(mesh.face(i).vertex(2));
faces[i] = new CellFace(cells[idx1],cells[idx2],cells[idx3],mesh.face(i),mesh);
}

return cells;
}

class Cell extends IParticle{
int growthDuration = 2000; //duration of growth and division
int growthInterval = 10;
int divisionInterval = 100;
double maxRadius = averageLength/2; // no bigger than half of average edge length
double growthSpeed=0.5*sizeRatio;

ArrayList< CellFace > faces; //store faces
boolean active = false;
double stayActiveProbability = 0; //activation control probability
boolean activateChild = false;
IVertex vertex;

IMesh mesh;

Cell(IVec pos, double rad, IMesh m){
super(pos, new IVec(0,0,0));
faces = new ArrayList< CellFace >();
vertex = new IVertex();
vertex.pos = pos; // set reference of particle pos directly
mesh = m;
fric(0.2);
}

Cell(IVertex vtx, double rad, IMesh m){
super(vtx.pos, new IVec(0,0,0));
faces = new ArrayList< CellFace >();
vertex = vtx;
mesh = m;
fric(0.2);
}

void interact(ArrayList< IDynamics > agents){
IVec neighborCenter = new IVec(0,0,0);
int neighborCount=0;
for(int i=0; i < agents.size(); i++){
if(agents.get(i) instanceof Cell){
Cell c = (Cell)agents.get(i);
if(c != this){
// push if closer than two radii
IVec dif = c.pos().dif(pos());
c.push(dif);
}
// count neighbors and calculate their center
neighborCount++;
}
}
}
}
if(neighborCount >= 1){ // push from center of neighbors
neighborCenter.div(neighborCount);
IVec dif = pos().dif(neighborCenter);
if(dif.len()>0){
this.push(dif.len(50*sizeRatio)); // constant force
}
}
this.push(nml().len(50*sizeRatio)); // pressure toward normal
}

void update(){
if(IG.time() < growthDuration){
if(time() > 0 && time()%divisionInterval==0){
if(active){ // divide when active flag is on
divide();
}
}
if(time()%growthInterval==0){
grow();
}
if(active){ clr(1.0,0.5,0); } // update color by being active
else{ clr(mesh.clr()); }
}
vertex.nml(nml()); // update mesh vertex normal
}

void grow(){ // growing cell size
}
}

void divide(){ // cell division
// pick two links and insert two faces between
if(linkIdx==null){ // no valid edge to split. skip division.
active = false;
return;
}

}
}
else{
}
}
IVec dir = new IVec();
}
dir.projectToPlane(nml());
Cell child = createChild(dir);
}
}

IVec nml(){ // calc vertex normal from face normal
if(faces.size()==0){ return new IVec(0,0,1); }
IVec n = new IVec(0,0,0);
for(int i=0; i < faces.size(); i++){ n.add(faces.get(i).nml()); }
return n.unit();
}

Cell createChild(IVec dir){
radius *= 0.5; //make both cell size half
Cell child = new Cell(pos().cp(dir), radius, mesh);
pos().sub(dir); // move to the other way
if(IRand.pct(stayActiveProbability)){ active=true; } //activation control
else{ active=false; }
if(activateChild){ child.active=true; }
else{ child.active=false; }
child.stayActiveProbability = stayActiveProbability; //activation behavior parameters are passed to child
child.activateChild = activateChild;
return child;
}

if(links.size() <= 2) return; // no need to sort
}
}
}

for(int i=0; i < links.size(); i++){
if(currentFace==null) break;
}
}
}

int[] findDividingLinks(){ // find splittable two links and return the index numbers

return null;
}

if(links.get(0).faces.size()==1 && IRand.pct(10)){ // if open edge, divide at one edge
for(int i=0; i < links.size(); i++){ // check all pairs
return new int[]{ (i+idx)%links.size() }; // index can be larger than links.size()
}
}
}
return null;
}

for(int i=0; i < links.size(); i++){ // check all pairs
return new int[]{ i+idx, i+idx+linkInterval }; // index can be larger than links.size()
}
}
return null;
}

ArrayList< CellFace > removingFaces = new ArrayList< CellFace >();
for(int j=0; j < linksToChild.size()-1; j++){
for(int i=0; i < faces.size(); i++){
break;
}
}
}

if(insertAtEdge){
for(int i=1; i < linksToChild.size(); i++){
}
for(int i=0; i < removingFaces.size(); i++){
removingFaces.get(i).del(); // replace face by deleting and recreating
Cell[] cells = removingFaces.get(i).oppositeCell(this);
new CellFace(child, cells[0], cells[1], mesh);
}
new CellFace(this, cell1, child, mesh);
}
else{
for(int i=1; i < linksToChild.size()-1; i++){
}

for(int i=0; i < removingFaces.size(); i++){
removingFaces.get(i).del(); // replace face by deleting and recreating
Cell[] cells = removingFaces.get(i).oppositeCell(this);
new CellFace(child, cells[0], cells[1], mesh);
}
new CellFace(this, cell1, child, mesh);
new CellFace(this, child, cell2, mesh);
}
}

for(int i=0; i < links.size() ; i++){
}
return null;
}
}

class CellLink extends IAgent{ // a link to connect 2 cells with spring force
double maxForce = 100;
Cell cell1, cell2;
ArrayList< CellFace > faces;

cell1 = c1; cell2 = c2;
faces = new ArrayList< CellFace >();
}

void interact(ArrayList< IDynamics > agents){
// spring force
IVec dif = cell1.pos().dif(cell2.pos());
if(force > maxForce){ force=maxForce; } // avoid putting too much force
else if(force < -maxForce){ force=-maxForce; }
dif.len(force);
cell1.pull(dif);
cell2.push(dif);
}

boolean contains(Cell c){ //check if link contains the cell
if(c==cell1 || c==cell2) return true;
return false;
}

void del(){
super.del(); // stop agent
}

CellFace oppositeFace(CellFace f){ // find other face on the link
if(faces.size()!=2) return null;
if(faces.get(0)==f) return faces.get(1);
return faces.get(0);
}

Cell oppositeCell(Cell c){ // find other cell on the link
if(cell1==c) return cell2;
if(cell2==c) return cell1;
IG.err("Link does not contain the input cell");
return null;
}

IVec oppositeDir(Cell c){ //calculate a vector to the other cell
return oppositeCell(c).pos().dif(c.pos());
}

return null;
}
}

class CellFace{ // a triangle grouping 3 cells and 3 links
Cell cell1, cell2, cell3;
IFace face;
IMesh mesh;

CellFace(Cell c1, Cell c2, Cell c3, IMesh m){
mesh = m;
// keep order of cells in right hand screw for consistent normal
cell1 = c2; cell2 = c1; cell3 = c3;
}
else{
cell1 = c1; cell2 = c2; cell3 = c3;
}
face = new IFace(cell1.vertex,cell2.vertex,cell3.vertex);
}

mesh = m;
Cell c1 = l3.sharedCell(l1);
Cell c2 = l1.sharedCell(l2);
Cell c3 = l2.sharedCell(l3);

// keep order of cells in right hand screw for consistent normal
if(l1.faces.size()==1 && l1.faces.get(0).cellOrder(c1,c2)){
cell1 = c2; cell2 = c1; cell3 = c3;
}
else{
cell1 = c1; cell2 = c2; cell3 = c3;
}
face = new IFace(cell1.vertex,cell2.vertex,cell3.vertex);
}

CellFace(Cell c1, Cell c2, Cell c3, IFace f, IMesh m){
mesh = m;
// keep order of cells in right hand screw for consistent normal
cell1 = c2; cell2 = c1; cell3 = c3;
}
else{
cell1 = c1; cell2 = c2; cell3 = c3;
}
face = f;
}

IVec center(){ // calc center of triangle
return IVec.center(cell1.pos(), cell2.pos(), cell3.pos());
}

IVec nml(){ // normal vector
return cell1.pos().nml(cell2.pos(), cell3.pos()).unit();
}

void flipNml(){ // flip normal
Cell tmp = cell2;
cell2 = cell3;
cell3 = tmp;
}

// true if ordr of c1-c2 follows the order of cell1-cell2-cell3, otherwise false
boolean cellOrder(Cell c1, Cell c2){
if(cell1==c1 && cell2==c2 || cell2==c1 && cell3==c2 || cell3==c1 && cell1==c2 ){ return true; }
return false;
}

}

for(int i=0; i < c1.links.size() ; i++){
}
return null;
}

Cell[] oppositeCell(Cell c){ // find 2 cells opposite of a cell in a triangle
if(cell1==c){ return  new Cell[]{ cell2, cell3 }; }
if(cell2==c){ return  new Cell[]{ cell3, cell1 };  }
if(cell3==c){ return  new Cell[]{ cell1, cell2 };  }
return null;
}

return null;
}

void del(){
cell1.faces.remove(this); // unregister this from cells and links
cell2.faces.remove(this);
cell3.faces.remove(this);