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.

add_library('igeo')

class MyParticle(IAgent) : 
    
    def __init__(self, p, v) : 
        self.pos = p
        self.vel = v
        self.mass = 1.0

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.

add_library('igeo')

class MyParticle(IAgent) : 
    
    def __init__(self, p, v) : 
        self.pos = p
        self.vel = v
        self.mass = 1.0
        self.frc = IVec()
    
    def push(self, f) : 
        self.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.

add_library('igeo')

class MyParticle(IAgent) : 
    dt = IConfig.updateRate; #static variable,  second per frame
    
    def __init__(self, p, v) : 
        self.pos = p
        self.vel = v
        self.mass = 1.0
        self.frc = IVec()

    def push(self, f) : 
        self.frc.add(f)
    
    def update(self) : 
        acc = self.frc.div(self.mass)
        self.vel.add(acc.mul(MyParticle.dt))
        self.pos.add(self.vel.cp().mul(MyParticle.dt)); #vel itself can't be changed
        self.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.

add_library('igeo')

def setup() : 
    size(480, 360, IG.GL)
    IG.duration(200)
    MyParticle(IG.v(0,0,0), IG.v(10,20,0))

class MyParticle(IAgent) : 
    dt = IConfig.updateRate#static variable, second per frame
    
    def __init__(self, p, v) : 
        self.pos = p
        self.vel = v
        self.mass = 1.0
        self.frc = IVec()
    
    def push(self, f) : 
        self.frc.add(f)
    
    def update(self) : 
        self.push(IG.v(0,-10,0))
        IPoint(self.pos.cp())
        
        acc = self.frc.div(self.mass)
        self.vel.add(acc.mul(MyParticle.dt))
        self.pos.add(self.vel.cp().mul(MyParticle.dt)); #vel itself can't be changed
        self.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.

add_library('igeo')

def setup() : 
    size(480, 360, IG.GL)
    IG.duration(200)
    MyParticle(IG.v(0,0,0), IG.v(10,20,0))

class MyParticleBase(IAgent) : 
    dt = IConfig.updateRate#static variable, second per frame
    
    def __init__(self, p, v) : 
        self.pos = p
        self.vel = v
        self.mass = 1.0
        self.frc = IVec()
    
    def push(self, f) : 
        self.frc.add(f)
    
    def update(self) : 
        acc = self.frc.div(self.mass)
        self.vel.add(acc.mul(MyParticleBase.dt))
        self.pos.add(self.vel.cp().mul(MyParticleBase.dt)); #vel itself can't be changed
        self.frc.zero() #reset frc by setting zero

class MyParticle(MyParticleBase) : 
    def __init__(self, p, v) :
        MyParticleBase.__init__(self,p,v)
    
    def update(self) :
        self.push(IG.v(0,-10,0))
        IPoint(self.pos.cp())
        MyParticleBase.update(self)

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.

add_library('igeo')

def setup() : 
    size(480, 360, IG.GL)
    IG.duration(200)
    MyParticle(IG.v(0,0,0), IG.v(10,20,0))

class MyParticle(IParticle) : 
    def __init__(self, p, v) : 
        IParticle.__init__(self,p,v)
    
    def update(self) : 
        self.push(IG.v(0,-10,0))
        IPoint(self.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.

add_library('igeo')

def setup() : 
    size(480, 360, IG.GL)
    MyParticle(IG.v(0,0,0), IG.v(0,0,0))
    MyParticle(IG.v(0,0,0), IG.v(0,0,0))
    MyParticle(IG.v(0,0,0), IG.v(0,0,0))


class MyParticle(IParticle) : 
    def __init__(self, p, v) :
        IParticle.__init__(self,p,v)
    
    def update(self) : 
        force = IRand.pt(-100,-100,0,100,100,0)
        self.push(force)
        IPoint(self.pos().cp())

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

add_library('igeo')

def setup() : 
    size(480, 360, IG.GL)
    MyParticle(IG.v(0,0,0), IG.v(10,0,0))

class MyParticle(IParticle) : 
    def __init__(self, p, v) : 
        IParticle.__init__(self,p,v)
    
    def update(self) : 
        force = self.vel().cp().mul(2).rot(PI/2)
        self.push(force);
        IPoint(self.pos().cp())

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

add_library('igeo')

def setup() : 
    size(480, 360, IG.GL)
    for i in range(10) : 
        MyParticle(IG.v(0,i*10,0), IG.v(10,0,0)).clr(0.1*i,0,0)

class MyParticle(IParticle) : 
    def __init__(self, p, v) : 
        IParticle.__init__(self,p,v)
    
    def interact(self, agents) : 
        for agent in agents : 
            if isinstance(agent, MyParticle) : 
                if agent is not self : 
                    #distance threshold
                    if agent.dist(self) < 20 : 
                        force = agent.pos().dif(self.pos())
                        force.mul(0.01)
                        self.push(force)
    
    def update(self) : 
        IPoint(self.pos().cp()).clr(self.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.

add_library('igeo')

def setup() : 
    size(480, 360, IG.GL)
    for i in range(10) : 
        MyParticle(IRand.pt(0,0,100,100), IG.v(10,0,0)).clr(IRand.clr())

class MyParticle(IParticle) : 
    def __init__(self, p, v) :
        IParticle.__init__(self,p,v)
    
    def interact(self, agents) : 
        for agent in agents : 
            if isinstance(agent, MyParticle) : 
                if agent is not self :
                    if agent.dist(self) < 20 :
                        force = agent.pos().dif(self.pos())
                        force.mul(1.2)
                        self.push(force)
    
    def update(self) :
        force = self.vel().cp().mul(2).rot(PI/2)
        self.push(force)
        IPoint(self.pos().cp()).clr(self.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.

add_library('igeo')

def setup() : 
    size(480, 360, IG.GL)
    for i in range(30) : 
        MyParticle(IRand.pt(0,0,100,100), IG.v(10,0,0)).clr(IRand.clr())

class MyParticle(IParticle) : 
    def __init__(self, p, v) : 
        IParticle.__init__(self,p,v)
    
    def interact(self, agents) : 
        for agent in agents : 
            if isinstance(agent, MyParticle) : 
                if agent is not self : 
                    if(agent.dist(self) < 20) : 
                        force = agent.pos().dif(self.pos())
                        force.mul(2.0)
                        self.push(force)
    
    def update(self) : 
        force = IRand.pt(-100,-100,0,100,100,0)
        self.push(force)
        IPoint(self.pos().cp()).clr(self.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.

add_library('igeo')

def setup() : 
    size(480, 360, IG.GL)
    IG.duration(200)
    MyParticle(IG.v(0,0,0), IG.v(10,20,0))
    MyGravity()

class MyParticle(IParticle) : 
    def __init__(self, p, v) :
        IParticle.__init__(self,p,v)
    
    def update(self) :
        IPoint(self.pos().cp())

class MyGravity(IAgent) : 
    
    def interact(self, agents) : 
        for agent in agents : 
            if isinstance(agent, MyParticle) : 
                agent.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.

add_library('igeo')

def setup() : 
    size(480, 360, IG.GL)
    IG.duration(200)
    MyParticle(IG.v(0,0,0), IG.v(10,20,0))
    MyGravity(IG.v(0,-10,0))

class MyParticle(IParticle) : 
    def __init__(self, p, v) : 
        IParticle.__init__(self,p,v)
    
    def update(self) :
        IPoint(self.pos().cp())

class MyGravity(IAgent) : 
    
    def __init__(self, g) : 
        self.gravity = g
    
    def interact(self, agents) : 
        for agent in agents : 
            if isinstance(agent, MyParticle) : 
                agent.push(self.gravity)

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

add_library('igeo')

def setup() : 
    size(480, 360, IG.GL)
    IG.duration(200)
    num = 40
    for i in range(num) : 
        a = i*2*PI/num
        MyParticle(IG.v(20,0,0).rot(a),IG.v(10,0,i).rot(a))
    MyGravity(IG.v(0,0,-10))

class MyParticle(IParticle) : 
    def __init__(self, p, v) : 
        IParticle.__init__(self,p,v)
    
    def update(self) :
        IPoint(self.pos().cp())

class MyGravity(IAgent) :
    def __init__(self, g) : 
        self.gravity = g
  
    def interact(self, agents) : 
        for agent in agents : 
            if isinstance(agent, MyParticle) : 
                agent.push(self.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);

add_library('igeo')

def setup() : 
    size(480, 360, IG.GL)
    IG.duration(200)
    num = 40
    for i in range(num) : 
        a = i*2*PI/num
        MyParticle(IG.v(20,0,0).rot(a), IG.v(10,0,i).rot(a))

    MyGravity(IG.v(0,0,-10))
    BouncePlane(-5)

class MyParticle(IParticle) : 
    def __init__(self, p, v) : 
        IParticle.__init__(self,p,v)
    
    def update(self) :
        IPoint(self.pos().cp())

class MyGravity(IAgent) :
    def __init__(self, g) : 
        self.gravity = g
  
    def interact(self, agents) : 
        for agent in agents : 
            if isinstance(agent, MyParticle) : 
                agent.push(self.gravity)

class BouncePlane(IAgent) : 
    def __init__(self, zval) : 
        self.z = zval
    
    def interact(self, agents) : 
        for agent in agents : 
            if isinstance(agent, MyParticle) : 
                if agent.pos().z() < self.z : 
                    agent.pos().add(0,0,self.z-agent.pos().z())
                    agent.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.

add_library('igeo')

def setup() : 
    size(480, 360, IG.GL)
    IG.duration(200)
    num = 40
    for i in range(num) : 
        a = i*2*PI/num
        MyParticle(IG.v(20,0,0).rot(a), IG.v(10,0,i).rot(a))
    MyGravity(IG.v(0,0,-10))
    BouncePlane(-5)

class MyParticle(IParticle) : 
    def __init__(self, p, v) : 
        IParticle.__init__(self,p,v)
        self.prevPos = None
    
    def update(self) :
        curPos = self.pos().cp()
        if self.prevPos is not None : 
            ICurve(self.prevPos, curPos).clr(0)
        self.prevPos = curPos

class MyGravity(IAgent) :
    def __init__(self, g) : 
        self.gravity = g
  
    def interact(self, agents) : 
        for agent in agents : 
            if isinstance(agent, MyParticle) : 
                agent.push(self.gravity)

class BouncePlane(IAgent) : 
    def __init__(self, zval) : 
        self.z = zval
    
    def interact(self, agents) : 
        for agent in agents : 
            if isinstance(agent, MyParticle) : 
                if agent.pos().z() < self.z : 
                    agent.pos().add(0,0,self.z-agent.pos().z())
                    agent.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.

add_library('igeo')

def setup() : 
    size(480, 360, IG.GL)
    IG.duration(200)
    num = 40
    particle = []
    for i in range(num) : 
        a = i*2*PI/num
        particle.append(MyParticle(IG.v(20,0,0).rot(a),IG.v(10,0,i).rot(a)))
    for i in range(1, num) : 
        GeometryAgent(particle[i-1], particle[i])
    MyGravity(IG.v(0,0,-10))
    BouncePlane(-5)

class MyParticle(IParticle) : 
    def __init__(self, p, v) : 
        IParticle.__init__(self,p,v)
        self.prevPos = None
    
    def update(self) :
        curPos = self.pos().cp()
        if self.prevPos is not None : 
            ICurve(self.prevPos, curPos).clr(0)
        self.prevPos = curPos

class MyGravity(IAgent) :
    def __init__(self, g) : 
        self.gravity = g
  
    def interact(self, agents) : 
        for agent in agents : 
            if isinstance(agent, MyParticle) : 
                agent.push(self.gravity)

class BouncePlane(IAgent) : 
    def __init__(self, zval) : 
        self.z = zval
    
    def interact(self, agents) : 
        for agent in agents : 
            if isinstance(agent, MyParticle) : 
                if agent.pos().z() < self.z : 
                    agent.pos().add(0,0,self.z-agent.pos().z())
                    agent.vel().ref(IG.zaxis)

class GeometryAgent(IAgent) : 
    def __init__(self, p1, p2) : 
        self.particle1 = p1
        self.particle2 = p2
        self.prevPos1 = None
        self.prevPos2 = None
    
    def update(self) : 
        curPos1 = self.particle1.pos().cp()
        curPos2 = self.particle2.pos().cp()
        if self.prevPos1 is not None and self.prevPos2 is not None :
            ISurface(self.prevPos1, curPos1, curPos2, self.prevPos2).clr(self.clr())
        self.prevPos1 = curPos1
        self.prevPos2 = curPos2
        # update color
        r = self.red() + IRand.get(-0.01,0.01)
        g = self.green() + IRand.get(-0.02,0.02)
        self.clr(r, g, 1.0)