KickJava   Java API By Example, From Geeks To Geeks.

Java > Open Source Codes > prefuse > util > force > NBodyForce


1 package prefuse.util.force;
2
3 import java.util.ArrayList JavaDoc;
4 import java.util.Arrays JavaDoc;
5 import java.util.Iterator JavaDoc;
6 import java.util.Random JavaDoc;
7
8 /**
9  * <p>Force function which computes an n-body force such as gravity,
10  * anti-gravity, or the results of electric charges. This function implements
11  * the the Barnes-Hut algorithm for efficient n-body force simulations,
12  * using a quad-tree with aggregated mass values to compute the n-body
13  * force in O(N log N) time, where N is the number of ForceItems.</p>
14  *
15  * <p>The algorithm used is that of J. Barnes and P. Hut, in their research
16  * paper <i>A Hierarchical O(n log n) force calculation algorithm</i>, Nature,
17  * v.324, December 1986. For more details on the algorithm, see one of
18  * the following links --
19  * <ul>
20  * <li><a HREF="http://www.cs.berkeley.edu/~demmel/cs267/lecture26/lecture26.html">James Demmel's UC Berkeley lecture notes</a>
21  * <li><a HREF="http://www.physics.gmu.edu/~large/lr_forces/desc/bh/bhdesc.html">Description of the Barnes-Hut algorithm</a>
22  * <li><a HREF="http://www.ifa.hawaii.edu/~barnes/treecode/treeguide.html">Joshua Barnes' recent implementation</a>
23  * </ul></p>
24  *
25  * @author <a HREF="http://jheer.org">jeffrey heer</a>
26  */
27 public class NBodyForce extends AbstractForce {
28
29     /*
30      * The indexing scheme for quadtree child nodes goes row by row.
31      * 0 | 1 0 -> top left, 1 -> top right
32      * -------
33      * 2 | 3 2 -> bottom left, 3 -> bottom right
34      */
35
36     private static String JavaDoc[] pnames = new String JavaDoc[] { "GravitationalConstant",
37             "Distance", "BarnesHutTheta" };
38     
39     public static final float DEFAULT_GRAV_CONSTANT = -1.0f;
40     public static final float DEFAULT_MIN_GRAV_CONSTANT = -10f;
41     public static final float DEFAULT_MAX_GRAV_CONSTANT = 10f;
42     
43     public static final float DEFAULT_DISTANCE = -1f;
44     public static final float DEFAULT_MIN_DISTANCE = -1f;
45     public static final float DEFAULT_MAX_DISTANCE = 500f;
46     
47     public static final float DEFAULT_THETA = 0.9f;
48     public static final float DEFAULT_MIN_THETA = 0.0f;
49     public static final float DEFAULT_MAX_THETA = 1.0f;
50     
51     public static final int GRAVITATIONAL_CONST = 0;
52     public static final int MIN_DISTANCE = 1;
53     public static final int BARNES_HUT_THETA = 2;
54     
55     private float xMin, xMax, yMin, yMax;
56     private QuadTreeNodeFactory factory = new QuadTreeNodeFactory();
57     private QuadTreeNode root;
58     
59     private Random JavaDoc rand = new Random JavaDoc(12345678L); // deterministic randomness
60
61     /**
62      * Create a new NBodyForce with default parameters.
63      */
64     public NBodyForce() {
65         this(DEFAULT_GRAV_CONSTANT, DEFAULT_DISTANCE, DEFAULT_THETA);
66     }
67     
68     /**
69      * Create a new NBodyForce.
70      * @param gravConstant the gravitational constant to use. Nodes will
71      * attract each other if this value is positive, and will repel each
72      * other if it is negative.
73      * @param minDistance the distance within which two particles will
74      * interact. If -1, the value is treated as infinite.
75      * @param theta the Barnes-Hut parameter theta, which controls when
76      * an aggregated mass is used rather than drilling down to individual
77      * item mass values.
78      */
79     public NBodyForce(float gravConstant, float minDistance, float theta) {
80         params = new float[] { gravConstant, minDistance, theta };
81         minValues = new float[] { DEFAULT_MIN_GRAV_CONSTANT,
82             DEFAULT_MIN_DISTANCE, DEFAULT_MIN_THETA };
83         maxValues = new float[] { DEFAULT_MAX_GRAV_CONSTANT,
84             DEFAULT_MAX_DISTANCE, DEFAULT_MAX_THETA };
85         root = factory.getQuadTreeNode();
86     }
87
88     /**
89      * Returns true.
90      * @see prefuse.util.force.Force#isItemForce()
91      */
92     public boolean isItemForce() {
93         return true;
94     }
95     
96     /**
97      * @see prefuse.util.force.AbstractForce#getParameterNames()
98      */
99     protected String JavaDoc[] getParameterNames() {
100         return pnames;
101     }
102     
103     /**
104      * Set the bounds of the region for which to compute the n-body simulation
105      * @param xMin the minimum x-coordinate
106      * @param yMin the minimum y-coordinate
107      * @param xMax the maximum x-coordinate
108      * @param yMax the maximum y-coordinate
109      */
110     private void setBounds(float xMin, float yMin, float xMax, float yMax) {
111         this.xMin = xMin;
112         this.yMin = yMin;
113         this.xMax = xMax;
114         this.yMax = yMax;
115     }
116
117     /**
118      * Clears the quadtree of all entries.
119      */
120     public void clear() {
121         clearHelper(root);
122         root = factory.getQuadTreeNode();
123     }
124     
125     private void clearHelper(QuadTreeNode n) {
126         for ( int i = 0; i < n.children.length; i++ ) {
127             if ( n.children[i] != null )
128                 clearHelper(n.children[i]);
129         }
130         factory.reclaim(n);
131     }
132
133     /**
134      * Initialize the simulation with the provided enclosing simulation. After
135      * this call has been made, the simulation can be queried for the
136      * n-body force acting on a given item.
137      * @param fsim the enclosing ForceSimulator
138      */
139     public void init(ForceSimulator fsim) {
140         clear(); // clear internal state
141
142         // compute and squarify bounds of quadtree
143 float x1 = Float.MAX_VALUE, y1 = Float.MAX_VALUE;
144         float x2 = Float.MIN_VALUE, y2 = Float.MIN_VALUE;
145         Iterator JavaDoc itemIter = fsim.getItems();
146         while ( itemIter.hasNext() ) {
147             ForceItem item = (ForceItem)itemIter.next();
148             float x = item.location[0];
149             float y = item.location[1];
150             if ( x < x1 ) x1 = x;
151             if ( y < y1 ) y1 = y;
152             if ( x > x2 ) x2 = x;
153             if ( y > y2 ) y2 = y;
154         }
155         float dx = x2-x1, dy = y2-y1;
156         if ( dx > dy ) { y2 = y1 + dx; } else { x2 = x1 + dy; }
157         setBounds(x1,y1,x2,y2);
158         
159         // insert items into quadtree
160 itemIter = fsim.getItems();
161         while ( itemIter.hasNext() ) {
162             ForceItem item = (ForceItem)itemIter.next();
163             insert(item);
164         }
165         
166         // calculate magnitudes and centers of mass
167 calcMass(root);
168     }
169
170     /**
171      * Inserts an item into the quadtree.
172      * @param item the ForceItem to add.
173      * @throws IllegalStateException if the current location of the item is
174      * outside the bounds of the quadtree
175      */
176     public void insert(ForceItem item) {
177         // insert item into the quadtrees
178 try {
179             insert(item, root, xMin, yMin, xMax, yMax);
180         } catch ( StackOverflowError JavaDoc e ) {
181             // TODO: safe to remove?
182 e.printStackTrace();
183         }
184     }
185
186     private void insert(ForceItem p, QuadTreeNode n,
187                         float x1, float y1, float x2, float y2)
188     {
189         // try to insert particle p at node n in the quadtree
190 // by construction, each leaf will contain either 1 or 0 particles
191 if ( n.hasChildren ) {
192             // n contains more than 1 particle
193 insertHelper(p,n,x1,y1,x2,y2);
194         } else if ( n.value != null ) {
195             // n contains 1 particle
196 if ( isSameLocation(n.value, p) ) {
197                 insertHelper(p,n,x1,y1,x2,y2);
198             } else {
199                 ForceItem v = n.value; n.value = null;
200                 insertHelper(v,n,x1,y1,x2,y2);
201                 insertHelper(p,n,x1,y1,x2,y2);
202             }
203         } else {
204             // n is empty, so is a leaf
205 n.value = p;
206         }
207     }
208     
209     private static boolean isSameLocation(ForceItem f1, ForceItem f2) {
210         float dx = Math.abs(f1.location[0]-f2.location[0]);
211         float dy = Math.abs(f1.location[1]-f2.location[1]);
212         return ( dx < 0.01 && dy < 0.01 );
213     }
214     
215     private void insertHelper(ForceItem p, QuadTreeNode n,
216                               float x1, float y1, float x2, float y2)
217     {
218         float x = p.location[0], y = p.location[1];
219         float splitx = (x1+x2)/2;
220         float splity = (y1+y2)/2;
221         int i = (x>=splitx ? 1 : 0) + (y>=splity ? 2 : 0);
222         // create new child node, if necessary
223 if ( n.children[i] == null ) {
224             n.children[i] = factory.getQuadTreeNode();
225             n.hasChildren = true;
226         }
227         // update bounds
228 if ( i==1 || i==3 ) x1 = splitx; else x2 = splitx;
229         if ( i > 1 ) y1 = splity; else y2 = splity;
230         // recurse
231 insert(p,n.children[i],x1,y1,x2,y2);
232     }
233
234     private void calcMass(QuadTreeNode n) {
235         float xcom = 0, ycom = 0;
236         n.mass = 0;
237         if ( n.hasChildren ) {
238             for ( int i=0; i < n.children.length; i++ ) {
239                 if ( n.children[i] != null ) {
240                     calcMass(n.children[i]);
241                     n.mass += n.children[i].mass;
242                     xcom += n.children[i].mass * n.children[i].com[0];
243                     ycom += n.children[i].mass * n.children[i].com[1];
244                 }
245             }
246         }
247         if ( n.value != null ) {
248             n.mass += n.value.mass;
249             xcom += n.value.mass * n.value.location[0];
250             ycom += n.value.mass * n.value.location[1];
251         }
252         n.com[0] = xcom / n.mass;
253         n.com[1] = ycom / n.mass;
254     }
255
256     /**
257      * Calculates the force vector acting on the given item.
258      * @param item the ForceItem for which to compute the force
259      */
260     public void getForce(ForceItem item) {
261         try {
262             forceHelper(item,root,xMin,yMin,xMax,yMax);
263         } catch ( StackOverflowError JavaDoc e ) {
264             // TODO: safe to remove?
265 e.printStackTrace();
266         }
267     }
268     
269     private void forceHelper(ForceItem item, QuadTreeNode n,
270                              float x1, float y1, float x2, float y2)
271     {
272         float dx = n.com[0] - item.location[0];
273         float dy = n.com[1] - item.location[1];
274         float r = (float)Math.sqrt(dx*dx+dy*dy);
275         boolean same = false;
276         if ( r == 0.0f ) {
277             // if items are in the exact same place, add some noise
278 dx = (rand.nextFloat()-0.5f) / 50.0f;
279             dy = (rand.nextFloat()-0.5f) / 50.0f;
280             r = (float)Math.sqrt(dx*dx+dy*dy);
281             same = true;
282         }
283         boolean minDist = params[MIN_DISTANCE]>0f && r>params[MIN_DISTANCE];
284         
285         // the Barnes-Hut approximation criteria is if the ratio of the
286 // size of the quadtree box to the distance between the point and
287 // the box's center of mass is beneath some threshold theta.
288 if ( (!n.hasChildren && n.value != item) ||
289              (!same && (x2-x1)/r < params[BARNES_HUT_THETA]) )
290         {
291             if ( minDist ) return;
292             // either only 1 particle or we meet criteria
293 // for Barnes-Hut approximation, so calc force
294 float v = params[GRAVITATIONAL_CONST]*item.mass*n.mass
295                         / (r*r*r);
296             item.force[0] += v*dx;
297             item.force[1] += v*dy;
298         } else if ( n.hasChildren ) {
299             // recurse for more accurate calculation
300 float splitx = (x1+x2)/2;
301             float splity = (y1+y2)/2;
302             for ( int i=0; i<n.children.length; i++ ) {
303                 if ( n.children[i] != null ) {
304                     forceHelper(item, n.children[i],
305                         (i==1||i==3?splitx:x1), (i>1?splity:y1),
306                         (i==1||i==3?x2:splitx), (i>1?y2:splity));
307                 }
308             }
309             if ( minDist ) return;
310             if ( n.value != null && n.value != item ) {
311                 float v = params[GRAVITATIONAL_CONST]*item.mass*n.value.mass
312                             / (r*r*r);
313                 item.force[0] += v*dx;
314                 item.force[1] += v*dy;
315             }
316         }
317     }
318
319     /**
320      * Represents a node in the quadtree.
321      */
322     public static final class QuadTreeNode {
323         public QuadTreeNode() {
324             com = new float[] {0.0f, 0.0f};
325             children = new QuadTreeNode[4];
326         } //
327 boolean hasChildren = false;
328         float mass; // total mass held by this node
329 float[] com; // center of mass of this node
330 ForceItem value; // ForceItem in this node, null if node has children
331 QuadTreeNode[] children; // children nodes
332 } // end of inner class QuadTreeNode
333
334     /**
335      * Helper class to minimize number of object creations across multiple
336      * uses of the quadtree.
337      */
338     public static final class QuadTreeNodeFactory {
339         private int maxNodes = 50000;
340         private ArrayList JavaDoc nodes = new ArrayList JavaDoc();
341         
342         public QuadTreeNode getQuadTreeNode() {
343             if ( nodes.size() > 0 ) {
344                 return (QuadTreeNode)nodes.remove(nodes.size()-1);
345             } else {
346                 return new QuadTreeNode();
347             }
348         }
349         public void reclaim(QuadTreeNode n) {
350             n.mass = 0;
351             n.com[0] = 0.0f; n.com[1] = 0.0f;
352             n.value = null;
353             n.hasChildren = false;
354             Arrays.fill(n.children, null);
355             if ( nodes.size() < maxNodes )
356                 nodes.add(n);
357         }
358     } // end of inner class QuadTreeNodeFactory
359
360 } // end of class NBodyForce
361
Popular Tags