Your description and ideas seem much more advanced than your code.

While it is often a very good idea to use classes/objects I get a strong 3D matrix feel from your schematic. Perhaps 3D arrays could be useful here (see example below) to keep track of existing cubes, colors etc. Perhaps you need multiple kind of systems/objects/objectManagers.

With regard to libraries I'm thinking Hemesh and Toxiclibs for their strengths and convenient classes/methods in the areas of math, mesh and interconnectivity information.

The example below has a similar look to your schematic. Click the mouse at different places on the x-axis to change the density of the cubes.

Code Example

1. import processing.opengl.*;
2.  
3. int dim = 10;
4. int halfdim = dim/2;
5. boolean[][][] matrix = new boolean[dim][dim][dim];
6. color[][][] colorMatrix = new color[dim][dim][dim];
7. int boxSize = 30;
8.  
9. void setup() {
10.   size(500, 500, OPENGL);
11.   setMatrix(0.5);
12.   smooth();
13. }
14.  
15. void draw() {
16.   background(255);
17.   translate(width/2, height/2, -400);
18.   rotateX(radians(frameCount*0.5));
19.   rotateY(radians(frameCount*0.2));
20.   rotateZ(radians(frameCount*0.7));
21.   for (int z=0; z<dim; z++) {
22.     for (int y=0; y<dim; y++) {
23.       for (int x=0; x<dim; x++) {
24.         if (matrix[x][y][z]) {
25.           pushMatrix();
26.           translate((x-halfdim)*boxSize, (y-halfdim)*boxSize, (z-halfdim)*boxSize);
27.           fill(colorMatrix[x][y][z]);
28.           box(boxSize);
29.           popMatrix();
30.         }
31.       }
32.     }
33.   }
34. }
35.  
36. void mousePressed() {
37.   setMatrix((float) mouseX/width);
38. }
39.  
40. void setMatrix(float range) {
41.   for (int z=0; z<dim; z++) {
42.     for (int y=0; y<dim; y++) {
43.       for (int x=0; x<dim; x++) {
44.         matrix[x][y][z] = random(1)>range?false:true;
45.         colorMatrix[x][y][z] = color(random(255), random(255), random(255));
46.       }
47.     }
48.   }
49. }

Thank you so much for the reply. It is just the beginning of the project. I am looking at toxiclibs at the moment. I am targeting to get the array and the growing working first. Then the environment and the decay. 

just a short suggestion. I think you should seperate the logic and the drawing into seperate classes/couple of classes. so you have a logic manager which stores the cube matrix. in there are the cubes, each cube described by a 3D Vector with a drawing independent coordinate system, a color, etc.. so you have a cube at 0,0,0 at next to it a cube with 1,0,0. with this system you can lookup more easily if a cube ist next to another.

There are many routes and possible solutions. Just start sketching and develop incrementally.

With regard to your question, yes you can use 3D arrays to store Cube's as well. For me the reason to use a 3D array here is that it fits the positioning of the objects and it will (once you need it) make it easier to optimize the neighbour searching through simple +/- 1 operations. The current class doesn't do any neighbour stuff, but showcases the concept of Cube's in a 3D array.

Code Example

1. import processing.opengl.*;
2.  
3. int dim = 10;
4. int halfdim = dim/2;
5. int maxSize = 50;
6. Cube[][][] matrix = new Cube[dim][dim][dim];
7.  
8. void setup() {
9.   size(800, 800, OPENGL);
10.   for (int z=0; z<dim; z++) {
11.     for (int y=0; y<dim; y++) {
12.       for (int x=0; x<dim; x++) {
13.         matrix[x][y][z] = new Cube(x, y, z, maxSize);
14.       }
15.     }
16.   }
17.   smooth();
18.   stroke(255);
19. }
20.  
21. void draw() {
22.   background(255);
23.   translate(width/2, height/2, sin(frameCount*0.01)*200);
24.   rotateX(radians(frameCount*0.3));
25.   rotateY(radians(frameCount*0.2));
26.   rotateZ(radians(frameCount*0.1));
27.   for (int z=0; z<dim; z++) {
28.     for (int y=0; y<dim; y++) {
29.       for (int x=0; x<dim; x++) {
30.         Cube currentCube = matrix[x][y][z];
31.         currentCube.update();
32.         currentCube.display();
33.       }
34.     }
35.   }
36. }
37.  
38. class Cube extends PVector {
39.   float cubeSize, maxCubeSize;
40.   color c;
41.  
42.   Cube(float x, float y, float z, float maxSize) {
43.     super((x-halfdim)*maxSize, (y-halfdim)*maxSize, (z-halfdim)*maxSize);
44.     cubeSize = maxCubeSize = maxSize;
45.     c = color(random(255), random(255), random(255));
46.   }
47.  
48.   void update() {
49.     cubeSize = constrain((sin(z+frameCount*0.01)+cos(x+y+frameCount*0.02)),0,1)*maxCubeSize;
50.   }
51.  
52.   void display() {
53.     pushMatrix();
54.     translate(x, y, z);
55.     fill(c);
56.     box(cubeSize);
57.     popMatrix();
58.   }
59. }