# SCM Repository

# View of /tests/vis15-bugs/circle.diderot

Parent Directory | Revision Log

Revision

File size: 3880 byte(s)

**4729**- (**download**) (**annotate**)*Tue Oct 11 21:32:41 2016 UTC*(2 years, 10 months ago) by*glk*File size: 3880 byte(s)

adding script for (fixed) bug with saving first snapshot

/* running this under valgrind produces: Syscall param write(buf) points to uninitialised byte(s)" that seems to be about the memory in the nrrd not being set due to saving out the first (post-initialization, before first iteration) snapshot */ input vec2{} ipos ("initial positions for all particles") = load("vec2.nrrd"); input real rad ("radius of particle's potential energy support") = 0.25; input real eps ("system convergence threshold, computed as the coefficient-of-variation of distances to nearest neighbors") = 0.05; real hhInit = 1; // initial step size real stepMax = rad; // limit on distance to travel per iter int iter = 0; // which iteration we're on // phi(0) and phi'(0) are bounded function real phi(real r) = (1 - r)^4; function real phi'(real r) = -4*(1 - r)^3; function real enr(vec2 x) = phi(|x|/rad); function vec2 frc(vec2 x) = phi'(|x|/rad) * (1/rad) * x/|x|; // chain rule /* The particle is initialized at position pos0 */ strand particle (vec2 pos0) { /* These strand variables are visible to the global update and to other strands making spatial queries. Any variable inside the scope of update{} will not be visible in this way. NOTE: "pos" is a special variable name that *must* be used to enable spatial queries of neighboring particles via sphere(). */ output vec2 pos = pos0/|pos0|; real hh = hhInit; vec2 step = [0,0]; // step along force real closest = 0; // distance to closest neighbor int ncount = 0; // how many neighbors did we have update { // Compute energy and forces on us from neighbors real energyLast = 0; vec2 force = [0,0]; ncount = 0; foreach (particle P in sphere(rad)) { energyLast += enr(P.pos - pos); force += frc(P.pos - pos); ncount += 1; } if (0 == ncount) { /* no neighbors; do nothing to do but set closest to big value */ closest = rad; continue; } // Else we have interacting neighbors vec2 norm = normalize(pos); // surface normal for unit circle // project force onto tangent surface force = (identity[2] - norm⊗norm)•force; /* Limiting the step size (even before testing for sufficient decrease, below) helps keep particles near the surface they are supposed to be sampling; this precaution matters more with a non-trivial (data-driven) constraint surface */ step = hh*force; // compute step along force if (|step| > stepMax) { // decrease hh by factor by which step was too big hh *= stepMax/|step|; // and find smaller step (of length stepMax) step = hh*force; } // take step and re-apply constraint vec2 posLast = pos; pos = normalize(pos + step); // find energy at new location, and distance to closest neighbor real energy = 0; closest = rad; foreach (particle P in sphere(rad)) { energy += enr(P.pos - pos); closest = min(closest, |P.pos - pos|); } // save actual step taken step = pos - posLast; if (energy - energyLast > 0.5*(pos - posLast)•(-force)) { hh *= 0.5; // backtrack pos = posLast; } else { hh *= 1.02; // opportunistically increase step size } } } global { 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; print("(iter ", iter, ") mean(hh)=", mean { P.hh | P in particle.all}, "; COV(closest) = ", covcl, "\n"); if (covcl < eps) { print("Stabilizing with COV(closest) ", covcl, " < ", eps, " (iter ", iter, ")\n"); stabilize; } iter += 1; } initially { particle(ipos[ii]) | ii in 0 .. length(ipos)-1 };

root@smlnj-gforge.cs.uchicago.edu | ViewVC Help |

Powered by ViewVC 1.0.0 |