Colisiones y Partículas con lenguaje de programación Processing

Los sistemas de partículas vienen empleándose en disciplinas como la Animación, los Gráficos 3D y la Realidad Virtual desde tiempos inmemoriales. Permiten simular modelos deformables como por ejemplo el agua, fenómenos atmosféricos como un tornado o las nubes e incluso comportamientos sociales distribuidos entre “mini seres” con algo de IA, como los peces de una pecera virtual o los pájaros de una bandada “sintética” que vuela por los cielos virtuales de un Videojuego. Sólo tenéis que mirar lo que Ron Fedkiw hace con ellas o pasaros por los tutoriales de David Baraff, por cierto éste último trabajando en Pixar.

Vamos a implementar un sencillo simulador en el que una partícula rebotará en el suelo. Partirá de una posición elevada y de una cierta velocidad inicial “horizontal” y perderá energía cada vez que impacte en el suelo, rebotando acorde a lo esperado. Todo fruto de la fuerza de gravedad. Lo podemos observar en la figura 1.

Figura 1. Nuestra partícula rebotará en el suelo gracias al cálculo de colisiones.

En la figura 1 observamos como nuestra partícula (esferas roja) parte de una posición elevada con una cierta velocidad inicial. Cae debido a la fuerza de la gravedad que terminará frenándola lo suficiente como para que pare.

Necesitamos que nuestra partícula disponga de masa (gramos o kilogramos, o ejemplo), posición XY y también velocidad XY. Además tiene que verse afectada por fuerzas. En nuestro caso, la de la gravedad, que según los preceptos de la Física:

Fuerza de Gravedad = Masa x g

Dónde g es la constante de gravitación universal para nuestro planeta Tierra, de valor 9.81 metros por segundo al cuadrado. Es curioso como en algunos Videojuegos en los que se desea tener en cuenta este efecto, por ejemplo en uno de futbol en el que la pelota bota por el campo, a veces se varía la constante (haciéndola valer 10 por ejemplo) para conseguir exagerar el efecto un poco y, aunque no siendo entonces físicamente correcto, aparezca como visualmente más plausible.

Bueno, lo último que nos queda por recordar en esta mini lección de Física es la segunda ley de Newton, por la que:

Fuerza = Masa x Aceleración

Que de hecho nos recuerda notablemente a la fórmula anterior. Y es que si un objeto con masa recibe una fuerza…se acelera! Y si eso ocurre a lo largo del tiempo, que transcurre, varía su velocidad…y por ende su posición! Ya tenemos el mecanismo por el cual se desplazará nuestra partícula pues!.

En un ordenador hace falta un procedimiento por el cual, vía el lenguaje de programación, implementemos ese movimiento. Necesitamos que el tiempo transcurra y que a la par, calculemos las aceleraciones que afectan a las partículas y sus variaciones en velocidad y posición. A estos procesos se les llama motores de inferencia numérica o “solvers”. Y existen distintos tipos (Euler, Runge-Kutta, Midpoint, Verlet, LepFrog, etc.) según su velocidad de ejecución y eficiencia.

Como en nuestro caso vamos a realizar una simulación bien simple emplearemos el más sencillo de todos. El “solver Euler”. Tenemos “las fórmulas” en la figura 2.

Figura 2. Fórmulas de cálculo de nuevas velocidades y posiciones gracias al solver Euler.

Tengamos en cuenta que con X nos referimos a la posición y que con V a la velocidad. Supongamos que nuestra partícula (círculo rojo) se encuentra en una posición actual Xn con una velocidad Vn (en el dibujo vector azul D). Lo que deseamos es desplazarla desde el momento de tiempo actual (t) al siguiente (t + h dónde h es el incremento de tiempo que transcurre). Lo que haremos será calcular las nuevas posición y velocidad Xn+1 y Vn+1 (en “el ramo” se denomina a este “par” como el “espacio de fase”) a partir de las anteriores Xn y Vn. Y eso lo hacemos gracias a las fórmulas indicadas en la figura 2:

• La nueva velocidad Vn+1 se obtiene como la anterior Vn más un factor que no es más que un incremento o decremento de velocidad. Este incremento o decremento se obtiene multiplicando el tiempo que transcurre (Dt o h) por la aceleración…y ésta se obtiene según la segunda Ley de Newton, como la fuerza f que afecta a la partícula dividida entre su masa m.
• La nueva posición Xn+1 se obtiene como la anterior Xn más un factor que no es más que un incremento de posición. Este incremento o decremento se obtiene multiplicando el tiempo que transcurre (Dt o h) por la velocidad.

Imaginad que dentro de nuestro bucle de iteración Processing, es decir dentro de la función draw(), calculamos continuamente la posición de nuestra partícula y a continuación la pintamos. El resultado será una animación en la que “se mueve sola” acelerándose y frenándose por acción de la fuerza de la gravedad. Vamos a ello pues y para realizarlo utilizaremos objetos de Processing. Veamos el código 1.

// A particle system with Processing
// One circular particle will fall from sky
// due to Gravity. Will Bounce if touching the floor
// until it lacks energy and stops!

// We begin by creating our particle object
particle part;

// We create several accumulators
// as variables to store forces,
// positions and velocities
float[] force = new float[2];
float[] oldPosition = new float[2];
float[] oldVelocity = new float[2];
// We also create a time variable for the system to evolve
float time;

// Our initialization setup function
void setup() {
// a 500 x 500 window
size(500,500);
// We initialize our particle object
// by calling the constructor of its class
// and passing several parameters toit
// which are color, radius and the initial
// 2D position and 2D velocity
part = new particle(color(255,0,0),50,50,450,15,0,1);
// The force accumulator in X equals 0 (there's no gravitational component)
force[0]=0.0;
// The force accumulator in Y equals -9.8 (this is Gravity)
force[1]=-9.8;
// Time will pass as 0.1 seconds per iteration
time=0.1;
}

Código 1. Objetos, variables e inicialización de nuestras partículas, 5 en total.

Observamos que lo primero es crear nuestra partícula como un objeto de un tipo que crearemos para tal ocasión (particle). Además preparamos 4 variables, 3 de ellas arrays, que necesitaremos (force, oldPosition, oldVelocity y time) para almacenar la fuerza, posiciones y velocidades antiguas de la partícula así como el tiempo, que transcurre, en nuestra simulación.

En la función de inicialización setup() creamos una ventana de 500 por 500 píxeles e inicializamos a nuestra partícula llamando a su constructor. En el código 3 observaremos como es este método. Básicamente crea a la partícula a partir de los parámetros de entrada que recibe, que son color, radio del círculo, posición y velocidad iniciales, en 2D por nuestra ventana y por tanto en coordenadas XY, y para terminar la masa. También definimos como deseamos que sea nuestra fuerza. Como sólo tenemos la gravedad y ésta actúa tan sólo verticalmente, ponemos a 0 la componente horizontal de la fuerza e igualamos a 9.8 m/s2 la vertical, tal y como ya hemos comentado anteriormente. Terminamos la función definiendo un incremento de tiempo (time) de 0.1 segundos para nuestra simulación. Esto quiere decir que para cada vez que iteremos y pintemos de nuevo, habrán transcurrido 0.1 segundos “para la partícula”. Se moverá por tanto, entre cuadro y cuadro generado, el espacio correspondiente a este tiempo.

En el código 2 podemos ver como nuestro bucle infinito, es decir la función draw(), refresca el color de fondo en negro y llama continuamente a los métodos asociados a nuestro objeto particle para a) Calcular la nueva posición ésta mediante el método particleMove(), b) Detectar la colisión con el suelo en caso de que la hubiera gracias a particleCollision () y c) La muestra en pantalla gracias al método particleRender(). En el código 3 veremos en qué consisten ambos métodos.

// The infinite loop
void draw() {
// Black background
background(0);
// We call our particle methods
// in order to 1) Move it
// 2) Evaluate collisions with the floor
// and 3) Render it
part.particleMove();
part.particleCollision();
part.particleRender();
}

Código 2. El bucle de simulación de nuestras partículas.

Vamos por tanto al código 3 de cara a tener claro como se define nuestro objeto “personalizado” de tipo particle así como sus variables y métodos internos.

// This is the particle class
// that defines our particular object
class particle {
// A particle has several attributes as:
// color, radius, 2D position, 2D velocity
// 2D acceleration and mass
color color_p;
float particleRadius;
float[] position = new float[2];
float[] velocity = new float[2];
float[] acceleration = new float[2];
float mass;
// The constructor for a particle
particle(color c,float r,float posX,float posY,float velX,float velY,float m) {
color_p=c;
particleRadius=r;
position[0]=posX;
position[1]=posY;
velocity[0]=velX;
velocity[1]=velY;
mass=m;
}
// The method that renders a particle
// as a circle with white stroke and its filling color
void particleRender() {
stroke(255);
fill(color_p);
// Remember that Y coords go top to down in Processing
// and therefore we need to invert the particle's position
// if we want our fireworks to go bottom-up
ellipse(position[0],height-position[1],particleRadius,particleRadius);
}
// The method that moves a particle
void particleMove() {
// First: save the old state
oldPosition[0]=position[0];
oldPosition[1]=position[1];
oldVelocity[0]=velocity[0];
oldVelocity[1]=velocity[1];
// Second: evaluate particle's acceleration as:
// acceleration = force/mass (this is one of Newton's laws)
acceleration[0]=force[0]/mass;
acceleration[1]=force[1]/mass;
// Third: evaluate new velocity from old velocity
// by using the Euler method that stands for:
// new velocity = old velocity + time * acceleration
velocity[0]=oldVelocity[0] + time*acceleration[0];
velocity[1]=oldVelocity[1] + time*acceleration[1];
// Fourth: evaluate new position from old position
// by using the Euler method that stands for:
// new position = old position + time * velocity
position[0]=oldPosition[0] + time*velocity[0];
position[1]=oldPosition[1] + time*velocity[1];
}
// The method that evaluates particle collision
// and restores its velocity to make it bounce
// If our particle reaches the bottom
// we invert its velocity in order to make it bounce
// and we also apply a restitution factor to
// provide with some energy loss (5%) as in reality
// We are implementing an inelastic simulation
void particleCollision(){
if (height-position[1]>=height){
velocity[1] = velocity[1] * -0.95;
}
}
}

Código 3. Definición de nuestros objetos Processing, de tipo particle.

Para crear un objeto en Processing utilizamos una estructura de tipo clase. Dentro de ella codificaremos el constructor que inicializa a cada uno de los objetos, en nuestro caso partículas, colocando los valores “que toque” a cada una de sus variables internas. Además se definen los métodos internos de los objetos, es decir aquellas funciones que se requieren para manejarlos. En definitiva:

• Un objeto particle se define por su color, particleRadius, position, velocity, acceleration y mass.
• El constructor (particle(color c,float r,float posX,float posY,float velX,float velY,float m)) recibe como parámetros estos valores y los introduce en cada una de las variables. Tanto la posición como la velocidad se refieren a las iniciales para cada partícula. Esencialmente centradas en la ventana y con velocidad vertical ascendente.
• El método que “pinta” la partícula es void particleRender(). Básicamente define color de borde blanco, de relleno el de la misma partícula y la pinta como un círculo de radio igual a particleRadius, centrado en la posición X = position[0] e Y = height-position[1]. Ésta última invertida dado que en Processing las Y’s crecen hacia abajo y en cambio queremos que nuestras partículas se muevan ascendentemente.
• El método que “mueve” la partícula y por tanto implementa al mencionado “Solver Euler” es void particleMove(). Tendrá que implementar las fórmulas que vimos en su momento en la figura 2. Así pues, empieza por almacenar tanto las posiciones XY como las velocidades XY actuales (oldPosition y oldVelocity), calcula la aceleración para la partícula actual según la segunda Ley de Newton (acceleration = force/mass) y por último calcula la posición y velocidad actuales, digamos “nuevas”, a partir de las “antiguas (velocity = oldVelocity + time * acceleration y position = oldPosition + time * velocity).
• Por último, el método que evalúa la colisión y reacciona ante ésta es particleCollision(). Lo que hacemos es realmente sencillo aunque puede resultar confuso si observamos el condicional. Recordemos que la coordenada Y crece hacia abajo y que nosotros la estamos tratando justamente al revés (botamos hacia arriba). Por eso operamos con height-position[1] y no con position[1]. Miramos si hemos llegado al “suelo· de la ventana (height) y si es así, invertimos la componente vertical de la velocidad (la multiplicamos por -1) para implementar el bote. Además conseguimos “disminuirla” multiplicándola por un coeficiente de restitución de 0.95 (pérdida del 0.05% de la energía a cada bote). Éste puede variarse de manera que si este coeficiente disminuye, se perderá más energía a cada bote en el suelo.

Por último, el callback de teclado keyPressed() se encarga de que el sistema vuelva a reiniciarse (llamamos de nuevo al constructor) cada vez que se pulse una tecla cualquiera.

// The key press event will call this callback
void keyPressed(){
// We reset the system everytime that the user presses a key
// by re-calling the particle's constructor
part = new particle(color(255,0,0),50,50,450,15,0,1);
}

Código 4. Callback de teclado. Cada vez que pulsemos una tecla la simulación se reiniciará.

// A particle system with Processing
// One circular particle will fall from sky
// due to Gravity. Will Bounce if touching the floor
// until it lacks energy and stops!

// We begin by creating our particle object
particle part;

// We create several accumulators
// as variables to store forces,
// positions and velocities
float[] force = new float[2];
float[] oldPosition = new float[2];
float[] oldVelocity = new float[2];
// We also create a time variable for the system to evolve
float time;

// Our initialization setup function
void setup() {
// a 500 x 500 window
size(500,500);
// We initialize our particle object
// by calling the constructor of its class
// and passing several parameters toit
// which are color, radius and the initial
// 2D position and 2D velocity
part = new particle(color(255,0,0),50,50,450,15,0,1);
// The force accumulator in X equals 0 (there's no gravitational component)
force[0]=0.0;
// The force accumulator in Y equals -9.8 (this is Gravity)
force[1]=-9.8;
// Time will pass as 0.1 seconds per iteration
time=0.1;
}

// The infinite loop
void draw() {
// Black background
background(0);
// We call our particle methods
// in order to 1) Move it
// 2) Evaluate collisions with the floor
// and 3) Render it
part.particleMove();
part.particleCollision();
part.particleRender();
}

// This is the particle class
// that defines our particular object
class particle {
// A particle has several attributes as:
// color, radius, 2D position, 2D velocity
// 2D acceleration and mass
color color_p;
float particleRadius;
float[] position = new float[2];
float[] velocity = new float[2];
float[] acceleration = new float[2];
float mass;
// The constructor for a particle
particle(color c,float r,float posX,float posY,float velX,float velY,float m) {
color_p=c;
particleRadius=r;
position[0]=posX;
position[1]=posY;
velocity[0]=velX;
velocity[1]=velY;
mass=m;
}
// The method that renders a particle
// as a circle with white stroke and its filling color
void particleRender() {
stroke(255);
fill(color_p);
// Remember that Y coords go top to down in Processing
// and therefore we need to invert the particle's position
// if we want our fireworks to go bottom-up
ellipse(position[0],height-position[1],particleRadius,particleRadius);
}
// The method that moves a particle
void particleMove() {
// First: save the old state
oldPosition[0]=position[0];
oldPosition[1]=position[1];
oldVelocity[0]=velocity[0];
oldVelocity[1]=velocity[1];
// Second: evaluate particle's acceleration as:
// acceleration = force/mass (this is one of Newton's laws)
acceleration[0]=force[0]/mass;
acceleration[1]=force[1]/mass;
// Third: evaluate new velocity from old velocity
// by using the Euler method that stands for:
// new velocity = old velocity + time * acceleration
velocity[0]=oldVelocity[0] + time*acceleration[0];
velocity[1]=oldVelocity[1] + time*acceleration[1];
// Fourth: evaluate new position from old position
// by using the Euler method that stands for:
// new position = old position + time * velocity
position[0]=oldPosition[0] + time*velocity[0];
position[1]=oldPosition[1] + time*velocity[1];
}
// The method that evaluates particle collision
// and restores its velocity to make it bounce
// If our particle reaches the bottom
// we invert its velocity in order to make it bounce
// and we also apply a restitution factor to
// provide with some energy loss (5%) as in reality
// We are implementing an inelastic simulation
void particleCollision(){
//if (height-position[1]>=height-0.01*height){
if (height-position[1]>=height){
velocity[1] = velocity[1] * -0.95;
}
}
}

// The key press event will call this callback
void keyPressed(){
// We reset the system everytime that the user presses a key
// by re-calling the particle's constructor
part = new particle(color(255,0,0),50,50,450,15,0,1);
}

Código 5. Partícula y colisiones. Programa completo.

Ahora sólo os queda copiar y pegar el código 5, que incluye el programa completo, y probar además de variar el ejemplo.