Possible Charged Particle Field

The N,S,E,W,F,B positions are not real positions in the world system. They are just calculation points in the math. We do not need to track them over time because they are just a way to determine the charge effects from different sides of the particle. I defined them to be set to the dimensions, not the particle, because I am using them to easily split the vectors.

By split the vectors, I mean to be able to extract their X, Y and Z components. For a given charge vector, the Y might be the linear velocity component and the X and Z are the spin components. This is determined by looking at which dimension points to the center of the particle and that will be the linear velocity dimension for that position.

With respect to what is considered the center, it actually varies during the calculations. For the algorithm described above, let's call it the Velocity algorithm (because it calculates velocity and spin), the current particle is considered center. However, for the Charge Density algorithm, the other particle is the center. This is just an artifact of the math. You could calculate it either way but I think it is easier if the source of the charge is center and the measurement point is relative to that.

The calculations are made, not for a moment, but for a time interval. I have not included any time math in that description but it will need to be handled. Essentially, time is handled by determining how much real time has passed since the last frame, looking at the units that you specify for things like velocity (m/s, km/h, etc) and converting them to the amount of time that has passed. So If we specified that a particle has a velocity of (0, 1, 0) which is defined as m/s, and the last frame was 1/60th of a second ago (60 frames per second), then the actual velocity to be used will be (0, 1/60, 0) per frame.

It is not the world system that is rotating, but each and every particle. It is possible that some will not have spin, of course. You have to find the right perspective to understand this stuff. The math above does perform a lot of rotations because it needs to arrange everything with respect to the current particle being calculated. That is just for the math though. It does not, in any way, actually rotate the particle. You will never see those rotations in the scene. You will see the rotations that are a result of the math, but not the rotations performed during the math.

When I say we want to orient a particle with respect to the current particle, it is just rotating the temporary copy of the position of that particle. Not the actual position of it. In fact, the positions and orientations used in these calculations are never written back to the particles. We read them in, but never write them back out because the result of these calculations are velocities and spins (1 each per particle). Velocity and spin are properties of the particle which are applied every frame. The application of those properties can happen at a different frequency than the calculations.

Suppose we find that these calculations are excruciatingly slow. If we were updating the position and orientation every time we calculate, then there is nothing we can do about it being slow. However, by disconnecting the calculation from the update, we can choose to calculate every 2nd frame, or 3rd, or what-ever. The velocities and spins are still applied every frame, but they themselves are not updated that often. This is one way to set the precision of the model.

I do have a rough outline of the Charge Density algorithm. It still needs work, and it is only an approximation, but it should be fairly good. The beauty of separating these 2 algorithms is that we can replace either one if we find a better way. I'll try to expand on the charge density algorithm and post it.

Charge Density Algorithm

Code:


// this array is used to approximate the charge density from equator to pole
// the values are percentages
// it is mirrored at the equator
// it handles the fact that the strongest emission is as 30° N/S, not at the equator itself
// the length of this array is arbitrary, you can include more values for more precision
CHARGE_DENSITY = [80, 90, 100, 90, 70 , 50, 30, 10, 5, 0];
// this value specifies the max distance from a particle that its charge field can be felt
// this is used for performance reasons
MAX_CHARGE_DISTANCE = ?
// this value is used to set the maximum amount of charge that a particle can emit
// it should be replaced by a 1/r^4 calculation at some stage
MAX_CHARGE_DENSITY = ?
// this value is used to specify the velocity applied by the max charge density
// this is a way to convert between a density and a velocity
MAX_CHARGE_STRENGTH = ?

// find vector to p2
v1 = p1 - p2
if( v1.length < p1Radius )
{
   // collision of actual particles
   return collision( p1, p2, ... );
}
// find distance as a factor of max charge distance
d = v1.length / MAX_CHARGE_DISTANCE;
if( d < 1 )
{
   // find equator to pole ratio
   v2 = [v1.x, 0, v1.z]
   l2 = v2.length / MAX_CHARGE_DISTANCE;
   y2 = v1.y / MAX_CHARGE_DISTANCE;
   if( y2 == 0 )
   {
      // perfectly on the equator
      epr = l2;
   }
   else
   {
      // somewhere between equator and pole, including perfectly at the pole
      epr = l2 / abs( y2 )
   }
   // find the charge density
   i = Integer.valueOf( epr * (CHARGE_DENSITY.length-1) );
   epf = CHARGE_DENSITY[i] / 100.0;
   df = d * MAX_CHARGE_DENSITY   // TODO 1/r^4
   v1.normalize();
   v1.scale( d * epf * MAX_CHARGE_STRENGTH )
   // TODO find resultant spin of the charge field, if we want to use that
}
else
{
   // outside of effective radius
   v1 = [0, 0, 0]; // = [0, 0, 0, 0] if using charge photon spin
}
return v1;


This one is a bit rougher than the last, but also closer to real code.

I have started to create a Javascript implementation of these algorithms. This will force me to look deeper at certain parts and get a better understanding of what data I have available at various stages. This will be the fallback implementation if shaders are not available on the running system. It does not use shaders at all, so I expect it to be a fair bit slower than the final form. However, I am designing it in such a way that it is shader friendly. I would do some things differently if I was not, but this will make it easier to port to the shader language or even OpenCL.

I can give you a better idea of how these functions look now. The Velocity algorithm is the driver. It is what the application must call to calculate the velocities and spins for all particles. The Velocity algorithm calls the Charge Density algorithm to perform its job. From the applications perspective, these algorithms use the position and orientation of all particles as read-only input, and the velocity and spin of all particles as write-only output. These 4 datasets are what we are going to transport through images (or textures) in the shader.

Each dataset is an array of floating point values. The values are grouped into a set number of items to form some larger entity. For example, the position array uses sets of 3 values to form a vector: X, Y, Z. The spin array uses 4 values to form a quaternion: X, Y, Z, W. When in texture form, each group of 3 or 4 values are called a texel or pixel. So the X, Y, Z(, W) are the same as R, G, B(, A) channels if it was real image data.

The function declaration for the Velocity algorithm might look like this:

Code:


void calcVelocityAndSpin(float[] positions, float[] orientations, float[] velocities, float[] spins);


It won't actually look like that in the shader, but it will equate to the same thing. It also doesn't look like that in my Javascript implementation but it does provide the data in the same way (as in 4 arrays) even though it is slower to do it like that in JS. An OpenCL implementation would look like that, but with some read-only and write-only flags added in.

Each item in these arrays, where an item is a group of 3 or 4 values, contains the data for a particular particle. The 1st item contains data for the 1st particle up to the Nth item for the Nth particle. To retrieve the correct data for a particular particle, you multiply its item index by the number of values in the item like this:

Code:


int posIndex = itemIndex * 3;
int orientIndex = itemIndex * 4;


So if we wanted to get the 3rd items position, we would go to the 3*3 = 9th value and take 3 of them to form a vector. If we wanted the 3rd items orientation, then we go to 3*4 = 12th value and take 4 of them to form a quaternion. This is generally how you pass around large amounts of data to the GPU.

The Charge Density function only needs the position and orientation of the emitting particle, the position of the measurement point, and a vector to store the result. It will look like this:

Code:


void calcChargeDensity(vec3 sourcePoint, vec4 sourceOrient, vec3 measurePoint, vec3* dst);


I have just discovered that GLSL does not supply any support for quaternions, so I will have to code up their math, much like I have for matrices. Not a big deal, but I also might see if I can avoid the quaternion math as much as possible because they are expensive operations. It might be better to convert the quaternions into matrices. I'm not sure at this point.

Nevyn wrote. I can give you a better idea of how these functions look now. The Velocity algorithm is the driver. It is what the application must call to calculate the velocities and spins for all particles. The Velocity algorithm calls the Charge Density algorithm to perform its job. From the applications perspective, these algorithms use the position and orientation of all particles as read-only input, and the velocity and spin of all particles as write-only output. These 4 datasets are what we are going to transport through images (or textures) in the shader.

Wow Nevyn, I love it when you talk data like that.

Thanks for the explanations. So it seems we will be building the Charged Particle Interaction Model (CPIM) from the existing Proto-Planet (aka webgl gpgpu debris) source after all. Currently CPIM has four shaders:
1. Fragment shader for position -->
2. Fragment shader for velocity -->
3. Particles vertex shader -->
4. Particles fragment shader -->
Makes perfect sense to add orientation to position; and spin to velocity. Spinning does add complexity.

CPIM includes one unusual Three.js source file. GPUComputationRenderer.js, author yomboprime https://github.com/yomboprime, https://github.com/yomboprime/GPGPU-threejs-demos. Did you happen to look at it? The author has four other related particle examples. I'm hoping it may have some bearing on your quaternion/matrix thinking.
.

I haven't looked at the existing shaders since you first mentioned this project, so I'm not sure how it will fit in, but I would expect it to be a new shader, not added to an existing one, or maybe added to the shader that contains the gravity math. We can't split the datasets that I have defined, they are all needed together in the same shader.

I'm not convinced that we need to use this example as the basis. We could create a new one from scratch, but there is nothing stopping us from doing so either. So I am happy to run with it for now. I do like the way it renders the spheres as glowing. I will probably attempt to add my charge shaders from AtomicViewer to that to show the equatorial emission.

I did see the GPUComputationRenderer and was pleased to do so. However, upon looking over it again, I'm not sure if I want to use it or not. It is a bit inefficient in its memory usage. Basically all textures are packed with 4 values per texel. Half of mine only need 3, so 1 byte is wasted per texel in those 2 datasets. That might not sound like much to worry about, but every little bit is worth worrying about at this level. That is 2 bytes per particle that will need to be written to GPU memory and copied back to main memory, 60 times a second, for no reason. However, I might get lazy and just go with it for convenience.

It won't help with quaternions or matrices though. That is all inside of the shader, so only GLSL can help with that and they don't have the routines available yet. The math isn't difficult, just expensive computationally. Lots of multiplications, and a multiply operation is about 9 times more expensive than an add operation at the processing level (CPU or GPU). However, converting quaternions into matrices can be almost as expensive so I need to make sure that it is worth the effort.

My first approach is usually to do it the way I see it and then optimize later. Once you have an algorithm working, it is easier to see how to change it.

Since I am about to commit the Javascript implementation of these algorithms to GIT, I thought I might go over how I setup my JS files.

In Javascript, when you define a variable, it has some scope. If you are inside of a function, then that function is the scope. If you are inside of a code block (surrounded by { and }), then that block is the scope. If you are outside of any function, then it is in the global scope. The variable does not exist outside of that scope. You can not use it if you can't see it. If the variable is in the global scope, then you can use it anywhere. If it is in a function, then you can use it anywhere in that function. If it is in a code block, then you can only use it in that code block and not anywhere else, even in the function that contains the block.

Now, global variables all get put into the same place (they are actually stored as properties of the window object that is supplied by the browser). Because of this, you can get name collisions and not realise that you are actually effecting other code just because they both declare the same variable. I have spent days tracking down such bugs.

To get around this, I create my own module that contains all of the code for that module and keeps it in its own scope. You have already used such modules. THREE is the module for ThreeJS. That is why you have to specify THREE in front of any class that you use from it, such as THREE.Vector3. The module object itself is still stored in the global scope, so you can access it anywhere, but now anything within that module is kept separate from the global scope.

This could be as simple as declaring an object like this:

Code:


var MyModule = {
  myVar: 1.0,
  myOtherVar: 'is a string'
};


However, I take it a step further so that I can create a private space so that no other code can access certain parts of that module. To accomplish that I use a function. Javascript is a bit of a strange language and it has a feature called Closures. Closures are basically a hierarchy of scope. If you are inside of a function and use a variable that you have not declared in that function, then it will look at the scope that that function was declared in, if not found it will look at the scope of that scope, etc, etc. If it still hasn't found it, then it will look at the global variables.

Here's an example:

Code:


var globalVariable = 1.0;
function myFunction()
{
  var functionVariable = 2.0;
  if( someCondition )
  {
    var blockVariable = 3.0;
  }
  // note that we do not need an if, or while, or any control structure to create a block
  {
    var anotherBlockVariable = 4.0;
  }
}


Now, those variables all have unique names so there is no problem. Let's look at another example where they all use the same name.

Code:


var variable = 1.0;
function myFunction()
{
  var variable = 2.0;
  if( someCondition )
  {
    var variable = 3.0;
  }
  // note that we do not need an if, or while, or any control structure to create a block
  {
    var variable = 4.0;
  }
}


That is perfectly usable code. I don't recommend doing that as it is confusing and prone to errors, but it will work and all of those variables are unique and will not collide at all. But if we forgot to use the var declaration, then things get weird.

Code:


variable = 1.0;
function myFunction()
{
  variable = 2.0;
  if( someCondition )
  {
    variable = 3.0;
  }
  // note that we do not need an if, or while, or any control structure to create a block
  {
    variable = 4.0;
  }
}


The first usage of variable will actually create it, even though it has no var declaration. All of the other uses will use that global variable though, so we are getting collisions. You can use this to your own advantage and I do so sometimes. What I want to show now though, is how to create a private space.

We accomplish that by creating our module inside of a function. Therefore, we can create variables in that function that are private to it. Any code in the module can access them, but nothing outside of that module can.

Code:


var MyModule = (function()
{
  var module = {
    myVar: 1.0,
    myOtherVar: 'is a string'
  };

  return module;
}());


That is the outline, but it does not actually create any private variables. The module variable is actually private, but we are returning it from the function to make it accessible as a global variable called MyModule. So in other code we could use it like this: MyModule.myVar = 2.0; Locally though, we must refer to it as module because MyModule may not exist yet. You actually can use MyModule in some places, but it is best to use the module local variable because then we can easily change the name of the module without having to replace its usage inside of its own code.

The private scope is created by declaring the module inside of a function, but we also have to immediately execute that function. To do that, we just add the parameter declaration to the end of the function. This module does not take any parameters but if we wanted to use THREE, then we would pass it in to this function like this:

Code:


var MyModule = (function( THREE )
{
  var module = {
    myVar: 1.0,
    myOtherVar: 'is a string'
  };

  return module;
}( THREE ));


Notice that I had to add THREE to the argument list in the function declaration (at the top) and also in the parameter list to call the function (at the bottom). We could actually rename THREE at this point if we wanted to. It would only be renamed within this function.

Code:


var MyModule = (function( FOUR )
{
  var module = {
    myVar: 1.0,
    myOtherVar: 'is a string'
  };

  var myVector = new FOUR.Vector3();

  return module;
}( THREE ));


I don't recommend doing that though, unless you absolutely have to.

The whole function and call declaration is wrapped in ( and ) just to make sure that it is executed.

Note that myVector is a private variable and can not be accessed from anywhere outside of this module. This can be very handy since that variable will exist for the life of the application. We can create objects and then use them within the functions and classes exposed by the module. They act like static variables in other languages. We might use them as constants, but they are actually variables so you can change their values if you want to.

I am getting a bit too advanced here, but I just want you to be able to recognize the structure and have at least some notion of why it is being used.

I should probably talk a bit about Javascript Classes as well. The first thing you need to know is that JS does not have classes. The second thing you need to know is that I use classes in JS a lot.

Confused yet?

Javascript is what is called a Prototyping Language. It doesn't have strict requirements on what any object can store within itself. You can even add and remove them as you go. In a language like Java, you define what a class can store and it will only ever store that data and you can not remove it. You can set it to null, so it has no value, but the variable itself will always exist on every object of that class.

So how do we use classes in JS? Well, JS uses prototypes to accomplish the same thing (with differences). You declare a function, and then attach a prototype to that function and the prototype becomes an extension of that function. To use the function as a class, you just use the new keyword to create an instance. e.g. var myVector = new THREE.Vector3();

THREE.Vector3 is actually a function and it has a prototype and that prototype contains all of the functions, which we now call methods, for that class. There is only 1 prototype, but every new object from that class contains a reference to it and can use the shared data and methods in it.

Let's look at an example for some clarity:

Code:


...

module.MyClass = function( name )
{
  this.name = name;
};

module.MyClass.prototype = Object.create( {} );

module.MyClass.prototype.getName = function()
{
  return this.name;
};

...


We could then use it like this, assuming it was declared in a module called MyModule:

Code:


var cls = new MyModule.MyClass( 'Nevyn' );
var name = cls.getName();
// of course, we can also access the name field directly like this:
name = cls.name;
// but I wanted to show classes, not field access


So what's the big deal? Well, it is more efficient than the alternative and we can set up class hierarchies by extending one class with another.

Firstly, the old way, before prototypes were added to JS:

Code:


...

module.MyClass = function( name )
{
  this.name = name;
  this.getName = function()
  {
    return this.name;
  };
};

...


While it is quite a bit less code to write, the problem is that every object created by this function will create its own function called getName. If you create 1000 instances of MyClass, then there will be 1000 instances of the getName function. However, with a prototype, there is only 1 instance of the getName function and all 1000 objects use it. This is much better on memory usage and much faster to execute because we don't have to create 999 functions that all do the same thing.

We can also create sub-classes. This allows you to extend one class with another. You do this by supplying the prototype of the class you want to extend when you create the prototype of the class that is extending.

Code:


...

module.MySubClass = function( firstName, lastName )
{
  module.MyClass.call( this, firstName );
  this.lastName = lastName;
};

module.MySubClass.prototype = Object.create( module.MyClass.prototype );

module.MySubClass.prototype.getName = function()
{
  return this.name + ' ' + this.lastName;
};

...


We can override methods of the parent class in this child. Instances of the parent class still use its own version of the method, but instances of the child class use their own version. Any methods not overridden are still available in the sub-class.

Again, it is not too important to be able to use this stuff, but you will need to recognize it to understand my code.

I'm getting worried about you, Airman! It's been a few days since I pushed my code to GIT, and I haven't heard from you. Are you still OK?

By the way, I had a look at the existing code and found a few things that I had not thought about. This will change the way the code is called and I still need to figure out exactly what the BufferGeometry object is going to store. It has to have something in it, but I don't really need it.

I am thinking that I will create another app along side this one for testing my implementation. I don't think this example is going to work for us but I will revisit that at a later date. I still need to get a good understanding of how to fit it into a shader, so my perception may change as I figure that out.

The way I setup my code, without looking at the existing code, expected a list of Object3D instances where each one is a sphere, for example. The actual visual representation is irrelevant, but the fact that it is its own Object3D is critical. I then take the position and orientation from each one, encode it into a data structure (array for me, texture for you), and send it to the code to calculate the velocities and spins of all particles, decode the data structure and write the velocity and spin back to each particle.

The existing code works differently, as it puts all of the positions into the BufferGeometry and operates on them in the shader. This will not work for me, as I want data to come back out of the shader, not just to change the position passed in to it.

The other problem I found is that when you add the data as a shader Attribute, you don't get to see all of the values for all of the other particles in the shader. You just see the value for the current vertex being processed. This is very different from OpenCL, where you get arrays which contain all data for all vertices, and have to find the correct data for the vertex being processed. So I need to change the way I think about this (which has been very OpenCL centric so far).

I'm getting worried about you, Airman! It's been a few days since I pushed my code to GIT, and I haven't heard from you. Are you still OK?

Hey Nevyn, thanks for asking, I'm fine. I was getting worried about you. I was occupied with the Book of Shaders tutorial. Silly me, I assumed bitbucket would send me an e-mail message when you changed CPIM.

Ok, I did a pull. I'll be re-reading, trying to catch up.
.

I wasn't sure that it wouldn't send you that email. I don't get them at work when others push to the same project, but I do get them when you push to this one. Maybe it is because I am the owner of the project and have given you write access to it.

Please note that there isn't really anything for you to do at the moment. Just try to understand what I have done so far, but remember that it is not for shaders, so some things will be different than the tutorials you are looking through. I tried to make sure that I added comments to explain various parts, but if any are confusing then just ask.

Charge Density Model

Given a particle, determine the charge field effects of other particles on it.

Conditions and setup:

Each particle is treated as if it was in a set orientation such that its North pole is +Y, East is +X and Front is +Z.
We want all calculations to be performed in this orientation so that we can split some vectors easily.
We will need to transform all other particles into this virtual coordinate system.
We do that by taking the orientation transform of this particle, inverting it, and applying to each other particle.
This will not actually change those particles, but will change the values in the temporary data we have for them.

Algorithm:

Code:

   Calculate sum of charge densities from other particles at the cardinal positions on this particle (N,S,E,W,F,B).
   For each position:
     Split charge vector into 3 components: 1 linear velocity component that points to the center of the particle; 2 spin components that point in the other 2 dimensions.
     The actual dimensions used are different for each position:
       North: V=-Y, S=X,Z
       South: V=+Y, S=X,Z
       East:  V=-X, S=Y,Z
       West:  V=+X, S=Y,Z
       Front: V=+Z, S=X,Y
       Back:  V=-Z, S=X,Y
     Where a vector with the same sign as the linear velocity dimension means a repulsion and the opposite means attraction (or do we ignore attraction and just let the other side apply repulsion?).
     The spin components work in both directions and just changes the angle from positive to negative.

///\\\///\\\///\\\///\\\///\\\///\\\

// Sorry Nevyn, I’m hung up. I either don’t agree or don't understand the “calculating sum of charge   
// densities” algorithm above. I bolded two spots which I’ll cover in detail below. For the time being  
// I guess it's probably better to discuss this here rather than through pushes and pulls, a future work
// routine. Please tell me where I'm going wrong.

We want all calculations to be performed in this orientation so that we can split some vectors easily.
// For the “charged density received” calculation we have at least two particles:
// 1. a single charge receiving particle which we will fix with 6 cardinal locations; and
// 2. n-1, all other remaining particles will be charge emitters. Unless I’m mistaken, I believe
// all “charge density received” calculations described are performed with the receiving particle
// in a single fixed orientation so that all incoming charge from all other particles present will be
// summed at the same six cardinal points simultaneously. On 10 July above you included a
// charge emission function. How is that function used in this algorithm?  

///\\\///\\\///\\\///\\\///\\\///\\\

     The actual dimensions used are different for each position:
       North: V=-Y, S=X,Z
       South: V=+Y, S=X,Z
       East:  V=-X, S=Y,Z
       West:  V=+X, S=Y,Z
       Front: V=+Z, S=X,Y
       Back:  V=-Z, S=X,Y
// How can V = +/-X, +/-Y, or +/-Z?
// Each set of “charge density received” computations begin by “fixing into position” the
// receiving particle’s cardinal positions, which must be in accordance with the receiving
// particle’s spin axis. The direction of V is determined by the direction from the center of
// the emitting particle to the location of charge received (at the cardinal position).
// I believe those directional lines for one emitter and one receiver can be found by:
// 1. (RN-EC)X + (RN-EC)Y + (RN-EC)Z.
// 2. (RS-EC)X + (RS-EC)Y + (RS-EC)Z.
// 3. (RE-EC)X + (RE-EC)Y + (RE-EC)Z.
// 4. (RW-EC)X + (RW-EC)Y + (RW-EC)Z.
// 5. (RF-EC)X + (RF-EC)Y + (RF-EC)Z.
// 6. (RB-EC)X + (RB-EC)Y + (RB-EC)Z.
// For example, RN means the receiving particle’s North position. EC is the charge
// emitting particle's center.
// Note, none of those lines are parallel, they diverge. When the emitter and receiving
// particles are close together the V directions to the opposite cardinal location pairs,
// N-S, E-W, or F-B emission lines approach right angles. The further apart the two
// particles are separated, the more parallel the six lines become. This doesn't seem
// to be the same as you've described it.
// What do you think?
.

We definitely do not want to use GIT for discussions. GIT is for source code and source code only.

The charge emission function is used to calculate the charge emission of 1 particle upon 1 other particle, or at least some point that we define. It is used in this algorithm to calculate the charge effects from all other particles, from the 6 cardinal positions. Therefore, each other particle is measured 6 times for that 1 particle.

The directional lines can easily be found without re-orienting the particles, but what you can't do so easily is split the vectors up into linear velocity and spin components. If we did not re-orient, then we would have to rotate all of the vectors (or use sin and cos which is the same thing) in order to split them which amounts to a lot more math and hence, a lot more time to compute. I am already doing more computations than I need to by re-orienting within the charge emission function. I will more-than-likely change that to perform all re-orientations in the other algorithm and pass in the already oriented particle. That way, I only need to re-orient once per other particle per emission receiver. I'll look into that when I have it working (or close to it).

I'm not sure if that covers your problems or not. This is not a simple algorithm, so please keep asking until we both understand.

To give you an idea of why I am splitting vectors, you can look at the math for the trajectory of a ball thrown horizontally and vertically, say at 45°. You split the vector into the horizontal and vertical components so that you can apply gravity to the vertical component without touching the horizontal. I am doing a similar thing by saying that the vector from the cardinal point to the center of the particle (of that cardinal point, not the emitter) represents the linear velocity and the other 2 dimensions represent spin or I should say tangential velocity. Since the cardinal points are in different directions from the center, they each have a different direction for the linear velocity component. Even when 2 sides use the same dimension, they differ in whether positive is repulsion or attraction.

With respect to summing the charge vectors at each cardinal point, we are trying to find the net charge influence from that point upon the particle, so we must sum them. I do have reservations about the relative position of the measurement point because we should be ignoring emitters on the other side of the particle. I have thought about it and will implement it because those forces are already handled by the cardinal point on that side. I just have larger issues to think about before that, so please remind me if I forget. That's why I put a comment in my notes about ignoring attractions.

Airman wrote:I believe all “charge density received” calculations described are performed with the receiving particle in a single fixed orientation so that all incoming charge from all other particles present will be summed at the same six cardinal points simultaneously.

Well that is an interesting sentence because it is full of correct and yet wrong.

Both parts of that sentence are correct, but one does not lead to the other.

Airman wrote:all “charge density received” calculations described are performed with the receiving particle in a single fixed orientation

is correct, and

Airman wrote:all incoming charge from all other particles present will be summed at the same six cardinal points simultaneously

is correct, but...

You can't join them together because they are unrelated to each other. We use the fixed orientation to split the vectors. Summing the charge would be both required, and the same math, even if we calculated the charge some other way.

.
Given my full of correct - yet wrong statements, I'll try again.



Here is my understanding of a single iteration of the calculating sum of charge densities algorithm. The algorithm sums the amount of charge received by one of the charged particles (cyan) from all the other charge emitting particles (yellow) present. The three frames show the charge density received at the receiving particle’s cardinal positions from each of the three emitting particles. Note that the receiving particle is fixed in space – its orientation will remain the same as shown during this algorithm iteration. Each line - drawn from the center of the emitting particles to each of the cardinal positions of the receiving particle - shows the linear direction of charge emission, V; we do not see the two orthogonal spin components also present at each cardinal position.

Is that correct? Or does the algorithm reorient the receiver’s cardinal points for each emitting particle?  
.

That shows the lines connecting the emission point to the cardinal points. But it does not show the velocity vectors. The velocity vectors are along the same direction as those lines, but are based at the cardinal point. This is because the vectors are a representation of the force applied at that cardinal point, from the other particle. Therefore, each velocity vector (or charge vector as I have called them before) starts at its cardinal point and points away from the center of the emitting particle with a length proportional to the force applied.

Take your center image above and I will number the lines 1 to 6 from top to bottom of the image so the order will be magenta, green, magenta, blue, green, magenta (is it meant to be red?). They all come from the same emission point but you can see that the force applied by line 3, the middle line, is almost completely linear velocity. However, lines 1 and 6 are almost completely tangential velocity, which will cause spin (in most cases but not if the charge is even on both sides, as they will cancel each other out and actually provide linear velocity).

So my algorithm measures those charge vectors for all other particles and stacks them together at those cardinal points. Then they are all added together to find the net charge influence at that point. That vector is then split into linear velocity and tangential velocity components.

Re-orienting Particles

The point of re-orienting particles is so that we can treat the receiver particle as if it was in a known orientation. That is, as if it did not actually have any orientation of its own. It is all about the receiver, not any of the other particles. However, it is the other particles that we have to re-orient because the orientation does not effect the location of the receiver particle, because the orientation is applied from the location of the receiver. It is effectively at (0, 0, 0). It does effect all other particle locations though, because they are definitely not at (0, 0, 0) with respect to the receiver particle. It is equivalent to taking your images above and rotating them with the cyan particle at the center of rotation.

It is basically saying 'I am the receiver and I am the center of everything and I define what is up and what is down.'. It puts everything else into the receivers local coordinate system or frame of reference, if you want to call it that.

Airman wrote:Each line - drawn from the center of the emitting particles to each of the cardinal positions of the receiving particle - shows the linear direction of charge emission, V; we do not see the two orthogonal spin components also present at each cardinal position.

While it does show the direction of charge emission being measured, that direction already contains both the linear and tangential components. You are thinking about a head-on collision but that is not what I am trying to accomplish. At times, there will be head-on collisions, which just means that all force is being applied along the linear velocity dimension for that cardinal point. The other 2 dimensions will be 0.

Let me paint a picture. I'd love to create it but it would take me too long. I need a vector scratch pad app.

We have a receiving particle at (0, 0, 0) with a radius of 1 and it has no orientation values so its north pole is at (0, 1, 0). The cardinal points are located here:

N: (0, 1, 0)
S: (0, -1, 0)
E: (1, 0, 0)
W: (-1, 0, 0)
F: (0, 0, 1)
B: (0, 0, -1)

Suppose we have an emitter particle at (5, 0, 0) and it also has no orientation so they are equator to equator.

From the East cardinal point we are going to get all charge along the X dimension. That is the linear velocity dimension for the East position so it is all linear velocity. The Y and Z components of the charge vector will be 0.

Suppose we have another emitter at (1, 5, 0) and it does have an orientation that rotates it 90° about the Z axis, which means that it is emitting down straight at the East position of the receiver.

From the East point, all charge will be along the Y axis and none along the X or Z. You can clearly see that this is a tangential velocity and must cause spin, not linear motion.

Suppose we had yet another emitter that was at (5, 4, 0) and oriented so that it also emits straight at the East cardinal point of the receiver, which would be at a 45° angle.

From the East position, half of the charge vector will be along the X axis and half will be along the Y axis (none on the Z) so we have both linear velocity and tangential velocity in the same charge vector.

That is why I want to split them and that is why I want them oriented along the axes. I can just grab a single value (from whichever dimension I care about at the time) instead of extracting the angles from absolute positions which would involve all 3 values in a calculation for a single dimension.

There is method to my madness!

Let's look at the charge vectors from above to see the larger algorithm.

R = receiver
E1 = emitter at (5, 0, 0)
E2 = emitter at (1, 4, 0), I changed this one a little bit
E3 = emitter at (5, 4, 0)

Relative to the East cardinal point each emitter is at:
E1 = (4, 0, 0);
E2 = (0, 4, 0);
E3 = (4, 4, 0).

The distance from East position to each one is:
E1 = 4
E2 = 4
E3 = sqrt( 4*4 + 4*4 ) = 5.65

Therefore the E1 and E2 particles have equal influence, but in different dimensions and E3 has a slightly smaller influence in both dimensions that E1 and E2 have.

For the sake of simplicity, let's just make the charge vector length equal to: 10 - distance. This is not how it would be in the real algorithm, but it doesn't really matter for this example.

The charge vectors would then have the following lengths:

E1 = 6
E2 = 6
E3 = 4.35

I just needed to get E3 as a smaller value than the others.

The charge vectors would be:

E1 = (-6, 0, 0)
E2 = (0, -6, 0)
E3 = (-3.07, -3.07, 0)

We sum these vectors to find the net charge vector.

   (-6, 0, 0)
+  (0, -6, 0)
+  (-3.07, -3.07, 0)
=  (-9.07, -9.07, 0)

For the East cardinal point, the linear velocity is in the -X dimension, so we take the X value: -9.07; as the linear velocity and put it into its own vector:

V = (-9.07, 0, 0).

The spin components are in the Y and Z dimensions so we take them and put them into their own vectors:

SY = (0, -9.07, 0),
SZ = (0, 0, 0).

At this stage, we need to merge all of the cardinal points together into a single linear velocity and a single tangential velocity which becomes the new velocity and spin for the receiver.

.
Thanks Nevyn, I think I've got it. I was still thinking of head-on collisions. Are all your ideas so difficult?

I think pictures are essential - for me anyway. I’m more than happy to provide diagrams. I may even overdo it. Please feel free to reject, mark-up, change or make diagram requests.

The following simple examples make perfect sense, including your previous text as explanation. As with my previous diagram, I show the emission/CardinalPosition (CP) directional lines only. I would need to zoom in to show the V components at each cardinal position. All particles lie in the z=0 plane. Note that my Cad system UCS Y and Z axis colors are opposite to the CP coloring choices I used here.







///////////////////////////////////////////////

Everyone, anyone, please feel free to add ideas.

///////////////////////////////////////////////

P.S. Corrected diagrams to agree with particle position change from (1,5,0) to (1,4,0)
.

Pictures are essential for me too, I just create them in my mind though. I don't have any easy-to-use tools to make them into images. The tools I do have are not vector based, so they are a little bit clumsy to use.

Are all of my ideas difficult? Probably. The ones worth mentioning are. If I thought they were simple, I'd probably just assume that they are obvious to everyone and not go into detail about it. Which is a really bad assumption because others don't always see what I see and vice-versa.

Image Request: Create an image of the receiver particle, take one cardinal point and draw the axes that are used for the linear and tangential component at that point. The axes would be rendered like the axis marker in the above images (the arrows with X, Y and Z labels) but also extended into the negative regions of those dimensions, at least for the tangential components, the linear is not so important unless the negative is considered the repulsion direction for that cardinal point (such as East in my example).

Then create a second image, starting with the above image, and add in an emitter particle and a line from the center of that emitter to the cardinal point, just like you have above.

Then create a third image, starting with the first image above, not the second, add in the emitter from the second image and instead of drawing the line from the center to the cardinal point, draw the charge vector for that line with its base at the cardinal point.

What this shows is how the linear and tangential components are part of the charge vector. You will be able to see how to split that vector along those dimensions because it will show the angles between the charge vector and those axes. You will be able to see how the charge vector points more along one dimension than another, if applicable based on their relative positions.

I think we should diagram all of this algorithm, well, as much as possible, and I might write a paper about it: How to Approximate the Charge Field of a Charged Particle. We can flesh out everything here, which we have done a lot of already, and see if it is worth it or not. Of course, we should really make sure that it works before we go too far down that path.

I have been thinking about the precision of this model and how to improve it. It isn't really the time to be doing that, since we have no idea of how this algorithm will perform, but I thought it worth mentioning that we are not limited to 6 cardinal points. Of course, they won't be cardinal if we have more (some might be, but not all) so we may have to find a new term for them. Charge Points? Whatever, it isn't important right now. I was thinking about using a dodecahedron.



Which is a spherical shape containing 20 vertices. It is easy to add even more vertices to make it even more spherical. However, with 20 vertices, we are already heading into slowness. That is 20 charge calculations per other particle per particle. This would also mean that I can't use my re-orientation trick. That only works for the cardinal points because they can be aligned with the axes. Well, I could re-orient for each Charge Point I guess. We also need to be able to determine the Charge Point that is on the opposite side from each other Charge Point. This is feasible with a dodecahedron, but not as easy as with the cardinal points.

That will become more feasible if I can find a more efficient way to implement this in a shader. I am thinking about splitting it into 2 shaders if I can. Rather than calculating the 6 (or 20) Charge Points in series (a for loop), I want to be able to calculate them in parallel (which is what shaders do) and then combine them in the next shader. It will use more memory but should run a lot faster. If it is fast enough, then we can think about using more Charge Points for more precision.

But I'm probably just muddying the waters. Don't worry too much about improvements like that. There is plenty of time for that later, once we have a working prototype.

.
The following example begins with emission direction lines but also includes projecting the charge density vectors received into their perpendicular and tangential axii at each of the cardinal positions.








.

There's too much going on in those images. Focus on a single cardinal point and don't worry about text. You can put the text into your post to explain them, but there won't be any need to put it in the image. Don't be afraid to zoom in on the cardinal point when showing something that only requires the area around that point to be in view.

The boxes around the charge vectors just hide them. I can't tell if those charge vectors are pointing in the right direction or not.

I'm not sure what those discs are trying to show, but they just get in the way. I assume it is the possible spin component directions. They might be better if they were transparent but I don't think they are needed.

Maybe make the charge vectors use a thicker line to render than the dimension markers on the cardinal points. We want the charge vectors to stand out as they are the important part.

I have just pushed some commits to GIT that have a working algorithm!

Use the test.html file to see it in action. I have created 2 charged particles that sit above and below each other so they are pole to pole. I also changed the charge density array to contain negative coefficients in the polar region (simulating attraction) and got some very interesting results.

The 2 particles are attracted to each other because they are in the polar region of the other particle. They move towards each other until they get deep enough that the charge repulsion from the equatorial region pushes them apart again. They move apart until the repulsion is gone and the attraction brings them back together again. They just bob up and down. It's so simple but it just blows me away! It means the algorithms are working as intended!

Then I took one of the particles and oriented it to emit down into the pole of the other. It pushes the other particle away but the other particles pole attracts it back, so the re-oriented particle chases the other one!

I also did the standard test of 2 particles, equator to equator, and they boringly moved apart.

You can play with it by setting the initial state of the particles. You can even add more if you want. I haven't tested that yet but the principle is the same. Go to line 119 of test.html to the initParticles function. You can see how to set the positions and there is a commented out line that shows how to set the orientation.

Have fun!

I just added a third particle (and fixed some bugs) and put them all one above the other. The middle particle gets bounced between the other 2, who are also bouncing up and down. Its a cool little dance that isn't exactly smooth as the middle particle is feeling forces from the other 2.

I then offset the middle particle in the X dimension and was amazed again as all particles pulled each other inwards (the middle particle sliding perfectly along the X axis) and then repulsed each other, causing them all to spin away.

I really wasn't expecting it to work this good without some serious tweaks. I felt the theory was sound, but not so sure that I had it all figured out in the math. I can't stop smiling! I'm really, really impressed with it. That doesn't mean it is perfect. I'm sure there are still problems and improvements to be made. The whole charge density algorithm can be replaced for one thing, but even with this fairly crude approximation it is working beautifully.

I'm sorry I can't show everyone else. You probably think I've finally gone mad, going crazy over some bouncing balls.

I added the charge field shader from AtomicViewer to show the charge emission of the particles, which helps a lot visually. Here's a screenshot of it, but the fun is in the live action.

After my excitement died down a bit, I set about cleaning up the initialization code so that it was easier to choose between various models for testing. The initParticles function now has a list of other functions (which all have a function name like init*Body*, eg initThreeBody02) that will create a ParticleArray object with various particle configurations. You choose which one you want, uncomment it, comment out any of the others and then refresh the page. I need to setup a mechanism to do that in the app, rather than in the code. At that point I could upload a copy to my site for everyone else to play with.

.
Polar aligned particles cycling their separation distance due to gravitational attraction and charge repulsion – of course! What a wonderful unexpected surprise. I shouted when the odd proton in particleArray = initThreeBody05() caused all three particles to tumble apart. Unless I’m mistaken, I don’t believe these particles are spinning yet. Dang Nevyn, your smile muscles must be hurting. Don’t worry about great charged particle discoveries like bouncing balls, I still love my marbles.





These two diagrams are better. I hope you approve, I need to get with the program.
.

They are spinning, but most of those configurations don't show it because they are aligned too well for direct collisions or have even force on both sides that cancel out. Have a look at initThreeBody05, and compare it to initThreeBody04 which is the same configuration but 04 has a Neutron in the middle and 05 has a Proton. You will see that the Neutron does not receive any spin, it just shoots off to the side (because it is offset on that side from the vertical alignment of the 2 Protons above and below). However, when there is a Proton in the middle (and offset), the charge fields interact and cause spin.

There is no gravity in this model either. The attraction I added was to represent the external ambient field pushing in from the outside, but not between the particles. This may not be needed once we implement gravity but I am also trying to represent the hole through the center of the Proton that will not resist any charge entering that area. This will need to be formalized into code too at some point.

Those images look better. I believe the charge vectors are pointing the right way. I'll give them a better look a bit later on.

That last image is showing what I want, but not really in the right way to be the most expressive. It is correct, but I think a change in the direction of the force will help to contrast the charge vectors with the cardinal point axes.

You have the receiver on the left and the emitter on the right and it will work better swapped around. Basically, I want the emitter to be negative, in all 3 dimensions, relative to the receiver. It will also help to make each dimension the same, so if the X axis has a difference of 5, then so should the Y and Z. This will cause the charge vectors to be around a 45° angle at each cardinal point, which helps to show how to break them up into linear and tangential components.

If the receiver is at (0, 0, 0), then put the emitter at (-5, -5, -5) and see how it looks. What I'm looking for is the charge vectors to be in the positive dimensions, as that is where the axis markers are pointing. Therefore, it is easier to see how those two things relate to each other.

.
Sure Nevyn, I'll be happy to redo the drawings as you've requested. If you'd like, please provide a desired set of charge density magnitudes/directions or I'll just make them up again.

///\\///\\\///\\\///\\\///\\\

Reviewing cdm-notes.txt. I wasn't able to give a proper review of the Charge Density Model last month. Now that I’ve given the subject some thought, in our discussion and the images I’ve made, here's my review.

During the charge density received calculation, the center of the receiving particle is the center of the world coordinate system - (0, 0, 0). The x, y, and z components - the linear and spin components of the charge density vector - at each cardinal position are required. All accelerations are received from some point in space, from that direction and at some magnitude. Although I may be mistaken, my thinking matches my images.

QUOTE
      Calculate sum of charge densities from other particles at the cardinal positions on this particle (N,S,E,W,F,B).
      For each position:
           Split charge vector into 3 components:
           A linear velocity component that points to the center of the particle; //points to the center of the - emitting? - particle.
UNQUOTE
Recommend changing the last sentence to –
         “A linear velocity component that may point either toward or away from the receiver center. //originating from the center of the emitter, the charge density vector may thought of as entering the cp’s on one side of the receiver and then exiting the receiver at the opposite cp".
 
Comment. Recently you’ve mentioned the idea of allowing through charge. When, say, the north or south charge density vector is well aligned with both the N/S axis direction and location, we might expect the N/S linear charge acceleration to be greatly reduced, but in those cases I would also expect some sort of trade off – like increased pair stability. It’s on the to-do list, too early for now.

Each orthogonal projection (3 at each cp) involves sines and cosines and a distance or two. How much charge received would we expect to actually penetrate to the relative opposite side of the receive particle? None really, but the cp idea is a great simplification for charge field approximations. I suppose there may be a way to simplify the charge density calculation down to just the three charge density vectors entering the receiver.

QUOTE
 “Two spin components that point in the other 2 dimensions. //point in the other 2 –orthogonal? – dimensions”.
UNQUOTE
Recommend changing that sentence to -
“Two components, both orthogonal to the linear velocity and tangential to the cardinal position. These two component directions cause particle spin".
 
QUOTE
     The actual dimensions used are different for each position:
       North: V=-Y, S=X,Z
       South: V=+Y, S=X,Z
       East:  V=-X, S=Y,Z
       West:  V=+X, S=Y,Z
       Front: V=+Z, S=X,Y
       Back:  V=-Z, S=X,Y
UNQUOTE
Comment. The table is incomplete. We know the north cp is located at (0,1,0). The charge density vector intercepting the north cp may come from any direction, with linear component not limited to -Y as you’ve indicated. Recommend changing table to.
       North: V=+/-Y, S=+/-X, +/-Z
       South: V=+/-Y, S=+/-X, +/-Z
       East:  V=+/-X, S=+/-Y, +/-Z
       West:  V=+/-X, S=+/-Y, +/-Z
       Front: V=+/-Z, S=+/-X, +/-Y
       Back:  V=+/-Z, S=+/-X, +/-Y

QUOTE
      “Where a vector with the same sign as the linear velocity dimension means a repulsion and the opposite means attraction (or do we ignore attraction and just let the other side apply repulsion?)”.
UNQUOTE
Comment. I would remove or change the sentence. As above, we will see the linear velocity components pointing both in and out of the receiver particle. I don’t believe the linear components of the charge density vector will ever indicate any attractions – only repulsions. The repulsion is felt in the directions of the linear components - in the same general direction as the six individual charge density vectors.

QUOTE
For linear velocity:

Code:

   We now have 6 vectors positioned at 6 positions.
   We reduce that to 1 vector by summing them all to find the resultant velocity.
   May need to ignore vectors that point away from the center of the particle (attraction).
   Add to the particles current velocity.
UNQUOTE
Comment. Strike out May need to ignore vectors that point away from the center of the particle (attraction). I don’t know how the algorithm is actually being performed. At present, half the linear components will point into the receiver and the other half will point out.

The spin stuff looks ok to me.

Let me know if you disagree with the above, or I'll proceed to clean it up a bit more, then change, commit and push those changes to you.
.

Nah, don't touch those notes. They were only meant to be a rough guide to remind myself of various thoughts. I will flesh them out a bit more to reflect the implementation and some realisations I have had since doing so.

Airman wrote:A linear velocity component that points to the center of the particle; //points to the center of the - emitting? - particle.

No, to the center of the receiver which is the owner of the charge (cardinal) point. You start at the charge point and look towards the center of the particle and that is the linear velocity direction for that charge point. The plane orthogonal to that direction is the spin.

Airman wrote:Each orthogonal projection (3 at each cp) involves sines and cosines and a distance or two. How much charge received would we expect to actually penetrate to the relative opposite side of the receive particle? None really, but the cp idea is a great simplification for charge field approximations. I suppose there may be a way to simplify the charge density calculation down to just the three charge density vectors entering the receiver.

Yes, I have been thinking about that. I haven't found a nice way to determine which charge points should be measured for a given emitter. It is not limited to 3 at a time though. There could be up to 5 for an emitter lined up with one of the cardinal points.

Airman wrote:The table is incomplete. We know the north cp is located at (0,1,0). The charge density vector intercepting the north cp may come from any direction, with linear component not limited to -Y as you’ve indicated

You are misunderstanding what the negative sign means. It does not mean that we ignore negative velocities. It just means that a velocity that matches the sign is considered a positive velocity. i.e. a repulsion. If the velocity (or the component of the velocity that we are interested in) is the opposite sign, then it is an attraction. It is really representing a transformation into the charge points local coordinate system.

.
We're talking a bit past one another, for example, the first quote you attribute to me in your preceding post is actually yours. I don't see a problem, we both understand the charge density vector and its cardinal position components better than your original notes. As per your instructions, I'll leave them alone, you'll clean them up at your leisure.

The following two drawings are as requested, with receiver at (0, 0, 0) and emitter at (-5, -5, -5).




.

It was really your comment after the statement that I was addressing.

The images now put the charge vectors in the right octet, but they are covering the axis markers. Maybe take away those added border lines on them. Can you make the axis markers longer, without making them thicker? About half the length of those charge vectors would be good. We need to be able to see the angles between the charge vector and axis marker.

Have you had a look over the algorithm in cdm.js? I know the code will be a bit confusing because I am using some things that aren't normally required, but if you find the 2 main functions then you can just focus on them as they provide the algorithm. The rest is just transferring data around.

.
I've looked at everything in the GIT folder at least once, of course I'm far from understanding any of it. I can't make any additions or changes till I do. I may as well concentrate on cdm.js next.


.

I have an idea about how to determine if a charge point should be measured or not. The premise is that a charge point on the far side of a receiver, from the emitter, should not receive any charge. To determine this, we look at each opposite charge point pair (the ones that are the same color in the images above).

Find the lines from each charge point to the center of the emitter.
Create a triangle by adding a third line between the 2 charge points (length will be 2 * radius).
Calculate the angle created at each charge point.
If angles are equal, then use both charge points, otherwise, use charge point with the larger angle.

The angles are calculated using the Law of Cosines which allows us to determine the angles when we only have the length of each side. If the charge points are labelled CP1 and CP2, and the center of the emitter is labelled E, the line from CP1 to E is labelled L1, from CP2 to E is L2 and from CP1 to CP2 is D (for diameter), then the angles at each charge point are given by:

cos(CP1 angle) = (L1.length * L1.length + D.length * D.length - L2.length * L2.length) / (2 * L1.length * D.length)
cos(CP2 angle) = (L2.length * L2.length + D.length * D.length - L1.length * L1.length) / (2 * L2.length * D.length)

That is an expensive operation but it is required. I can make it a little bit better by calculating both angles at the same time, and therefore saving some multiplications since we use the square of the lengths in each calculation. I can also avoid the arccos function (the cos moves to the other side of the equation and becomes an arccos), which is rudely expensive, because I don't care about the actual values, just that one angle is larger than the other. Therefore, if the right side of those equations are compared to each other, then the smaller value equates to the larger angle (arccos will invert that relationship).

I also tried to combine them into one equation by putting them into a ratio (dividing one by the other) and if the answer to that is above or below one, then we know the relationship between angles. But the resulting equation seemed to be even more complicated and expensive.

I'll give this a try tonight and see how it goes.

.
Good stuff Nevyn. I’m agog, almost overwhelmed by this wonderful learning opportunity. Thanks.

Sounds like you want diagrams with the opposite cardinal positions joined together. Will do.
 
I made a few additions to the particleArray list you provided, but I haven’t committed or pushed any of them yet. The most interesting result I found was when I started with particleArray = initThreeBody05(), then changed the top and bottom protons to neutrons.

Code:

function initThreeBody06()
 {
 var particleArray = new ChargeInteractionModel.ParticleArray( 3 );
 particleArray.items[0].object3D.position.set( 0, 5, 0 );
 particleArray.items[0].type = ChargeInteractionModel.NEUTRON;
 particleArray.items[0].recreateNode();
 particleArray.items[1].object3D.position.set( 1, 0, 0 );
 particleArray.items[2].object3D.position.set( 0, -5, 0 );
 particleArray.items[2].type = ChargeInteractionModel.NEUTRON;
 particleArray.items[2].recreateNode();
 return particleArray;
 }

From the browser, it’s plain to see the neutrons - initially stationary - begin to spin and bob. An odd spin reversal occurs at about 20sec in. Here’s a small gif - looking up at the neutron, proton and neutron, rotated to a smaller aspect, showing the apparent anomaly.

The site wont let me post the 1 Meg gif, turned it into a poor jpg image. Here's a slightly better image.


Turns out, I like Quaternions. I accept the fact that the protons are spinning. I had gotten the impression from the fact that I never saw the proton wireframes rotate.

Might there be an easy way to slow down action?
.

That is a nice little test to perform. I didn't think about that one. I was mostly focused on Protons rather than the Neutrons. I will check it out tonight.

I have noticed some strange spins that seem to suddenly kick it the opposite way. I assumed it was complex charge interactions going on because I had Protons in there. Neutrons remove that possibility, so you may have found evidence of a bug or it may be perfectly explainable. Of course, once I make the changes mentioned in my last post, things will act differently. Hopefully, it will fix it at the same time and save me the effort of tracking it down. Of course, tracking down bugs like that is how you truly understand an algorithm. You have to get right inside of it and make sure that every step is doing what you want it to do. I often find new ways to do things as I look at it closely like that.


There are a few ways to slow it down.

1. You can reduce the force applied from the charge field. This can be found in cdm.js file, find the calculateChargeDensity function and above that you will see a constant called something like MAX_CHARGE_STRENGTH = 0.2;, you can change that to a smaller value to reduce the force.

2. You can introduce a frame count in between calculations. That is, you can tell it to not calculate on every frame, which it is doing now, and it will ignore N frames before calculating again. This will need to wrap both the calculations and the updating of particles (which applies the calculated velocity and spin), otherwise it will just use the previously calculated velocity and spin during those ignored frames.

To do this you need to find the render function in test.html (near the bottom of the file) and do a little magic do stop it calculating when you don't want it to. It is a little more difficult to do than option 1, and even harder to explain. I will add that in so it will be easy to change the number of frames to ignore just by changing a single variable.

3. Run it on an old computer! OK, that is not as feasible as the other 2, but it would work.


I have created a new test function that creates a random number of particles, giving them a random position, orientation, velocity and spin, with a 25% chance of being a Neutron. You can easily change the max number of particles, the max bounds to be positioned and the chance of being a Neutron. This allows you to just refresh the page and get a new scenario. Not so good for testing the algorithm, but good for seeing how it handles lots of different interactions.

Actually, it's just cool to watch!

I have found that I can get about 200 - 300 particles and still reach 60 frames per second running on an RX580 graphics card and Intel i5 processor at around 3.5GHz (I think). The card doesn't really have much to do with it, since a pure JS implementation is doing it all in the CPU, not the GPU. Still, I am pretty happy with the performance so far. A shader implementation should be able to handle 1000 particles easily on a decent graphics card.

I haven't committed that yet, but I will tonight.

.

I have an idea about how to determine if a charge point should be measured or not. The premise is that a charge point on the far side of a receiver, from the emitter, should not receive any charge. To determine this, we look at each opposite charge point pair (the ones that are the same color in the images above).
…
I'll give this a try tonight and see how it goes.

Good luck, I guess.

I think your idea of sampling the charge received at each cardinal point is perfectly valid. Without such a simplification computations may become untenable. Sorry I made mention of charge entering or exiting, into or out of the receiver, can we make due with just three cardinal points? I didn’t intend those comments to detract from your cardinal points idea.

All the necessary computations to determine whether the receiver’s cardinal point is on the front or back side the receiver sound like a huge load. We know that the distance to the center of the receiver is the halfway point. All these lengths, any chance we can convert distances into distances^2?

I haven’t figured out how you’ve done it yet, but given the discussion, I can't wait. With respect to a possible reduced set of cardinal points, i.e. 2 equatorial cp’s and a N or S pole. Have you considered that direct charge density received at a pole will result in a gyroscopic spin axis precession?
.

You didn't detract from the idea of charge points, only showed that it is currently receiving forces that it should not. You are perfectly correct and I was already aware of it, but hadn't mentioned it. It's great that you saw it. Shows that you are understanding the algorithm.

Yes, determining which charge points to measure can get expensive. Initially, I was looking for a cheap way to do it with the most obvious being (which I mentioned in my notes) to ignore any force that has a linear velocity component (with respect to that charge point) that has a sign opposite to that charge points linear velocity component (the +'s and -'s in my table of cardinal points). I actually did implement that the other day, but I quickly realised that it would remove the attractions that I had put into the charge vectors.

In the end, we may not need those attractions, but I have a little feeling that they might still be needed. Maybe reduced from what they are now (which is equal and opposite to the equatorial charge strength at every angle from equator to pole, that is, 10° from the equator is the same, but opposite value, as 10° from the pole).

However, I am not performing this math for performance reasons. I am doing it for physical reasons, and that means that it is not negotiable. The functionality is required, but we may be able to find a cheaper way to determine the same thing. If we can remove those attractions, then I can do it very cheaply (basically just a less than or greater than comparison).

In the end, no matter how slow it is, you have to make sure that the code matches the expected outcomes. There is no point having lightning fast code if it doesn't give you the correct answer.

While nearly all scenarios will only involve 3 charge points, there is still the possibility of 5 being used at once. Put the receiver at (0, 0, 0) and the emitter at (5, 0, 0) and the receiver will get a direct hit on the East charge point. The West should not receive any charge. However, all of the others are in the same, relative, position as each other, so either they all receive charge or none of them do. This is a tricky scenario because all charge points, except E and W, are receiving charge that points away from the receiver, but still hits at a tangent. Currently, each of those points (N, S, F, B) will receive an attraction along the line from the charge point to the center of the receiver, which is just plain wrong. Maybe in this case, they should only use the East charge point. If I went with the idea of not accepting charge vectors that point away from the receiver, as these would, then it would mean that only the East point is used in this scenario. You may be right. I'll have a play with it tonight.

We need to implement gravity too. I went looking through Miles' papers the other day to find whatever I could on the strength of gravity at this level, but didn't find what I was looking for.

I also want to implement the effect that the receiver has on the charge that it receives. Is the charge reduced because it had to pass through the receivers emission? Is it ignored because it goes through the poles and straight through the central hole in the receiver?

I just had a new idea that could replace all of this! We can forget about charge points if we look at where the actual BPhoton of the receiver is in its spin cycle. This becomes the one and only charge point and we take the line back to the center of the receiver to determine which are the linear and spin components. This would involve using sines and cosines to figure out those components, but since we are only looking at a single charge point, it becomes feasible to do it that way. Of course, I would need to perform the math to determine where that BPhoton is, which is quite complicated.

That idea needs a lot of fleshing out, but is probably more accurate in the long run. Probably more complicated too. Calculating where the BPhoton is is not too bad since I already have that code in one form or another. The charge calculations might be better or they might be worse, I'm not sure at this stage.

Do you see how one thing leads to another? If I hadn't gone through this algorithm using various charge points, I might not have seen how to use it for the real BPhoton.

I have implemented the above algorithm to exclude charge points based on the angle between them and the emitter. I think it works better now. More smooth. Most configurations remain the same, but some have changed. I was shocked to see that initThreeBody04 stopped pushing the Neutron but figured out why and created a new version of that configuration called initThreeBody04_02. See if you can figure out why they behave the way they do. Very interesting configuration, the second version shows acceleration of the Neutron quite nicely.

Then I ran your configuration and my jaw hit the ground! The Neutrons were being corralled into the poles of the Proton between them. They moved about, but remained within the pole area. They gained all sorts of spins too. Amazing to watch. Loved it!

These changes have been pushed to GIT. Have a look, especially at initThreeBody04 and 04_02 but also your own configuration. Let me know if it is acting differently to how it did before these changes.

I think I just figured out what that sudden spin reversal is. It occurs when the particle is moving from within the equatorial charge zone of an emitter and into the polar charge zone or back the other way. It is the attraction I added that is causing it. The attraction should only operate on the linear velocity component and not the spin, but it is applied to both. I'm not sure how to fix it yet, either.

I fixed that, but it wasn't the problem. It was a problem, though, so it is fixed now.

I then thought that it was the quaternions overloading. When a binary number gets too large, it suddenly goes negative. If you keep adding 1 to an integer, you will eventually get back to the value that you started with. So I normalized the spin quaternions, since a quaternion can store a rotation of 370° just as well as 10°, but visually, they are the same thing. But that wasn't the problem either.

I really don't know what it is now.

I have put up a version on my website. It uses Airmans configuration with 2 Neutrons, positioned above and below a Proton. Just watch them dance around. It is quite mesmerizing.

https://www.nevyns-lab.com/mathis/app/cpim/test.html

.

Have a look, especially at initThreeBody04 and 04_02 but also your own configuration. Let me know if it is acting differently to how it did before these changes.

Gosh Nevyn, it's good to see you rolling.

Great, a pop quiz. I’m afraid I knew the answer before I read your question. I’d already reviewed your GIT changes in Sourcetree.

The main change appears to be the 3X reduction in the proton emission field density, from 3000 to 1000.

Code:

//effect.maxParticles = 3000;
effect.emissionRate = 1000;

That slowed emission reactions down 3X, which means I won't have to move back to my old computer in the spare bedroom in order to slow things down further. In “my NPN” configuration, I see the reverse spin anomaly now takes three times as long before it occurs. The reduced proton emission field strength makes the neutron positions above and below the proton appear more stable – more bobs and cycling neutron spins occur before the neutrons leave the proton.

I see the slight increase in offset positions from 1 to 1.05 as a way to compensate for the reduced emission field strength, to reduce the resulting waiting time necessary to “eject” the offset particles to the sides. The 3X reduction and 1 to 1.05 offset changes together result in what appears to be no change at all.

I thought I saw a radius increase that could suggest you’re using test.html as a geometry/sandbox while working on the alternative cardinal position methods - but that’s a stretch.

I am studying, and will continue to study cdm.js for some time to come. I see your latest changes to the charge density received calculation, reducing the number of cardinal position used. I'm not at all certain I see any "smoothing", less calculations has got to be good.

What am I missing? I have several staged/uncommitted file changes you made to cdm-notes.txt and cdm.js txt. I thought they would update at the last pull, I don’t know how to incorporate or reconcile them at present. Everything else is just wonderful.
.

LongtimeAirman wrote:The main change appears to be the 3X reduction in the proton emission field density, from 3000 to 1000.

Code:

//effect.maxParticles = 3000;
effect.emissionRate = 1000;

Nope, that's not it. That only effects the charge shader I imported from AV. It is only a visual representation of the charge emission and does not play any part in the calculations. But I don't blame you for thinking that caused the slowdown.

The slowness is caused by what I mentioned yesterday. I reduced the value of MAX_CHARGE_STRENGTH from 0.2 to 0.1.

LongtimeAirman wrote:I see the slight increase in offset positions from 1 to 1.05 as a way to compensate for the reduced emission field strength, to reduce the resulting waiting time necessary to “eject” the offset particles to the sides. The 3X reduction and 1 to 1.05 offset changes together result in what appears to be no change at all.

Not quite. The initThreeBody04 function puts the Neutron at X=1 while the Protons above and below it are at X=0 (they have a different Y value, obviously). The radius of a Proton or Neutron is 1, so the Protons were emitting straight up and down at the West charge point on the Neutron. The North and South charge points also received charge. This means that all charge received by that Neutron is balanced. The West charge point only receives spin since it receives purely tangential vectors. The North and South charge points receive both tangential and linear velocity, but again, it is equal and opposite so no net force. The Neutron does not move.

In initThreeBody04_02, I moved the Neutron over to X=1.05 so that the West charge point would receive a little bit of linear velocity. Because I used such a small increase, you can see how the Protons give it a little bit of force the first time they move towards the Neutron but less as they move back up. The second time they come down they give it another push and the Neutron shoots away at a quicker velocity.

LongtimeAirman wrote:I thought I saw a radius increase that could suggest you’re using test.html as a geometry/sandbox while working on the alternative cardinal position methods - but that’s a stretch.

Not sure what you saw. Maybe the line that looks like var D = RADIUS * 2;, but that is not actually changing the radius, only calculating the diameter to create the third side of the triangle for the charge point exclusion code.

LongtimeAirman wrote:What am I missing? I have several staged/uncommitted file changes you made to cdm-notes.txt and cdm.js txt. I thought they would update at the last pull, I don’t know how to incorporate or reconcile them at present. Everything else is just wonderful.

If you are prepared to lose all of your own changes, then you can do a hard reset to the top commit of mine. You do this by right-clicking on the commit you want to reset to and click on the 'Reset to this commit' menu item. It will display a dialog box with some options and one of them is a dropdown menu with options like Soft, Medium and Hard. A hard reset will discard all of your own changes and the source will reflect the GIT structure and content.

Alternatively, in the File Status tab or Log/History tab, click on the 'Unsaved changes' line at the top of the commit tree and the box at the bottom will show all of the changed files. You can right-click on each one and select 'Discard' to get rid of your changes to that file. You can probably select multiple files and discard them all at once.

If you want to keep your changes but they aren't merging, then it gets more complicated. Let me know if this is the path you need to take.

.
I was so sure I had the answer; but no, wrong again. I’ll get you yet Red Barron!  

But then, it may take a lot longer than expected, as of a few hours ago, it appears you’ve made at least eight commits today. Good gosh Nevyn. I’ll try to give them a proper review tomorrow.  
.

Well, I was separating various little changes into individual commits just in case we wanted to take some out later. Each commit is an 'All or Nothing' type of thing. It also helps you to see the related changes for a given update. If you read the comment against each commit, you should be able to tell which ones relate to which problem/solution.

.
It’s difficult making my slow progress in this project interesting to the average reader – even charge field devotees. Don’t feel sorry for me, I think its fun using this fancy new platform.

Thanks for including initRandomBody01 - making initial configurations of large numbers of random particles - Nevyn. Like everything else you do, it’s got more puzzling details than I can currently count, it’s also given me a few new ideas.

Today I added a basic lattice configuration. The lattice is 4 particle layers thick. The top and bottom layers are protons, and the two middle layers are neutrons. I posted two versions: 1. initLatticeBody01. Proton separations (15) greater than the emission radius (10) let this lattice just sit there – like a brick.



And 2. initLatticeBody02. Proton separations are slightly less than the emission radius – the brick blows apart.



Each proton is emitting charge. The breaking lattice makes a larger computational load, see the slower 29 FPS (frames per second) here, vice 55 FPS for “the brick”. Note that the two larger neutrons below and to the right of the image center are overlapping – passing through one another. Collisions seem so far away right now.

I’ve got plenty of little questions: like how did the proton appear in the center of the brick? I tried incorporating little random changes to the neutron’s y positions with no luck. Several comments indicate notional tasks, not yet operational, and not worthy of being called bugs.
 
Still, baby steps. My immediate goal – nothing is certain in this life – is to launch a modified initOneBody01_02 at the wall. That should be fun.

Nevyn, you’ve mentioned needing to implement gravity a time or two, any chance you might be able to incorporate expansion gravity?
.

I haven't seen the code yet, so I don't know why that proton sits in the middle. You may be creating more particles than you need. Remember that array indices are 0-based, so the last item in the array is Array.length - 1. If you created a ParticleArray with 3 particles, then you reference them as 0, 1 and 2. Not 1, 2, 3. Check the line that looks like:

Code:


var particleArray = new ChargeInteractionModel.ParticleArray( <number of items> );


and the following code that will be setting the positions of those particles. Make sure that the highest item index that you use is 1 less than the above. My guess is that you are not changing one of the items, so it remains in the center. It could also be the first item that you are missing, which should be addressed as 0. You can also make sure that you are addressing all items in the array, and have not missed one. That is, all of the array indices that you use are sequential.

There is no particle collision code yet. I have a method in there to handle it, and the algorithm checks for it and calls that method, but it does nothing (apart from setting the return vector to 0, 0, 0, which it must do). Maybe you would like to implement that method? Just a standard inelastic collision. Use a mass of 1 (which means it can probably be ignored in the equation, but don't, as we will set mass and radius soon enough). Find the RADIUS declaration and add in another one for MASS. That makes it easier to change it, but eventually I want that as passed in data. Then we can deal with electrons as well. Remember that the collision method is called for 1 particle only. Not both particles in the collision. The other particle will call the method again for itself. It is inefficient, but I'm not too concerned about that at the moment.

However, the way the collision method is called will mean that it will not work for neutrons. I will need to change the driver algorithm to handle neutrons better. It was written with charged particles in mind and I added in neutrons later and haven't revisited that section.

It is quite amazing how much action we get with just the charge field alone. Kind of shows how it is more important than gravity at this scale.

I'm not sure how to implement gravity yet. I need to read through a few of Miles' papers again to get a feel for how to add it in along side of the charge field code, rather than it including that charge field stuff itself, as the mainstream math does. If I can do it via expansion, without anything actually expanding, then I will do it that way. That is, I don't want the radii of the particles to change at all, but I may be able to calculate how they would change and then use that to determine the gravity vector.

What I need to know is how to balance the charge field with gravity. What is their relative strengths at this scale? As long as we stick to that ratio, then we can make the values whatever we want.

I think a more serious change is required to implement the collision method. We need to know the current velocity to calculate the collision and we don't have that available at the moment. You can still start implementing it and any data that you need, but don't have, then just add a new method parameter for it and I will make sure that it is available and passed in.

.
I was innocent. My picture showed a new initial configuration with a breakdown resulting in converging and overlapping neutrons; a point worth mentioning. I don’t believe anyone expected any collisions yet, although there’s no denying you could have easily done so. At least the neutrons didn’t merge like the particles in Proto-planet.

I repeat, I’m happy for this opportunity to work with you – helping build a charge field collision model. You're the boss. I have plenty of appreciation for your charge field knowledge and math, software and engineering skills. Believe me, I’m anxious to contribute, I'm devoting time to it every day.

For maximum efficiency, I would assume collision detection should be the first task performed in the “charge density received” calculation since a collision can override the calculation. It also makes sense to calculate the accelerations due to gravity alongside charge in the same “charge density received” calculation – 1/r^2 for gravity and r^4 for charge. If you feel we must establish some scale of radius/tangential velocity, then the proton scale is fine. Macro sizes take longer, but I don't believe our model can do the difference justice. We certainly don’t need a gravity calculation for every cardinal position, just a single calculation, from this.particle’s position to the other particle’s position (center to center).

I agree, from what I understand, "a more serious change is required to implement the collision method". Your current collision constraints seem unworkable. Any collision requires being able to reset both particles’ motion set: position, orientation, velocity and spin.

For us, a “simple” collision might mean that we can allow this.particle’s position to remain the same; however, this.particle’s orientation should be reset. I will argue that the “other particle” colliding with this.particle must at least be pushed back to the position at which contact would occur, with a new rebound velocity – based on it’s previous velocity. Even better, we should include the appropriate rebound distance.
.