In classical mechanics, motion of an object is described by vector variables of position, velocity, acceleration and force and a scalar variable of mass. A velocity is defined as a derivative of the position of an object with respect to time. An acceleration is defined as a derivative of the velocity with respect to time and which is a second degree derivative of the position.

Then a force, a mass and an acceleration are related by this equation of Newton's laws of motion, and therefore a velocity and a position are also related through derivation with respect to time.

To implement a particle with motion behavior following Newton's laws, you'd start with putting fields of a position, a velocity, an acceleration, a force and a mass in an agent class. The type of variable of mass is scalar (double) and others are vectors (IVec). Because a mass of an object usually doesn't change, the variable of a mass can be seen as a constant value and is initialized at the declaration of the variable. Usually in simulation of motion, an object is initialized with an initial position and velocity, and the constructor of the class takes two input vectors to initialize the instance's position and velocity.

import processing.opengl.*;
import igeo.*;

class MyParticle extends IAgent{
  IVec pos;
  IVec vel;
  IVec acc;
  IVec frc;
  double mass = 1.0; // default mass
  
  MyParticle(IVec p, IVec v){
    pos = p;
    vel = v;
  }
}

Whereas a position and a velocity are internal states of an object, a force is something which can be applied to the object from the outside of the object. In each time frame, a force is applied to the object and updates its velocity and position. An object could receive multiple vector forces and those forces are calculated into a resultant force vector as described in the following formula.

A resultant vector is shown graphically below and it's just a summation vector of all received vectors in one time frame.

To have the mechanism to receive multiple forces, from the outside of the class or the inside, push() method is defined as the following.

import processing.opengl.*;
import igeo.*;

class MyParticle extends IAgent{
  IVec pos, vel, acc, frc;
  double mass = 1.0; // default mass
  
  MyParticle(IVec p, IVec v){
    pos = p;
    vel = v;
    frc = new IVec();
  }

  void push(IVec f){
    frc.add(f);
  }
}

To update an object's acceleration, velocity and position, you only need to use the resultant force vector, not each individual force vector, because the equations in Newton's laws are all linear-independent. According to the Newton's laws and the definition of an acceleration and a velocity, the agent's acceleration, velocity and position can be calculated as the following.

To implement those formulas shown above, the code of update() method can be written as the following.

import processing.opengl.*;
import igeo.*;

class MyParticle extends IAgent{
  final double dt = IConfig.updateRate; //second per frame
  IVec pos, vel, acc, frc;
  double mass = 1.0; // default mass
  
  MyParticle(IVec p, IVec v){
    pos = p;
    vel = v;
    frc = new IVec();
  }

  void push(IVec f){
    frc.add(f);
  }

  void update(){
    acc = frc.div(mass);
    vel.add(acc.mul(dt));
    pos.add(vel.cp().mul(dt)); //vel itself can't be changed
    frc.zero(); //reset frc by setting zero
  }
}

The variable dt is a time in second between each time frame to update agents. It has the keyword final to use it as a constant variable and the value of the time between frames of the iGeo system is stored in IConfig.updateRate. This value is set by the method IG.updateRate(double) or its alias method IG.rate(double) or method IG.speed(double). In the update method, the acceleration (acc) is defined by dividing the force (frc) by the mass. Then the velocity (vel) is updated by adding the acceleration, which is an increment (or decrement) of velocity per unit time, multiplied by the interval time of frames (dt). In the similar way, the position (pos) is updated with the velocity, which is an increment of position per unit time, but because the value of the velocity needs to be preserved (on the other hand, the acceleration and the force don't need to be preserved because they are only for the duration of each time frame), the content of vel is copied with cp() method before multiplying dt. And then the force (frc) is reset to zero for the next time frame.

To use this particle class, you'd add the code to apply some force to it, and also to create some geometries. The following example creates one particle instance having a force to push towards negative y-direction and creating a point geometry on each time frame.

import processing.opengl.*;
import igeo.*;

void setup(){
  size(480, 360, IG.GL);
  IG.duration(200);
  new MyParticle(IG.v(0,0,0), IG.v(10,20,0));
}

class MyParticle extends IAgent{
  final double dt = IConfig.updateRate;//second per frame
  IVec pos, vel, acc, frc;
  double mass = 1.0; // default mass
  
  MyParticle(IVec p, IVec v){
    pos = p;
    vel = v;
    frc = new IVec();
  }

  void push(IVec f){
    frc.add(f);
  }

  void update(){
    push(IG.v(0,-10,0));
    new IPoint(pos.cp());
    
    acc = frc.div(mass);
    vel.add(acc.mul(dt));
    pos.add(vel.cp().mul(dt)); //vel itself can't be changed
    frc.zero(); //reset frc by setting zero
  }
}

To have a clearer organization of the code, the force application part and the geometry creation part can be separated from the base particle class as its child class. The following code shows MyParticleBase class as the base class and the MyParticle class as the separated child class inheriting MyParticleBase.

import processing.opengl.*;
import igeo.*;

void setup(){
  size(480, 360, IG.GL);
  IG.duration(200);
  new MyParticle(IG.v(0,0,0), IG.v(10,20,0));
}

class MyParticle extends MyParticleBase{
  MyParticle(IVec p, IVec v){ super(p,v); }

  void update(){
    push(IG.v(0,-10,0));
    new IPoint(pos.cp());
    super.update();
  }
}

class MyParticleBase extends IAgent{
  final double dt = IConfig.updateRate;//second per frame
  IVec pos, vel, acc, frc;
  double mass = 1.0; // default mass
  
  MyParticleBase(IVec p, IVec v){
    pos = p;
    vel = v;
    frc = new IVec();
  }

  void push(IVec f){ frc.add(f); }

  void update(){
    acc = frc.div(mass);
    vel.add(acc.mul(dt));
    pos.add(vel.cp().mul(dt)); //vel itself can't be changed
    frc.zero(); //reset frc by setting zero
  }
}

iGeo has a built-in particle class named IParticle which implements the same algorithm with the above. Using IParticle class, the code above can be re-written as the following replacing MyParticleBase with IParticle. Then you also need to replace pos with pos(), vel with vel(), acc with acc(), frc with frc() and mass with mass() when you access to them.

import processing.opengl.*;
import igeo.*;

void setup(){
  size(480, 360, IG.GL);
  IG.duration(200);
  new MyParticle(IG.v(0,0,0), IG.v(10,20,0));
}

class MyParticle extends IParticle{
  MyParticle(IVec p, IVec v){ super(p,v); }

  void update(){
    push(IG.v(0,-10,0));
    new IPoint(pos().cp());
  }
}

When you inherit IParticle class, you don't need to have a line of
super.update();
in the update method.

Here is another example of use of IParticle to define a random walk behavior by adding a random vector as force.

import processing.opengl.*;
import igeo.*;

void setup(){
  size(480, 360, IG.GL);
  new MyParticle(IG.v(0,0,0), IG.v(0,0,0));
  new MyParticle(IG.v(0,0,0), IG.v(0,0,0));
  new MyParticle(IG.v(0,0,0), IG.v(0,0,0));
}

class MyParticle extends IParticle{
  MyParticle(IVec p, IVec v){ super(p,v); }

  void update(){
    IVec force = IRand.pt(-100,-100,0,100,100,0);
    push(force);
    new IPoint(pos().cp());
  }
}

Yet another IParticle example adding force based on the particle's velocity.

import processing.opengl.*;
import igeo.*;

void setup(){
  size(480, 360, IG.GL);
  new MyParticle(IG.v(0,0,0), IG.v(10,0,0));
}

class MyParticle extends IParticle{
  MyParticle(IVec p, IVec v){ super(p,v); }

  void update(){
    IVec force = vel().cp().mul(2).rot(PI/2);
    push(force);
    new IPoint(pos().cp());
  }
}

The following code shows an example of interaction between the same particle class.

import processing.opengl.*;
import igeo.*;

void setup(){
  size(480, 360, IG.GL);
  for(int i=0; i < 10; i++){
    new MyParticle(IG.v(0,i*10,0), IG.v(10,0,0)).clr(0.1*i,0,0);
  }
}

class MyParticle extends IParticle{
  MyParticle(IVec p, IVec v){ super(p,v); }

  void interact(ArrayList < IDynamics > agents){
    for(int i=0; i < agents.size(); i++){
      if(agents.get(i) instanceof MyParticle){
        MyParticle p = (MyParticle)agents.get(i);
        if(p != this){
          if(p.dist(this) < 20){ //distance threshold
            IVec force = p.pos().dif(pos());
            force.mul(0.01);
            push(force);
          }
        }
      }
    }
  }
  
  void update(){
    new IPoint(pos().cp()).clr(clr());
  }
}

A particle of this particle class checks other particles if they are within a certain threshold. If so, it calculates a difference vector from the particle to another particle position as the force direction

IVec force = p.pos().dif(pos());

and adjust the intensity of the force by multiplying some number.

force.mul(0.01);

Then it push itself by feeding the force vector to the push method.

push(force);

As result, particles attract each other and gather closer.

If you put the rotating behavior in update method as shown in the previous section, keeping the attraction behavior in the interact method, you can see how the interaction behavior intervenes the rotation behavior as shown below.

import processing.opengl.*;
import igeo.*;

void setup(){
  size(480, 360, IG.GL);
  for(int i=0; i < 10; i++){
    new MyParticle(IRand.pt(0,0,100,100), IG.v(10,0,0)).clr(IRand.clr());
  }
}

class MyParticle extends IParticle{
  MyParticle(IVec p, IVec v){ super(p,v); }

  void interact(ArrayList < IDynamics > agents){
    for(int i=0; i < agents.size(); i++){
      if(agents.get(i) instanceof MyParticle){
        MyParticle p = (MyParticle)agents.get(i);
        if(p != this){
          if(p.dist(this) < 20){
            IVec force = p.pos().dif(pos());
            force.mul(1.2);
            push(force);
          }
        }
      }
    }
  }
  
  void update(){
    IVec force = vel().cp().mul(2).rot(PI/2);
    push(force);
    new IPoint(pos().cp()).clr(clr());
  }
}

The following code combines the random walk algorithm by adding a random vector as a force in the update method and the attraction algorithm in the interact method.

import processing.opengl.*;
import igeo.*;

void setup(){
  size(480, 360, IG.GL);
  for(int i=0; i < 30; i++){
    new MyParticle(IRand.pt(0,0,100,100), IG.v(10,0,0)).clr(IRand.clr());
  }
}

class MyParticle extends IParticle{
  MyParticle(IVec p, IVec v){ super(p,v); }

  void interact(ArrayList < IDynamics > agents){
    for(int i=0; i < agents.size(); i++){
      if(agents.get(i) instanceof MyParticle){
        MyParticle p = (MyParticle)agents.get(i);
        if(p != this){
          if(p.dist(this) < 20){
            IVec force = p.pos().dif(pos());
            force.mul(2.0);
            push(force);
          }
        }
      }
    }
  }
  
  void update(){
    IVec force = IRand.pt(-100,-100,0,100,100,0);
    push(force);
    new IPoint(pos().cp()).clr(clr());
  }
}

The particle agent above was applying a force to itself inside update() method, but if you think the force as gravity or some type of an external force, it's more natural to separate the code to apply the force outside of the particle class. You can do it by defining another agent which interacts with the particle agents through interact() method. In the following code, the gravity force is implemented as MyGravity class inheriting IAgent.

import processing.opengl.*;
import igeo.*;

void setup(){
  size(480, 360, IG.GL);
  IG.duration(200);
  new MyParticle(IG.v(0,0,0), IG.v(10,20,0));
  new MyGravity();
}

class MyParticle extends IParticle{
  MyParticle(IVec p, IVec v){ super(p,v); }
  void update(){ 
    new IPoint(pos().cp()); 
  }
}

class MyGravity extends IAgent{

  void interact(IDynamics agent){
    if(agent instanceof MyParticle){
      MyParticle particle = (MyParticle)agent;
      particle.push(IG.v(0,-10,0));
    }
  }
}

The code above use simpler definition of interact method
void interact(IDynamics agent)
rather than the longer and faster definition of interact method
void interact(ArrayList< IDynamics > agents)
to keep the code simple. The if-condition on the first line if( agent instanceof MyParticle ) in the interact method checks if the the current checking agent is an instance of MyParticle class. Then if so, it's casting the variable agent which is an instance of the parent class of IAgent, IDynamics, into a type of MyParticle. Then finally it apply a force to the particle agent via push() method. For the detail of interact() method, please see this tutorial page about interact method.

The following code shows a small revision of the previous code. For better reusability of code, it's better to bring a variable into MyGravity class to control the force vector. The vector variable gravity is defined and initialized at the constructor of the class. Actual direction and amount of the gravity force is given at setup() method.

import processing.opengl.*;
import igeo.*;

void setup(){
  size(480, 360, IG.GL);
  IG.duration(200);
  new MyParticle(IG.v(0,0,0), IG.v(10,20,0));
  new MyGravity(IG.v(0,-10,0));
}

class MyParticle extends IParticle{
  MyParticle(IVec p, IVec v){ super(p,v); }
  void update(){ new IPoint(pos().cp()); }
}

class MyGravity extends IAgent{
  IVec gravity;
  
  MyGravity(IVec g){ gravity=g; }

  void interact(IDynamics agent){
    if(agent instanceof MyParticle){
      MyParticle particle = (MyParticle)agent;
      particle.push(gravity);
    }
  }
}

The example code below shows multiple particles with different initial positions and velocities reacting to one gravity agent.

import processing.opengl.*;
import igeo.*;

void setup(){
  size(480, 360, IG.GL);
  IG.duration(200);
  int num = 40;
  for(int i=0; i < num; i++){
    double a = i*2*PI/num;
    new MyParticle(IG.v(20,0,0).rot(a),IG.v(10,0,i).rot(a));
  }
  new MyGravity(IG.v(0,0,-10));
}

class MyParticle extends IParticle{
  MyParticle(IVec p, IVec v){ super(p,v); }
  void update(){ new IPoint(pos().cp()); }
}

class MyGravity extends IAgent{
  IVec gravity;
  
  MyGravity(IVec g){ gravity=g; }

  void interact(IDynamics agent){
    if(agent instanceof MyParticle){
      MyParticle particle = (MyParticle)agent;
      particle.push(gravity);
    }
  }
}

As one of interacting agent with the particle agent, you can define bouncing behavior on a virtual xy-plane. The code below has a class BouncePlane and to create this bounce agent, you provide z height at the constructor. In the interact method, first it checks if the particle is below the z level and if so, it push the particle up to the z level and the flip the velocity of the particle by ref() method on this line.
particle.vel().ref(IG.zaxis);

import processing.opengl.*;
import igeo.*;

void setup(){
  size(480, 360, IG.GL);
  IG.duration(200);
  int num = 40;
  for(int i=0; i < num; i++){
    double a = i*2*PI/num;
    new MyParticle(IG.v(20,0,0).rot(a), IG.v(10,0,i).rot(a));
  }
  new MyGravity(IG.v(0,0,-10));
  new BouncePlane(-5);
}

class MyParticle extends IParticle{
  MyParticle(IVec p, IVec v){ super(p,v); }
  void update(){ new IPoint(pos().cp()); }
}

class MyGravity extends IAgent{
  IVec gravity;
  MyGravity(IVec g){ gravity=g; }
  void interact(IDynamics agent){
    if(agent instanceof MyParticle){
      MyParticle particle = (MyParticle)agent;
      particle.push(gravity);
    }
  }
}

class BouncePlane extends IAgent{
  double z;
  
  BouncePlane(double zval){ z=zval; }

  void interact(IDynamics agent){
    if(agent instanceof MyParticle){
      MyParticle particle = (MyParticle)agent;
      if(particle.pos().z() < z){
        particle.pos().add(0,0,z - particle.pos().z());
        particle.vel().ref(IG.zaxis);
      }
    }
  }
}

To create geometries on a trace of particle agents, one way is to create geometries inside update method. The previous two codes above create a point object inside the update method. If you want to create a line between the position of the particle and the previous position, first you need to let the agent remember the previous position by defining an instance field of prevPos. Then you put the copy of the current position into prevPos at the end of update method, because the content of pos() is constantly changing. In the very first execution of update method of an agent, there is no previous position and the content of prevPos is null. Because of this, the if-condition if(prevPos != null){ is used to skip the very first execution to create a line.

import processing.opengl.*;
import igeo.*;

void setup(){
  size(480, 360, IG.GL);
  IG.duration(200);
  int num = 40;
  for(int i=0; i < num; i++){
    double a = i*2*PI/num;
    new MyParticle(IG.v(20,0,0).rot(a), IG.v(10,0,i).rot(a));
  }
  new MyGravity(IG.v(0,0,-10));
  new BouncePlane(-5);
}

class MyParticle extends IParticle{
  IVec prevPos;
  MyParticle(IVec p, IVec v){ super(p,v); }
  void update(){
    IVec curPos = pos().cp();
    if(prevPos != null){
      new ICurve(prevPos, curPos).clr(0);
    }
    prevPos = curPos;
  }
}

class MyGravity extends IAgent{
  IVec gravity;
  MyGravity(IVec g){ gravity=g; }
  void interact(IDynamics agent){
    if(agent instanceof MyParticle){
      MyParticle particle = (MyParticle)agent;
      particle.push(gravity);
    }
  }
}

class BouncePlane extends IAgent{
  double z;
  
  BouncePlane(double zval){ z=zval; }

  void interact(IDynamics agent){
    if(agent instanceof MyParticle){
      MyParticle particle = (MyParticle)agent;
      if(particle.pos().z() < z){
        particle.pos().add(0,0,z - particle.pos().z());
        particle.vel().ref(IG.zaxis);
      }
    }
  }
}

You might notice there is a point geometry at the latest position of the particle agent although there is no code to create a point geometry inside MyParticle class. This is a point of IParticle class to show the current position of the particle. If you don't want to show it, you can hide by the method of IParticle hide(). The example to hide the point is shown in the next code.

There is another way to create geometries on a trace of particles. It is to create another independent agent class just to create geometries in every time frame. The next code has the class MyGeometry class. This class has two instance fields of MyParticle, particle1 and particle2 and it creates a surface in every time frame inside the update method. The instances of MyGeometry class are initialized in the setup method by feeding adjacent two particle agents.

import processing.opengl.*;
import igeo.*;

void setup(){
  size(480, 360, IG.GL);
  IG.duration(200);
  int num = 40;
  MyParticle[] particle = new MyParticle[num];
  for(int i=0; i < num; i++){
    double a = i*2*PI/num;
    particle[i] =
    new MyParticle(IG.v(20,0,0).rot(a), IG.v(10,0,i).rot(a));
  }
  for(int i=1; i < num; i++){
    new GeometryAgent(particle[i-1], particle[i]);
  }
  new MyGravity(IG.v(0,0,-10));
  new BouncePlane(-5);
}

class MyParticle extends IParticle{
  IVec prevPos;

  MyParticle(IVec p, IVec v){ 
    super(p,v);
    hide(); // hiding the point of IParticle
  }
  void update(){
    IVec curPos = pos().cp();
    if(prevPos != null){
      new ICurve(prevPos, curPos);
    }
    prevPos = curPos;
  }
}

class MyGravity extends IAgent{
  IVec gravity;
  MyGravity(IVec g){ gravity=g; }
  void interact(IDynamics agent){
    if(agent instanceof MyParticle){
      MyParticle particle = (MyParticle)agent;
      particle.push(gravity);
    }
  }
}

class BouncePlane extends IAgent{
  double z;
  
  BouncePlane(double zval){ z=zval; }

  void interact(IDynamics agent){
    if(agent instanceof MyParticle){
      MyParticle particle = (MyParticle)agent;
      if(particle.pos().z() < z){
        particle.pos().add(0,0,z - particle.pos().z());
        particle.vel().ref(IG.zaxis);
      }
    }
  }
}

class GeometryAgent extends IAgent{
  MyParticle particle1, particle2;
  IVec prevPos1, prevPos2;

  GeometryAgent(MyParticle p1, MyParticle p2){
    particle1 = p1;
    particle2 = p2;
  }

  void update(){
    IVec curPos1 = particle1.pos().cp();
    IVec curPos2 = particle2.pos().cp();
    if(prevPos1 != null && prevPos2 != null){
      new ISurface(prevPos1, curPos1, curPos2, prevPos2).clr(clr());
    }
    prevPos1 = curPos1;
    prevPos2 = curPos2;
    // update color
    double r = red() + IRand.get(-0.01,0.01);
    double g = green() + IRand.get(-0.02,0.02);
    clr(r, g, 1.0);
  }
}