6 steps for Long Range Interaction

Share your projects and Chipmunk enhancements.
Post Reply
User avatar
juanpi
Posts: 24
Joined: Wed Jan 30, 2008 3:52 pm
Contact:

6 steps for Long Range Interaction

Post by juanpi »

My friends,
Here is the implementation of the long range interaction.
1. We need to use the void* data that is in every cpBody structure
2. This pointer will point a new structure called Sing (from singularity). Here is an optional code you can do it as you wish.

Code: Select all

typedef struct paramActSingularity {

	/** @brief The type of singularity */ 
    char* type;
	/** @brief Scalar value of singularity */ 
    cpFloat value;
	/** @brief Position in the actor in reference to the CoM (body coordinates)*/ 
    vect position;

} acSingularity;

typedef struct ActorSingularity{
       	/** @brief Number of singularities*/ 
	int Nsing; 
       	/** @brief Array of singularities*/ 
	acSingularity* sing; 
       	/** @brief Global position of the singularities*/ 
	cpVect* Gpos; 
}Sing;
2b. Now we are ready to define long range interaction between two actors.

3. Add this code to your function run() or equivalent

Code: Select all

       //Reset forces
	cpSpaceEachBody(space, &eachBody, NULL);		

	// Long range interaction N*N/2 instead of N*N
	for(int j=0; j<bodies->num; j++)
	{
	  for(int k=j+1; k<bodies->num; k++)
	  {
	     LRangeForceApply((cpBody*)bodies->arr[j],(cpBody*)bodies->arr[k]);
	  }
	}

	cpSpaceStep(space, dt);
4. Define the LRangeForceApply(cpBody*, cpBody *) function

Code: Select all

void LRangeForceApply(cpBody *a, cpBody *b){
       
        Sing *aux = (Sing*)a->data;
        Sing *aux2 = (Sing*)b->data;
        cpVect delta;
        cpFloat r2=0.0;
        cpVect f;

        // Calculate the forces between the singularities of different bodies
        for (int i=0; i<aux->Nsing; i++)
        {
         for (int j=0; j<aux2->Nsing; j++)
         {
       	   delta = cpvsub(aux2->Gpos[j],aux->Gpos[i]);
	   r2 = cpvlengthsq(delta);

          // Coulomb interaction -k*q1*q2/r^2
	   f=cpvmult(cpvnormalize(delta),-aux->sing[i].value*aux2->sing[j].value/r2);
	   
	   //Force applied to body A
	   cpBodyApplyForce(a,f,aux->Gpos[i]);
	 
	   //Reaction on body B
	   cpBodyApplyForce(b,cpvneg(f),aux2->Gpos[j]);
	 }
        }
}
4b. The value of f has the vlaue of the force in the body. There you can put any force you can emulate as interaction of two point singularities. (for example magnets could be thought of two charges of the same sing very close together)
4c. Observe that the singularities of the same body do not produce forces.

5. To add the Singularities to the bodies do the following (for example)

Code: Select all

		Sing* aux;
                aux->Nsing=1; // Here the number of singularities per body

               acSingularity singularity;
               singularity.type="electrical";
               singularity.value=1.0e0;
               singularity.position=cpv(0.0e0,1.0e0); //Slightly above the center of mass

                cpArray* bodies=space->bodies;
		for (int i=0; i<bodies->num; i++)
		{
                 // Memory allocation
		  aux=(Sing*)malloc(sizeof(Sing));
   		  aux->sing=(acSingularity*)malloc(aux->Nsing*sizeof(acSingularity));
		  aux->Gpos=(cpVect*)malloc(aux->Nsing*sizeof(cpVect));
		  
		  cpBody *body = (cpBody *)bodies->arr[i];
		  
                  // Loop over the singularities of each body
                  for(int k=0;k<aux->Nsing;k++)
		  {
			aux->sing[k]=singularity[k];
			aux->Gpos[k]=cpvadd(body->p,cpv(aux->sing[k].position.x,aux->sing[k].position.y));
		  }
		  body->data=aux;
                  body->position_func=SingularBodyUpdatePosition;
	   }
6. The function to update the position is just simply modified by

Code: Select all

void ChargedBodyUpdatePosition(cpBody *body, cpFloat dt)
{
/** @todo Increase the order of the solvers */
	cpVect dp = cpvmult(cpvadd(body->v, body->v_bias), dt);
	body->p = cpvadd(body->p, dp);

	cpBodySetAngle(body, body->a + (body->w + body->w_bias)*dt);

	// Update position of the charges
	aux = (Sing*)body->data;
	for (int i=0; i<aux->Nsing; i++)
            aux->Gpos[i]=cpvadd(body->p,cpvrotate(cpv(aux->sing[i].position.x,
													  aux->sing[i].position.y), body->rot));
            
 	body->v_bias = cpvzero;
	body->w_bias = 0.0f;
}
There are some preliminaries videos on Youtube (search for Chipmunk Long Range interaction). I hope soon there will be more (I am programming Velocity Verlet algorithm to improve energy conservation (and possibly if somebody wants to do it, improve the calculations of some constraints)) http://en.wikipedia.org/wiki/Verlet_integration

A description of the basic forces you would like to apply can be found here
http://en.wikipedia.org/wiki/Force

Enjoy!
maximile
Posts: 157
Joined: Mon Aug 20, 2007 12:53 pm
Location: London, UK
Contact:

Re: 6 steps for Long Range Interaction

Post by maximile »

Looking very interesting. Would this be useful for, say, simulating planets orbiting a star?
User avatar
juanpi
Posts: 24
Joined: Wed Jan 30, 2008 3:52 pm
Contact:

Re: 6 steps for Long Range Interaction

Post by juanpi »

lll
Last edited by juanpi on Mon Aug 18, 2008 1:24 am, edited 1 time in total.
User avatar
juanpi
Posts: 24
Joined: Wed Jan 30, 2008 3:52 pm
Contact:

Re: 6 steps for Long Range Interaction

Post by juanpi »

Yes, of course. But you will have some trouble keeping planets in orbits for a long period of time with a low order solver like Euler.
I have implemnted a algorithm that conserves energy to higher order without too much calculation (check my other post Velocity Verlet). Anyway you can try and comment the results. If you want to implement a hack you could restart the orbit of the planets from time to time...but I am sure there are more intelligent hacks for this somewhere.

so...good luck!

PS:
Gravity force
F= - m1 * m2 / r^2
and applies in the direction of the vector r, going from one mass to the other.
Post Reply

Who is online

Users browsing this forum: No registered users and 6 guests