Home My Page Projects Code Snippets Project Openings diderot

SCM Repository

[diderot] View of /tests/vis15-bugs/halftone-bug.diderot
 [diderot] / tests / vis15-bugs / halftone-bug.diderot

View of /tests/vis15-bugs/halftone-bug.diderot

Thu Oct 6 18:24:09 2016 UTC (2 years, 10 months ago) by glk
File size: 6356 byte(s)
`something about the pcIter function`
```vec2{} ipos = {[0,0],[0,0.1]};
input vec2 radmm ("particle minimum > 0 and maximum radius") = [0.01, 0.1];
input real eps ("system convergence threshold, computed as the coefficient-of-variation of normalized distances to nearest neighbors") = 0.05;
input int pcp ("periodicity of population control (if greater than zero)") = 2;
input real hhInit ("initial step size for potential energy gradient descent") = 1;

real newDist = 0.55;  // how far away to birth, as multiple of radius
real stepMax = 1;     // limit on travel per iter, as multiple of radius
int iter = 0;         // which iteration we're on

// Univariate energy functions; see ../circle/circle.diderot for alternatives
// well radius 2/3, depth -0.02
function real phi(real r) =
1 + r*(-3.87 + r*(4.725 - r*1.8225))  if r < 0.666666666666 else
2.7 +  r*(-12.96 + r*(22.68 + r*(-17.28 + r*4.86)));
function real phi'(real r) =
-0.0225*(-2 + 3*r)*(-86 + 81*r)  if r < 0.666666666666 else
6.48*(r - 1)*(r - 1)*(-2 + r*3);
// Energy and force from particle (with radius rad) at vec2 x

// From a given vec2, find a random-ish value uniformly sampling [0,1)
function real posrnd(vec2 v, real rad) {
return fmod(|fmod(p[0],1)| + |fmod(p[1],1)|, 1);
}

/* this one compiles okay
function int pcIter() = 1 if (pcp > 0 && iter > 0 && 0 == iter % pcp) else 0;
*/
/* BUG this one crashes the compiler */
function int pcIter() {
if (!(pcp > 0 && iter > 0 && 0 == iter % pcp)) {
return 0;
}
return 1 if (iter/pcp) % 2 == 0 else -1;
}

/* The particle is initialized at position pos0, with initial stepsize hh0.
The first set of particles gets hhInit for the initial stepsize, but new
particles created by population control benefit from getting the stepsize
that was adaptively learned by the parent. */
strand particle (vec2 pos0, real hh0) {
output vec2 pos = pos0;
real hh = hh0;
vec2 step = [0,0];   // step along force
real rad = 0;        // zero isn't a valid value
real closest = 1;    // distance to closest neighbor (as mult of radius)
int ncount = 0;      // how many neighbors did we have
real undone = 1;     // same as in ../sphere
update {
if (|pos[0]| > 1 || |pos[1]| > 1) { die; }
// compute energy and forces on us from neighbors
real energyLast = 0;
vec2 force = [0,0];
ncount = 0;
foreach (particle P in sphere(rad)) {
vec2 x = P.pos - pos;
if (|x| == 0) {
die; // same rationale as in ../sphere
}
ncount += 1;
}
if (0 == ncount) {
if (pcIter() > 0) {  // no neighbors, so let's make one
vec2 npos = pos + newDist*rad*[cos(a),sin(a)];
new particle(npos, hh);
undone = 1;
} else {
undone *= pow(0.5, 1.0/(2*(pcp if pcp > 0 else 1)));
}
// set closest to something in case used in global update
closest = newDist;
continue;
}

/* Else have interacting neighbors; find step, limit step size */
step = hh*force;
step = hh*force;
}

// Take step, re-apply constraint, find new energy
vec2 posLast = pos;
pos += step;
real energy = 0;
closest = 1;
ncount = 0;
vec2 mnd = [0,0]; // mean neighbor difference
foreach (particle P in sphere(rad)) {
vec2 nd = P.pos - pos;
mnd += nd;
ncount += 1;
}

// Armijo-Goldstein sufficient decrease condition
if (energy - energyLast > 0.5*(pos - posLast)•(-force)) {
hh *= 0.5;  // energy didn't decrease as expected: backtrack
pos = posLast;  // try again next iteration
// no progress, so decrease of undone
} else {
hh *= 1.02; // opportunistically increase step size
// indicate progress; may be over-written below
undone *= pow(0.5, 1.0/(2*(pcp if pcp > 0 else 1)));
// try to have between 5 and 8 neighbors
if (pcIter() != 0) {
if (ncount < 6) {
if (pcIter() > 0) {
}
undone = 1;
} else if (ncount > 8 && energy > 0 && pcIter() < 0) {
/* the logic for using posrnd() this way is the same as in
../sphere, but here we can exploit how the potential
function has a (negative) well, so energy > 0 means that the
total system energy would be lower without this particle */
if (posrnd(pos, rad) < (ncount - 8.0)/ncount) {
die;
}
// else not done if too many neighbors, w/ population control
undone = 1;
}
}
}

// Record actual step taken, in case used in global update
step = pos - posLast;
}
}

global {
/* Compute coefficient-of-variation of distance to closest neighbor */
real meancl = mean { P.closest | P in particle.all};
real varicl = mean { (P.closest - meancl)^2 | P in particle.all};
real covcl = sqrt(varicl) / meancl;
real meanncount = mean { real(P.ncount) | P in particle.all};
real maxundone = max { P.undone | P in particle.all};
print("(iter ", iter, "; pcIter ", pcIter(), ") COV(cl)=", covcl,
"; mean(cl)=", meancl,
"; mean(ncount)=", meanncount,
"; max(undone)=", maxundone, "\n");

if (covcl < eps           // seem to be geometrically uniform
&& maxundone < 0.5) { // and no particle recently set undone=1
print("Stabilizing ", numActive(), " points with COV(closest) ", covcl,
" < ", eps, " and maxundone ", maxundone, " < 0.5 (iter ", iter, ")\n");
stabilize;
}
iter += 1;
}

initially { particle(ipos[ii], hhInit) | ii in 0 .. length(ipos)-1 };
```