particle simulation using fast multiple method

Discuss how polywell fusion works; share theoretical questions and answers.

Moderators: tonybarry, MSimon

happyjack27
Posts: 1435
Joined: Wed Jul 14, 2010 5:27 pm

particle simulation using fast multiple method

Post by happyjack27 »

============

EDIT: code now hosted on sourceforge: https://sourceforge.net/projects/octree ... particles/

============

i'm going to post in a new thread here as i should have done before. not sure i can move my other posts into here, but i wouldn't complain if an admin did :)

Note this is an O(N) time n-body simulator. That means that the simulation time only grows linearly with the number of particles. That means you can scale it to HUGEMONGOUS particle counts. and it fits well onto non-uniform memory architectures (NUMA), so hardware-wise, with a few modifications you can put it on massive computer clusters without being bottlenecked by i/o - the i/o pattern stays grid-like, and the i/o demand per channel grows very slowly to not-at-all with the simulation size.

Code for the Particle class:

Code: Select all


public class Particle {
	static final float ATOMIC_CHARGE = 1;
	static final float ELECTRON_MASS = 1;
	static final float PROTON_MASS = 1;
	static final float RC2 = 1; //reciprocal speed of light squared?
	
	float x,y,z;
	float dx,dy,dz;
	public float ddx,ddy,ddz;
	public float oddx,oddy,oddz;
	public int element_number = 0;
	public float mass = 0;
	public float charge = 0;
	public static float[] charge_density = new float[]{-ATOMIC_CHARGE/ELECTRON_MASS,ATOMIC_CHARGE/PROTON_MASS};
	boolean rigid = false;
	
	public Particle(int element_number2, float x2, float y2, float z2) {
		element_number = element_number2;
		x = x2;
		y = y2;
		z = z2;
	}

	public Particle(int element_number2, float x2, float y2, float z2, float mass, boolean rigid) {
		element_number = element_number2;
		x = x2;
		y = y2;
		z = z2;
		this.mass = mass;
		this.rigid = rigid;
		this.charge = mass*charge_density[element_number];
	}

	public float[] getNetForceOnParticle(Particle p) {
		float deltax = x-p.x;
		float deltay = y-p.y;
		float deltaz = z-p.z;
		float rds = (float)1.0/(deltax*deltax + deltay*deltay + deltaz*deltaz); //reciprocal distance squared
		float reciprocal_norm = (float)Math.sqrt(rds);
		
		float electrostatic_force = charge*p.charge*rds;

		//electric force vector
		float ex = electrostatic_force*deltax*reciprocal_norm;
		float ey = electrostatic_force*deltay*reciprocal_norm;
		float ez = electrostatic_force*deltaz*reciprocal_norm;
		
		//magnetic force vector
		float deltadx = dx-p.dx;
		float deltady = dy-p.dy;
		float deltadz = dz-p.dz;
		
		float bx = RC2*(deltady*ez + deltadz*ey);
		float by = RC2*(deltadz*ex + deltadx*ez);
		float bz = RC2*(deltadx*ey + deltady*ex);
		
		float mx = -(deltady*bz + deltadz*by);
		float my = -(deltadz*bx + deltadx*bz);
		float mz = -(deltadx*by + deltady*bx);
		
		return new float[]{
				ex+mx,
				ey+my,
				ez+mz,
				};
	}

	public void combine(Particle p) {
		float rsm = (float)1.0/(mass+p.mass); //reciprocal sum of masses;

		x = rsm*(x*mass+p.x*p.mass);
		y = rsm*(y*mass+p.y*p.mass);
		z = rsm*(z*mass+p.z*p.mass);

		dx = rsm*(dx*mass+p.dx*p.mass);
		dy = rsm*(dy*mass+p.dy*p.mass);
		dz = rsm*(dz*mass+p.dz*p.mass);
		
		mass += p.mass;
		charge += p.charge;
	}

	public void separate(Particle p) {
		float rsm = (float)1.0/(mass-p.mass); //reciprocal sum of masses;

		x = rsm*(x*mass-p.x*p.mass);
		y = rsm*(y*mass-p.y*p.mass);
		z = rsm*(z*mass-p.z*p.mass);

		dx = rsm*(dx*mass-p.dx*p.mass);
		dy = rsm*(dy*mass-p.dy*p.mass);
		dz = rsm*(dz*mass-p.dz*p.mass);
		
		mass -= p.mass;
		charge -= p.charge;
	}

	public void applyForcesToParticle(Particle p) {
		float[] forces = getNetForceOnParticle(p);
		
		p.ddx += forces[0];
		p.ddy += forces[1];
		p.ddz += forces[2];
	}

	public void integrateTimeStep(float dt) {
		if( rigid)
			return;
		
		//do initial setup calculations
		float dtt = dt*dt/2;
		float dttt = dt*dtt/3;
		float rdt = (float)1.0/dt;

		ddx /= mass;
		ddy /= mass;
		ddz /= mass;
		
		float dddx = (ddx-oddx)*rdt;
		float dddy = (ddy-oddy)*rdt;
		float dddz = (ddz-oddz)*rdt;
		
	    //adjustLastTimeStep
	    x += dddx*dttt;
	    y += dddy*dttt;
	    z += dddz*dttt;
	      
	    dx += dddx*dtt;
	    dy += dddy*dtt;
	    dz += dddz*dtt;
	      
		//doThisTimeStep
		x += dx*dt + ddx*dtt;
		y += dy*dt + ddy*dtt;
		z += dz*dt + ddz*dtt;
		
		dx += ddx;
		dy += ddy;
		dz += ddz;
		
		//store force for next time step
		oddx = ddx;
		oddy = ddy;
		oddz = ddz;
		
		ddx = 0;
		ddy = 0;
		ddz = 0;
	}
}
Code for the Octet class:

Code: Select all

import java.util.*;

public class Octet {
   static float threshold = (float)0.1;

   float x,y,z; //the pseudo-particles store the inertias (dx,dy,dz)   
   float delta; //radius
   float dv; //volume
   
   Octet parent = null;
   Octet[][][] children;
   Vector<Particle> particles;
   Vector<Particle> pseudo_particles;
   float num_particles = 0;
   
   public static void main() {
      Octet rootOctet = new Octet(0,0,0,128,(float)0.125,null);
      Vector<Particle> particles = new Vector<Particle>();
      
      //TODO: generate or load a particle list here
      
      iterate( particles,  rootOctet, (float)0.01);
   }

   public Octet(float x, float y, float z, float delta, float min_delta, Octet parent) {
      this.x = x;
      this.y = y;
      this.z = z;
      this.delta = delta;
      this.dv = 8*delta*delta*delta; //since delta is only half the length of a side, we have to multiply by 2*2*2.
      this.parent = parent;
      
      float half_delta = delta/2;
      if( delta > min_delta) {
         children = new Octet[2][2][2];
         for( int dz = 0; dz < 2; dz++) {
            for( int dy = 0; dy < 2; dy++) {
               for( int dx = 0; dx < 2; dx++) {
                  children[dz][dy][dx] = new Octet(
                        x-half_delta+delta*dx,
                        y-half_delta+delta*dy,
                        z-half_delta+delta*dz,
                        half_delta,
                        min_delta,
                        this
                        );
               }
            }               
         }
      } else {
         particles = new Vector<Particle>();
      }
   }
   
   public void addParticleDown(Particle p) {
      num_particles++;
      if( particles == null) {
         children[p.z > z ? 1 : 0][p.y > y ? 1 : 0][p.x > x ? 1 : 0].addParticleDown(p);
      } else {
         particles.add(p);
      }
      
      //now add it at this level
      for( Particle pp : pseudo_particles) {
         if( pp.element_number == p.element_number) {
            pp.combine(p);
            return;
         }
      }
      Particle pp = new Particle(p.element_number,x,y,z);
      pp.combine(p);
      pseudo_particles.add(p);
   }
   public void removeParticleDown(Particle p) {
      num_particles--;
      for( Particle pp : pseudo_particles) {
         if( pp.element_number == p.element_number) {
            pp.separate(p);
            if( pp.mass == 0)
               pseudo_particles.remove(pp);
            break;
         }   
      }
      if( particles == null) {
         children[p.z > z ? 1 : 0][p.y > y ? 1 : 0][p.x > x ? 1 : 0].removeParticleDown(p);
      } else {
         particles.remove(p);
      }
   }

   public static void iterate(Vector<Particle> particles, Octet rootOctet, float deltat) {
      for( Particle p : particles) {
         rootOctet.addParticleDown(p);
         p.ddx = 0;
         p.ddy = 0;
         p.ddz = 0;
      }
      for( Particle p : particles) {
         rootOctet.applyForcesToParticle(p);
      }
      for( Particle p : particles) {
         rootOctet.removeParticleDown(p);
      }
      for( Particle p : particles) {
         p.integrateTimeStep(deltat);
      }
   }

   public void applyForcesToParticle(Particle p) {
      
      //if at bottom, just apply forces
      if( particles != null) {
         for( Particle ps : particles) {
            ps.applyForcesToParticle(p);
         }
         return;
      }
      
      if( getNetPsuedoParticleForceOnParticle(p) > threshold) {
         for( int z = 0; z < 2; z++) {
            for( int y = 0; y < 2; y++) {
               for( int x = 0; x < 2; x++) {
                  children[z][y][x].applyForcesToParticle(p);
               }
            }               
         }
         return;
      } else {
         for( Particle ps : pseudo_particles) {
            ps.applyForcesToParticle(p);
         }
         return;
      }
   }

   float getDistanceSquared(Particle p) {
      float deltax = p.x-x;
      float deltay = p.y-y;
      float deltaz = p.z-z;
      return deltax*deltax + deltay*deltay + deltaz*deltaz;
   }
   
   //TODO: give this function a better implementation.
   public float getNetPsuedoParticleForceOnParticle(Particle p) {
      return num_particles/getDistanceSquared(p);
   }
}
Particle class + Octet class = DONE! An entire plasma simulation program. No visual display, though... or particle list importer... the later should be pretty trivial. could easily be multi-threaded via parallel-izing the loops in the iterate() function.
Last edited by happyjack27 on Fri Dec 13, 2013 3:44 am, edited 1 time in total.

happyjack27
Posts: 1435
Joined: Wed Jul 14, 2010 5:27 pm

Re: particle simulation using fast multiple method

Post by happyjack27 »


happyjack27
Posts: 1435
Joined: Wed Jul 14, 2010 5:27 pm

Re: particle simulation using fast multiple method

Post by happyjack27 »

revision 2, commited to repository:
* you can now load a particle list from a file.

i added the jsonobject class, and made the particle object extends it. then added a util class which extends it too. the util class is now the main entry class. the util class just deserializes itself from the a .json file you select, and this provides the particle list. then it iterates the simulation a bunch of times.

there is also a sample.json file there to show you the format for a particle file.

(note i still haven't set the atomic_charge constant yet. if anyone can tell me what that should be...)

also still no visual display. woudl probably pretty easy to make - just instantiate a JFrame, put a canvass on it, draw on the canvass. (might want to add a 3-d transform so it can be zoomed and panned and rotated.)

EDIT: filled in serialization part, and added save to file function in the Util class, so now the simulation can be saved to disk, to be continued later.

happyjack27
Posts: 1435
Joined: Wed Jul 14, 2010 5:27 pm

Re: particle simulation using fast multiple method

Post by happyjack27 »

version 7 now.

*multi-threading the main loop.
*added basic frame for user interface, with an image panel that shows the particles (no paning,rotating, etc. yet)
*reorganized a bit, into packages

kcdodd
Posts: 722
Joined: Tue Jun 03, 2008 3:36 am
Location: Austin, TX

Re: particle simulation using fast multiple method

Post by kcdodd »

Is this static or full EM?
Carter

happyjack27
Posts: 1435
Joined: Wed Jul 14, 2010 5:27 pm

Re: particle simulation using fast multiple method

Post by happyjack27 »

full em. dynamic particle-to-particle electric and magnetic

i got a few ideas to improve it, like automatically increasing or decreasing the depth of the octree (turning it into a sparse octree), adapting the threshold on a per-particle basis, writting code to shift loads and memory across computes in a cluster, etc. but first i want to make something that you can see and work with. i guess i just have to add some buttons and some pan/rotate/zoom.

kcdodd
Posts: 722
Joined: Tue Jun 03, 2008 3:36 am
Location: Austin, TX

Re: particle simulation using fast multiple method

Post by kcdodd »

What I meant to ask is does it use the electrostatic approximation (and magnetostatic). I.E. instantaneous field E1 = kq2(r1-r2)/(r1-r2)^3.
Carter

happyjack27
Posts: 1435
Joined: Wed Jul 14, 2010 5:27 pm

Re: particle simulation using fast multiple method

Post by happyjack27 »

it uses point-to-point lorentz force and coloumb force.

for interactions with far away particles, multiple particles are lumped together into a pseudo-particle, representing essentially the sum of the forces of the lumped together one. this pseudo-particle is then used in point-to-point lorentz & coloumb force. this little trick saves it from having to do a lot of calculations that ultimately have little to no impact on the accuracy.

more precisely it calulates the force like so (anyone please chime in if they see anything wrong):

Code: Select all


   public float[] getNetForceOnParticle(Particle p) {
      float deltax = x-p.x;
      float deltay = y-p.y;
      float deltaz = z-p.z;
      float rds = (float)1.0/(deltax*deltax + deltay*deltay + deltaz*deltaz); //reciprocal distance squared
      float reciprocal_norm = (float)Math.sqrt(rds);
      
      float electrostatic_force = charge*p.charge*rds;

      //electric force vector
      float ex = electrostatic_force*deltax*reciprocal_norm;
      float ey = electrostatic_force*deltay*reciprocal_norm;
      float ez = electrostatic_force*deltaz*reciprocal_norm;
      
      //magnetic force vector
      float deltadx = dx-p.dx;
      float deltady = dy-p.dy;
      float deltadz = dz-p.dz;
      
      float bx = RC2*(deltady*ez + deltadz*ey);
      float by = RC2*(deltadz*ex + deltadx*ez);
      float bz = RC2*(deltadx*ey + deltady*ex);
      
      float mx = -(deltady*bz + deltadz*by);
      float my = -(deltadz*bx + deltadx*bz);
      float mz = -(deltadx*by + deltady*bx);
      
      return new float[]{
            ex+mx,
            ey+my,
            ez+mz,
            };
   }


happyjack27
Posts: 1435
Joined: Wed Jul 14, 2010 5:27 pm

Re: particle simulation using fast multiple method

Post by happyjack27 »

revision 13;
-added data structures and serialization for elements, particle sinks (aka loss vectors), and particle sources (aka e-guns/ion guns)
note, their functionality isn't actually implemented yet, just the data structures and loading/saving to disk (or serializing/deserializing for network transport).

the long-term idea is you'll be able to control the color, charge, mass, and time-scale of the elements independentely in-situ (i.e. live) through sliders (or you can type in exact figures), and similarly you'll be able to control parameters of the particle sources in-situ. for sources, there's a group number, so you'll control a whole group at once. also for sinks, you'll be able to see the total kinetic energy loss rate per sink group.

sources and sinks are spheres. to represent the vacuum chamber, you just use an "inverted" particle sink. inverted means anything _outside_ of it is lost.

now the configuration file look something like this:

Code: Select all

particles : [
	{x:1,y:0,z:0,dx:0,dy:0,dz:0,mass:1,charge:-1,element_number:0,rigid:false},
	{x:0,y:1,z:0,dx:0,dy:0,dz:0,mass:1,charge:1,element_number:1,rigid:false},
	{x:0,y:0,z:1,dx:0,dy:0,dz:0,mass:1,charge:-1,element_number:0,rigid:false},
	{x:-1,y:0,z:0,dx:0,dy:0,dz:0,mass:1,charge:1,element_number:1,rigid:false},
	{x:0,y:-1,z:0,dx:0,dy:0,dz:0,mass:1,charge:-1,element_number:0,rigid:false},
	{x:0,y:0,z:-1,dx:0,dy:0,dz:0,mass:1,charge:1,element_number:1,rigid:false},
	{x:0,y:0,z:0, dx:0,dy:0,dz:2,mass:1,charge:-1,element_number:0,rigid:true},
],
elements : [
	{element_number:0,charge:-1,mass:1,r:255,g:0,b:0},
	{element_number:1,charge:1,mass:1,r:0,g:255,b:0},
],
sinks : [
	{x:0,y:0,z:0,radius:10,element_number:0,inverted:true}
	{x:0,y:0,z:0,radius:10,element_number:1,inverted:true}
],
sources : [
	{x:0,y:0,z:0,radius:5,element_number:0,inverted:false}
	{x:0,y:0,z:0,radius:5,element_number:1,inverted:false}
],

happyjack27
Posts: 1435
Joined: Wed Jul 14, 2010 5:27 pm

Re: particle simulation using fast multiple method

Post by happyjack27 »

-implemented sources and sinks.
-updated sample.json accordingly

note: particle source production rate is tied to the time_scale for the element. this is because the time_scale effects the movement speed, and thus the loss rate to the sinks. to make it so this doesn't effect particle balance, production rate has to be linked to the time_scale for the element, so that increase in loss rate is tied to a proportional increase in production rate. statistics on loss and production rate are tied to the time_scale, so even though visually the e guns/ion guns will appear to have higher voltage, this will be divided by the time_scale when reporting figures, to get a per-real-second figure.

note also that while particles at a higher time scale will appear to move faster than one at a lower timescale, that is only because their update step is assuming more time has passed. the magnitude of their velocity vector is not actually any greater. and thus they are NOT generating a proportionally higher current.

i think i'm going to skip adding rotation to the display for now. i know too little about quarternions...

prestonbarrows
Posts: 78
Joined: Sat Aug 03, 2013 4:41 pm

Re: particle simulation using fast multiple method

Post by prestonbarrows »

happyjack27 wrote:also still no visual display. woudl probably pretty easy to make - just instantiate a JFrame, put a canvass on it, draw on the canvass. (might want to add a 3-d transform so it can be zoomed and panned and rotated.)
I would suggest using HDF5 for dumping/loading/visualizing your results. It is an open source database format developed by the national labs for handling and visualizing huge data sets.

http://www.hdfgroup.org/HDF5/

VisIt is a free program for viewing the files. It allows for some really powerful functions for post-processing on datasets with millions or billions of points.

http://www.visitusers.org/

happyjack27
Posts: 1435
Joined: Wed Jul 14, 2010 5:27 pm

Re: particle simulation using fast multiple method

Post by happyjack27 »

prestonbarrows wrote:
happyjack27 wrote:also still no visual display. woudl probably pretty easy to make - just instantiate a JFrame, put a canvass on it, draw on the canvass. (might want to add a 3-d transform so it can be zoomed and panned and rotated.)
I would suggest using HDF5 for dumping/loading/visualizing your results. It is an open source database format developed by the national labs for handling and visualizing huge data sets.

http://www.hdfgroup.org/HDF5/

VisIt is a free program for viewing the files. It allows for some really powerful functions for post-processing on datasets with millions or billions of points.

http://www.visitusers.org/
thanks. looks like a lot of learning to do something pretty basic. if i could find a json schema or database definition for one of these, i might add an exporter, but i couldn't find anything of the sort in the documentation. (and it looks like there's a lot of wheel-re-invention going on, at least in the hdf5 one. for example: http://www.hdfgroup.org/HDF5/doc/HL/RM_ ... elds_index -- umm, how about just use odbc?)

i have a basic visual display up now, with perspective projection. i don't have z-sorting yet, but i don't anticipate any issues with that. and i don't have rotation yet because i want to use quarternions. (i might just give in and go with basic matrix multiplication.) other than that, it's all coded (though not tested). though it'd be nice to offload that to hardware rendering on a gpu in the long run. i was thinking jCUDA.

having said that, you can export a snapshot to .json. now that you mention post-processing and all that, i'm thinking i might want to include another export option that would produce more detailed .json, such magnetic and electric forces on particles, ke loss to particle sinks, etc, for potential import into visualizing software. (though you'd probably have to write a converter for the particular visualizer you want to use)

all in all i got a basic display coded, though it's not hardware accelerated, so it's probably slow. down the road i might be ambitious enough to improve upon it, but if anyone is feeling ambitious now - the code is open source, and the import/export format is too. (see "sample.json" https://sourceforge.net/p/octreemethodf ... ample.json)

DeltaV
Posts: 2245
Joined: Mon Oct 12, 2009 5:05 am

Re: particle simulation using fast multiple method

Post by DeltaV »

viewtopic.php?f=3&t=1003&p=52212&sid=c3 ... f59#p52212

Define an angle-about-axis rotational difference quaternion (normalized), which takes an (arbitrary magnitude) "old" position vector into a "new" position vector, relative to some reference frame (or, say, from the frame's x-axis (pos_vec_magnitude, 0, 0) to the "new" position, same number of multiplications),

[[d, D]] = [[cos(angle/2), axis_unit_vector sin(angle/2)]] .

Then, using the general multiplication rule for quaternions (associative but not commutative),

[[c, C]] = [[b, B]] [[a, A]] = [[ba - B.A, bA + aB + BxA]] ,

rotate the position vector (formally, an "axial" vector standing-in for the "polar" position vector) by employing a "similarity transform",

[[0, pos_vector_new]] = [[d, D]] [[0, pos_vector_old]] [[d, -D]]

where [[d, -D]] (note sign) is the "conjugate" (for normalized quaternions, the same as the "inverse").

happyjack27
Posts: 1435
Joined: Wed Jul 14, 2010 5:27 pm

Re: particle simulation using fast multiple method

Post by happyjack27 »

okay, looks simple, except i don't know what 90% of that stuff is. "the general rule for multiplying quarternions". yeah, that's what i'm looking for. if someone could just write down the algorithm for me in a standard programming language such as C or java, that would be all i need.

i have two angles, a rotation about the x axis, lets call it "xtheta", and a rotation about the y axis, lets call it "ytheta", and a WHOLE lot of points, (x,y,z). in the end i want a 3x3 matrix:
[
[xx,xy,xz],
[yx,yy,yz],
[zx,zy,zz],
]

such that i can take that that array of points, p, and go like:

Code: Select all

for( int i = 0; i < p.length; i++) {
   p[i].transformed_x = xx*p[i].x + xz*p[i].y + xz*p[i].z;
   p[i].transformed_y = yx*p[i].x + yz*p[i].y + yz*p[i].z;
   p[i].transformed_z = zx*p[i].x + zz*p[i].y + zz*p[i].z;
   drawPointAt(p[i].transformed_x,p[i].transformed_y,p[i].transformed_z);
}
(NOTE: it's important that in the end i get a simple transform like this, that takes very few computations, instead of having to run through the whole process for each point. because there will be millions of points. having to do one more operation per point means having too do a million more operations.)

i do not want a library of classes for dealing with quartenions and classes of vectors and matrices, made up of 50 files, each of which is 100s of lines of code. I can find those way too easily on the internet. i want a maximum of 100 lines of code, that does only what it needs to do - convert 2 numbers: xtheta and ytheta, to 9 numbers:xx,xy,xz,yx,yy,yz,zx,zy,zz.

i'm tempted to just apply a xtheta matrix ala sin and cosine, followed by a ytheta matrix. simple and effective. but i'm curiuos how much the order matters...they'd be dragging their mouse to rotate on x and y simultationsly.

DeltaV
Posts: 2245
Joined: Mon Oct 12, 2009 5:05 am

Re: particle simulation using fast multiple method

Post by DeltaV »

happyjack27 wrote:okay, looks simple, except i don't know what 90% of that stuff is.
Traditional (Gibbs-Heaviside) 3-D vector analysis, as was commonly taught to engineers of my generation.
The twist is the quaternionic pairing of a scalar and an "axial" vector ("pseudovector"), [[s, v_axial]], similar to the 2-D complex-plane pairing of a real and an imaginary number. Axial vectors involve rotation. Ordinary ("polar") vectors involve translation, a prime example being a radial position vector p. Historical accident (or ****-up, or conspiracy) that the same arrows/notation were used for both. Axial vectors embody not just magnitude and direction, but also handedness or "chirality" or "parity".
happyjack27 wrote:"the general rule for multiplying quarternions". yeah, that's what i'm looking for. if someone could just write down the algorithm for me in a standard programming language such as C or java, that would be all i need.
My source code is 1000 miles away at present and not accessible. It embodies several heuristic tweaks to improve the user experience. Maybe I can remember the core of it and post some pseudocode later. But it is getting late. The tweaks take up most of that routine. Variable mouse input scaling based on zoom factor, center offset, etc. The product rule itself is just a few lines of code. Cross product, dot product. Separate routine.
happyjack27 wrote:i'm tempted to just apply a xtheta matrix ala sin and cosine, followed by a ytheta matrix. simple and effective. but i'm curiuos how much the order matters...they'd be dragging their mouse to rotate on x and y simultationsly.
You shouldn't think of it as sequential rotations. Infinitesimal rotations commute, finite rotations do not. One display-update cycle is as "infinitesimal" as your code can get. You should think of it as rotation about some axis contained in the x-y plane of the screen, passing through the projected (onto the 2-D screen) center of your particle cloud.

Post Reply