1 /**
  2  * Constructs a new charge force, with an optional charge constant. The charge
  3  * constant can be negative for repulsion (e.g., particles with electrical
  4  * charge of equal sign), or positive for attraction (e.g., massive particles
  5  * with mutual gravity). The default charge constant is -40.
  6  *
  7  * @class An n-body force, as defined by Coulomb's law or Newton's law of
  8  * gravitation, inversely proportional to the square of the distance between
  9  * particles. Note that the force is independent of the <i>mass</i> of the
 10  * associated particles, and that the particles do not have charges of varying
 11  * magnitude; instead, the attraction or repulsion of all particles is globally
 12  * specified as the charge {@link #constant}.
 13  *
 14  * <p>This particular implementation uses the Barnes-Hut algorithm. For details,
 15  * see <a
 16  * href="http://www.nature.com/nature/journal/v324/n6096/abs/324446a0.html">"A
 17  * hierarchical O(N log N) force-calculation algorithm"</a>, J. Barnes &
 18  * P. Hut, <i>Nature</i> 1986.
 19  *
 20  * @name pv.Force.charge
 21  * @param {number} [k] the charge constant.
 22  */
 23 pv.Force.charge = function(k) {
 24   var min = 2, // minimum distance at which to observe forces
 25       min1 = 1 / min,
 26       max = 500, // maximum distance at which to observe forces
 27       max1 = 1 / max,
 28       theta = .9, // Barnes-Hut theta approximation constant
 29       force = {};
 30 
 31   if (!arguments.length) k = -40; // default charge constant (repulsion)
 32 
 33   /**
 34    * Sets or gets the charge constant. If an argument is specified, it is the
 35    * new charge constant. The charge constant can be negative for repulsion
 36    * (e.g., particles with electrical charge of equal sign), or positive for
 37    * attraction (e.g., massive particles with mutual gravity). The default
 38    * charge constant is -40.
 39    *
 40    * @function
 41    * @name pv.Force.charge.prototype.constant
 42    * @param {number} x the charge constant.
 43    * @returns {pv.Force.charge} this.
 44    */
 45   force.constant = function(x) {
 46     if (arguments.length) {
 47       k = Number(x);
 48       return force;
 49     }
 50     return k;
 51   };
 52 
 53   /**
 54    * Sets or gets the domain; specifies the minimum and maximum domain within
 55    * which charge forces are applied. A minimum distance threshold avoids
 56    * applying forces that are two strong (due to granularity of the simulation's
 57    * numeric integration). A maximum distance threshold improves performance by
 58    * skipping force calculations for particles that are far apart.
 59    *
 60    * <p>The default domain is [2, 500].
 61    *
 62    * @function
 63    * @name pv.Force.charge.prototype.domain
 64    * @param {number} a
 65    * @param {number} b
 66    * @returns {pv.Force.charge} this.
 67    */
 68   force.domain = function(a, b) {
 69     if (arguments.length) {
 70       min = Number(a);
 71       min1 = 1 / min;
 72       max = Number(b);
 73       max1 = 1 / max;
 74       return force;
 75     }
 76     return [min, max];
 77   };
 78 
 79   /**
 80    * Sets or gets the Barnes-Hut approximation factor. The Barnes-Hut
 81    * approximation criterion is the ratio of the size of the quadtree node to
 82    * the distance from the point to the node's center of mass is beneath some
 83    * threshold.
 84    *
 85    * @function
 86    * @name pv.Force.charge.prototype.theta
 87    * @param {number} x the new Barnes-Hut approximation factor.
 88    * @returns {pv.Force.charge} this.
 89    */
 90   force.theta = function(x) {
 91     if (arguments.length) {
 92       theta = Number(x);
 93       return force;
 94     }
 95     return theta;
 96   };
 97 
 98   /**
 99    * @ignore Recursively computes the center of charge for each node in the
100    * quadtree. This is equivalent to the center of mass, assuming that all
101    * particles have unit weight.
102    */
103   function accumulate(n) {
104     var cx = 0, cy = 0;
105     n.cn = 0;
106     function accumulateChild(c) {
107       accumulate(c);
108       n.cn += c.cn;
109       cx += c.cn * c.cx;
110       cy += c.cn * c.cy;
111     }
112     if (!n.leaf) {
113       if (n.c1) accumulateChild(n.c1);
114       if (n.c2) accumulateChild(n.c2);
115       if (n.c3) accumulateChild(n.c3);
116       if (n.c4) accumulateChild(n.c4);
117     }
118     if (n.p) {
119       n.cn += k;
120       cx += k * n.p.x;
121       cy += k * n.p.y;
122     }
123     n.cx = cx / n.cn;
124     n.cy = cy / n.cn;
125   }
126 
127   /**
128    * @ignore Recursively computes forces on the given particle using the given
129    * quadtree node. The Barnes-Hut approximation criterion is the ratio of the
130    * size of the quadtree node to the distance from the point to the node's
131    * center of mass is beneath some threshold.
132    */
133   function forces(n, p, x1, y1, x2, y2) {
134     var dx = n.cx - p.x,
135         dy = n.cy - p.y,
136         dn = 1 / Math.sqrt(dx * dx + dy * dy);
137 
138     /* Barnes-Hut criterion. */
139     if ((n.leaf && (n.p != p)) || ((x2 - x1) * dn < theta)) {
140       if (dn < max1) return;
141       if (dn > min1) dn = min1;
142       var kc = n.cn * dn * dn * dn,
143           fx = dx * kc,
144           fy = dy * kc;
145       p.fx += fx;
146       p.fy += fy;
147     } else if (!n.leaf) {
148       var sx = (x1 + x2) * .5, sy = (y1 + y2) * .5;
149       if (n.c1) forces(n.c1, p, x1, y1, sx, sy);
150       if (n.c2) forces(n.c2, p, sx, y1, x2, sy);
151       if (n.c3) forces(n.c3, p, x1, sy, sx, y2);
152       if (n.c4) forces(n.c4, p, sx, sy, x2, y2);
153       if (dn < max1) return;
154       if (dn > min1) dn = min1;
155       if (n.p && (n.p != p)) {
156         var kc = k * dn * dn * dn,
157             fx = dx * kc,
158             fy = dy * kc;
159         p.fx += fx;
160         p.fy += fy;
161       }
162     }
163   }
164 
165   /**
166    * Applies this force to the specified particles. The force is applied between
167    * all pairs of particles within the domain, using the specified quadtree to
168    * accelerate n-body force calculation using the Barnes-Hut approximation
169    * criterion.
170    *
171    * @function
172    * @name pv.Force.charge.prototype.apply
173    * @param {pv.Particle} particles particles to which to apply this force.
174    * @param {pv.Quadtree} q a quadtree for spatial acceleration.
175    */
176   force.apply = function(particles, q) {
177     accumulate(q.root);
178     for (var p = particles; p; p = p.next) {
179       forces(q.root, p, q.xMin, q.yMin, q.xMax, q.yMax);
180     }
181   };
182 
183   return force;
184 };
185