Home My Page Projects Code Snippets Project Openings diderot
Summary Activity Tracker Tasks SCM

SCM Repository

[diderot] View of /branches/lamont/test/implicit-surface/unit-circle.diderot
ViewVC logotype

View of /branches/lamont/test/implicit-surface/unit-circle.diderot

Parent Directory Parent Directory | Revision Log Revision Log


Revision 2225 - (download) (annotate)
Thu Feb 28 05:14:28 2013 UTC (7 years, 1 month ago) by lamonts
File size: 2851 byte(s)
Updated unit-circle
// unit-circle
//
// Demo of distributing particles on the unit circle 
//
real{} initPosns = load("../../data/implicit-surface/circlePosns.nrrd"); // file has 300 particle positions  
int numOfParticles = 100; // But I'm only look at 100 for right now. 
input real r = 0.08;  // particle radius 
input real d = 0.16;  // desired difference of radius (d < r && d > 0)
input real R = 0.15;  // neighbor radius 
real MAX_FORCE = 2; 
input real h = 0.09;  // integration step size to control how quickly iterations change position. 
input int maxi = 1000; // maximum number of steps in the program 
vec4 qWinDim = [-1,1,-1,1]; // i.e. [XMIN, XMAX, YMIN, YMAX]  (required for the query function)
vec2 qGridDim = [12,12];  // how many grid cells you want in each direction for the uniform grid (required for the query function) 
vec2 qCellDim = [r,r];    // the length in each direction for a cell (required for the query function) 
int steps = 1; 

strand Particle (int i, real posx, real posy) {
    vec2 pos = [posx,posy];
    output vec2 outPos  = pos;
    int id = i;  
    update {

        vec2 energy =[0.0,0.0];
        vec2 force = [0.0,0.0]; 

        foreach(Particle neighbor in sphere(R)){ 
            real neighborDist = dist (pos,neighbor.pos);   
            vec2 unitDirection = normalize(pos-neighbor.pos); 

            if(neighborDist >= 0 && neighborDist <= 1) {
                   force += ((-4.0 * ((1 - neighborDist)^3)) * unitDirection);     
                   energy += (((1 - neighborDist)^4) * unitDirection); 
             }
        }
        real angle = 0.0;

        //Looking at how strand #0 changes for 30 iterations. 
        /*if(id == 0) { 
           print("pos before = ", pos, "\n"); 
           print("force = ", force, "\n"); 
           angle = atan2(pos[1],pos[0])* 180/π; 
           angle = (angle + 360) if angle < 0 else angle; // atan2 rangs from (π,-π] so need to range from [0,2π)
           print("Angle = ", angle, "\n"); 
        }*/

        force *= h;
 
        real isZero = force • force; 

        // limit the force
 	    if ( isZero != 0 && |force| > MAX_FORCE) {
	         force = normalize(force); 
             force *= MAX_FORCE; 
	     }  
        
        pos += force;   
        pos = normalize(pos);

       /* if(id == 0) { 
           print("new pos = ", pos, "\n"); 
           angle = atan2(pos[1],pos[0])* 180/π; 
           angle = (angle + 360) if angle < 0 else angle; 
           print("new angle = ", angle, "\n"); 
           print("========= iteration done ==================\n"); 
        } */ 
        outPos = pos; 

        if(steps >= 30) { 
            stabilize;
        }
    }
}

// can add statements in here to do some global computation
global{
  steps+=1; 
} 

initially {Particle(i,initPosns{i*2},initPosns{i*2+1}) | i in 0 .. numOfParticles-1 };


root@smlnj-gforge.cs.uchicago.edu
ViewVC Help
Powered by ViewVC 1.0.0