Home My Page Projects Code Snippets Project Openings diderot

# SCM Repository

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

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

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 };

```