feat: Points (#82)
All checks were successful
/ test (push) Successful in 2m17s

Replaces the former sole Element entity by two, Sphere and Point, both implementing an Element trait. Adds Point display, uses the former Element display for Sphere. Adds a new "canned" configuration, and the ability to add, select, and nudge Point entities.

Co-authored-by: Aaron Fenyes <aaron.fenyes@fareycircles.ooo>
Reviewed-on: #82
Co-authored-by: Vectornaut <vectornaut@nobody@nowhere.net>
Co-committed-by: Vectornaut <vectornaut@nobody@nowhere.net>
This commit is contained in:
Vectornaut 2025-05-01 19:25:13 +00:00 committed by Glen Whitney
parent 360ce12d8b
commit a2478febc1
9 changed files with 815 additions and 330 deletions

View file

@ -1,15 +1,23 @@
use nalgebra::{DMatrix, DVector, DVectorView, Vector3};
use nalgebra::{DMatrix, DVector, DVectorView};
use rustc_hash::FxHashMap;
use slab::Slab;
use std::{collections::BTreeSet, rc::Rc, sync::atomic::{AtomicU64, Ordering}};
use std::{
any::{Any, TypeId},
cell::Cell,
collections::BTreeSet,
rc::Rc,
sync::atomic::{AtomicU64, Ordering}
};
use sycamore::prelude::*;
use web_sys::{console, wasm_bindgen::JsValue}; /* DEBUG */
use crate::{
display::DisplayItem,
engine::{
Q,
change_half_curvature,
local_unif_to_std,
point,
realize_gram,
sphere,
ConfigSubspace,
@ -33,31 +41,87 @@ pub type ElementColor = [f32; 3];
static NEXT_ELEMENT_SERIAL: AtomicU64 = AtomicU64::new(0);
pub trait ProblemPoser {
fn pose(&self, problem: &mut ConstraintProblem, elts: &Slab<Element>);
fn pose(&self, problem: &mut ConstraintProblem, elts: &Slab<Rc<dyn Element>>);
}
#[derive(Clone, PartialEq)]
pub struct Element {
pub trait Element: ProblemPoser + DisplayItem {
// the default identifier for an element of this type
fn default_id() -> String where Self: Sized;
// create the default example of an element of this type
fn default(id: String, id_num: u64) -> Self where Self: Sized;
// the regulators that should be created when an element of this type is
// inserted into the given assembly with the given storage key
/* KLUDGE */
// right now, this organization makes sense because regulators identify
// their subjects by storage key, so the element has to be inserted before
// its regulators can be created. if we change the way regulators identify
// their subjects, we should consider refactoring
fn default_regulators(_key: ElementKey, _assembly: &Assembly) -> Vec<Rc<dyn Regulator>> where Self: Sized {
Vec::new()
}
fn id(&self) -> &String;
fn label(&self) -> &String;
fn representation(&self) -> Signal<DVector<f64>>;
// the regulators the element is subject to. the assembly that owns the
// element is responsible for keeping this set up to date
fn regulators(&self) -> Signal<BTreeSet<RegulatorKey>>;
// a serial number that uniquely identifies this element
fn serial(&self) -> u64;
// take the next serial number, panicking if that was the last one left
fn next_serial() -> u64 where Self: Sized {
// the technique we use to panic on overflow is taken from _Rust Atomics
// and Locks_, by Mara Bos
//
// https://marabos.nl/atomics/atomics.html#example-handle-overflow
//
NEXT_ELEMENT_SERIAL.fetch_update(
Ordering::SeqCst, Ordering::SeqCst,
|serial| serial.checked_add(1)
).expect("Out of serial numbers for elements")
}
// the configuration matrix column index that was assigned to the element
// last time the assembly was realized, or `None` if the element has never
// been through a realization
fn column_index(&self) -> Option<usize>;
// assign the element a configuration matrix column index. this method must
// be used carefully to preserve invariant (1), described in the comment on
// the `tangent` field of the `Assembly` structure
fn set_column_index(&self, index: usize);
}
// the `Element` trait needs to be dyn-compatible, so its method signatures can
// only use `Self` in the type of the receiver. that means `Element` can't
// implement `PartialEq`. if you need partial equivalence for `Element` trait
// objects, use this wrapper
#[derive(Clone)]
pub struct ElementRc(pub Rc<dyn Element>);
impl PartialEq for ElementRc {
fn eq(&self, ElementRc(other): &Self) -> bool {
let ElementRc(rc) = self;
Rc::ptr_eq(rc, &other)
}
}
pub struct Sphere {
pub id: String,
pub label: String,
pub color: ElementColor,
pub representation: Signal<DVector<f64>>,
// the regulators this element is subject to. the assembly that owns the
// element is responsible for keeping this set up to date
pub regulators: Signal<BTreeSet<RegulatorKey>>,
// a serial number, assigned by `Element::new`, that uniquely identifies
// each element
pub serial: u64,
// the configuration matrix column index that was assigned to this element
// last time the assembly was realized, or `None` if the element has never
// been through a realization
column_index: Option<usize>
column_index: Cell<Option<usize>>
}
impl Element {
impl Sphere {
const CURVATURE_COMPONENT: usize = 3;
pub fn new(
@ -65,83 +129,161 @@ impl Element {
label: String,
color: ElementColor,
representation: DVector<f64>
) -> Element {
// take the next serial number, panicking if that was the last number we
// had left. the technique we use to panic on overflow is taken from
// _Rust Atomics and Locks_, by Mara Bos
//
// https://marabos.nl/atomics/atomics.html#example-handle-overflow
//
let serial = NEXT_ELEMENT_SERIAL.fetch_update(
Ordering::SeqCst, Ordering::SeqCst,
|serial| serial.checked_add(1)
).expect("Out of serial numbers for elements");
Element {
) -> Sphere {
Sphere {
id: id,
label: label,
color: color,
representation: create_signal(representation),
regulators: create_signal(BTreeSet::default()),
serial: serial,
column_index: None
}
}
// the smallest positive depth, represented as a multiple of `dir`, where
// the line generated by `dir` hits the element (which is assumed to be a
// sphere). returns `None` if the line misses the sphere. this function
// should be kept synchronized with `sphere_cast` in `inversive.frag`, which
// does essentially the same thing on the GPU side
pub fn cast(&self, dir: Vector3<f64>, assembly_to_world: &DMatrix<f64>) -> Option<f64> {
// if `a/b` is less than this threshold, we approximate
// `a*u^2 + b*u + c` by the linear function `b*u + c`
const DEG_THRESHOLD: f64 = 1e-9;
let rep = self.representation.with_untracked(|rep| assembly_to_world * rep);
let a = -rep[3] * dir.norm_squared();
let b = rep.rows_range(..3).dot(&dir);
let c = -rep[4];
let adjust = 4.0*a*c/(b*b);
if adjust < 1.0 {
// as long as `b` is non-zero, the linear approximation of
//
// a*u^2 + b*u + c
//
// at `u = 0` will reach zero at a finite depth `u_lin`. the root of
// the quadratic adjacent to `u_lin` is stored in `lin_root`. if
// both roots have the same sign, `lin_root` will be the one closer
// to `u = 0`
let square_rect_ratio = 1.0 + (1.0 - adjust).sqrt();
let lin_root = -(2.0*c)/b / square_rect_ratio;
if a.abs() > DEG_THRESHOLD * b.abs() {
if lin_root > 0.0 {
Some(lin_root)
} else {
let other_root = -b/(2.*a) * square_rect_ratio;
(other_root > 0.0).then_some(other_root)
}
} else {
(lin_root > 0.0).then_some(lin_root)
}
} else {
// the line through `dir` misses the sphere completely
None
serial: Self::next_serial(),
column_index: None.into()
}
}
}
impl ProblemPoser for Element {
fn pose(&self, problem: &mut ConstraintProblem, _elts: &Slab<Element>) {
let index = self.column_index.expect(
format!("Element \"{}\" should be indexed before writing problem data", self.id).as_str()
impl Element for Sphere {
fn default_id() -> String {
"sphere".to_string()
}
fn default(id: String, id_num: u64) -> Sphere {
Sphere::new(
id,
format!("Sphere {id_num}"),
[0.75_f32, 0.75_f32, 0.75_f32],
sphere(0.0, 0.0, 0.0, 1.0)
)
}
fn default_regulators(key: ElementKey, assembly: &Assembly) -> Vec<Rc<dyn Regulator>> {
vec![Rc::new(HalfCurvatureRegulator::new(key, assembly))]
}
fn id(&self) -> &String {
&self.id
}
fn label(&self) -> &String {
&self.label
}
fn representation(&self) -> Signal<DVector<f64>> {
self.representation
}
fn regulators(&self) -> Signal<BTreeSet<RegulatorKey>> {
self.regulators
}
fn serial(&self) -> u64 {
self.serial
}
fn column_index(&self) -> Option<usize> {
self.column_index.get()
}
fn set_column_index(&self, index: usize) {
self.column_index.set(Some(index));
}
}
impl ProblemPoser for Sphere {
fn pose(&self, problem: &mut ConstraintProblem, _elts: &Slab<Rc<dyn Element>>) {
let index = self.column_index().expect(
format!("Sphere \"{}\" should be indexed before writing problem data", self.id).as_str()
);
problem.gram.push_sym(index, index, 1.0);
problem.guess.set_column(index, &self.representation.get_clone_untracked());
}
}
pub struct Point {
pub id: String,
pub label: String,
pub color: ElementColor,
pub representation: Signal<DVector<f64>>,
pub regulators: Signal<BTreeSet<RegulatorKey>>,
pub serial: u64,
column_index: Cell<Option<usize>>
}
impl Point {
const WEIGHT_COMPONENT: usize = 3;
pub fn new(
id: String,
label: String,
color: ElementColor,
representation: DVector<f64>
) -> Point {
Point {
id,
label,
color,
representation: create_signal(representation),
regulators: create_signal(BTreeSet::default()),
serial: Self::next_serial(),
column_index: None.into()
}
}
}
impl Element for Point {
fn default_id() -> String {
"point".to_string()
}
fn default(id: String, id_num: u64) -> Point {
Point::new(
id,
format!("Point {id_num}"),
[0.75_f32, 0.75_f32, 0.75_f32],
point(0.0, 0.0, 0.0)
)
}
fn id(&self) -> &String {
&self.id
}
fn label(&self) -> &String {
&self.label
}
fn representation(&self) -> Signal<DVector<f64>> {
self.representation
}
fn regulators(&self) -> Signal<BTreeSet<RegulatorKey>> {
self.regulators
}
fn serial(&self) -> u64 {
self.serial
}
fn column_index(&self) -> Option<usize> {
self.column_index.get()
}
fn set_column_index(&self, index: usize) {
self.column_index.set(Some(index));
}
}
impl ProblemPoser for Point {
fn pose(&self, problem: &mut ConstraintProblem, _elts: &Slab<Rc<dyn Element>>) {
let index = self.column_index().expect(
format!("Point \"{}\" should be indexed before writing problem data", self.id).as_str()
);
problem.gram.push_sym(index, index, 0.0);
problem.frozen.push(Point::WEIGHT_COMPONENT, index, 0.5);
problem.guess.set_column(index, &self.representation.get_clone_untracked());
}
}
pub trait Regulator: ProblemPoser + OutlineItem {
fn subjects(&self) -> Vec<ElementKey>;
fn measurement(&self) -> ReadSignal<f64>;
@ -168,7 +310,7 @@ impl InversiveDistanceRegulator {
pub fn new(subjects: [ElementKey; 2], assembly: &Assembly) -> InversiveDistanceRegulator {
let measurement = assembly.elements.map(
move |elts| {
let representations = subjects.map(|subj| elts[subj].representation);
let representations = subjects.map(|subj| elts[subj].representation());
representations[0].with(|rep_0|
representations[1].with(|rep_1|
rep_0.dot(&(&*Q * rep_1))
@ -198,11 +340,11 @@ impl Regulator for InversiveDistanceRegulator {
}
impl ProblemPoser for InversiveDistanceRegulator {
fn pose(&self, problem: &mut ConstraintProblem, elts: &Slab<Element>) {
fn pose(&self, problem: &mut ConstraintProblem, elts: &Slab<Rc<dyn Element>>) {
self.set_point.with_untracked(|set_pt| {
if let Some(val) = set_pt.value {
let [row, col] = self.subjects.map(
|subj| elts[subj].column_index.expect(
|subj| elts[subj].column_index().expect(
"Subjects should be indexed before inversive distance regulator writes problem data"
)
);
@ -221,8 +363,8 @@ pub struct HalfCurvatureRegulator {
impl HalfCurvatureRegulator {
pub fn new(subject: ElementKey, assembly: &Assembly) -> HalfCurvatureRegulator {
let measurement = assembly.elements.map(
move |elts| elts[subject].representation.with(
|rep| rep[Element::CURVATURE_COMPONENT]
move |elts| elts[subject].representation().with(
|rep| rep[Sphere::CURVATURE_COMPONENT]
)
);
@ -249,7 +391,7 @@ impl Regulator for HalfCurvatureRegulator {
match self.set_point.with(|set_pt| set_pt.value) {
Some(half_curv) => {
let representation = assembly.elements.with_untracked(
|elts| elts[self.subject].representation
|elts| elts[self.subject].representation()
);
representation.update(
|rep| change_half_curvature(rep, half_curv)
@ -262,13 +404,13 @@ impl Regulator for HalfCurvatureRegulator {
}
impl ProblemPoser for HalfCurvatureRegulator {
fn pose(&self, problem: &mut ConstraintProblem, elts: &Slab<Element>) {
fn pose(&self, problem: &mut ConstraintProblem, elts: &Slab<Rc<dyn Element>>) {
self.set_point.with_untracked(|set_pt| {
if let Some(val) = set_pt.value {
let col = elts[self.subject].column_index.expect(
let col = elts[self.subject].column_index().expect(
"Subject should be indexed before half-curvature regulator writes problem data"
);
problem.frozen.push(Element::CURVATURE_COMPONENT, col, val);
problem.frozen.push(Sphere::CURVATURE_COMPONENT, col, val);
}
});
}
@ -286,7 +428,7 @@ type AssemblyMotion<'a> = Vec<ElementMotion<'a>>;
#[derive(Clone)]
pub struct Assembly {
// elements and regulators
pub elements: Signal<Slab<Element>>,
pub elements: Signal<Slab<Rc<dyn Element>>>,
pub regulators: Signal<Slab<Rc<dyn Regulator>>>,
// solution variety tangent space. the basis vectors are stored in
@ -317,66 +459,61 @@ impl Assembly {
// --- inserting elements and regulators ---
// insert a sphere into the assembly without checking whether we already
// insert an element into the assembly without checking whether we already
// have an element with the same identifier. any element that does have the
// same identifier will get kicked out of the `elements_by_id` index
fn insert_sphere_unchecked(&self, elt: Element) -> ElementKey {
// insert the sphere
let id = elt.id.clone();
let key = self.elements.update(|elts| elts.insert(elt));
fn insert_element_unchecked<T: Element + 'static>(&self, elt: T) -> ElementKey {
// insert the element
let id = elt.id().clone();
let key = self.elements.update(|elts| elts.insert(Rc::new(elt)));
self.elements_by_id.update(|elts_by_id| elts_by_id.insert(id, key));
// regulate the sphere's curvature
self.insert_regulator(HalfCurvatureRegulator::new(key, &self));
// create and insert the element's default regulators
for reg in T::default_regulators(key, &self) {
self.insert_regulator(reg);
}
key
}
pub fn try_insert_sphere(&self, elt: Element) -> Option<ElementKey> {
pub fn try_insert_element(&self, elt: impl Element + 'static) -> Option<ElementKey> {
let can_insert = self.elements_by_id.with_untracked(
|elts_by_id| !elts_by_id.contains_key(&elt.id)
|elts_by_id| !elts_by_id.contains_key(elt.id())
);
if can_insert {
Some(self.insert_sphere_unchecked(elt))
Some(self.insert_element_unchecked(elt))
} else {
None
}
}
pub fn insert_new_sphere(&self) {
pub fn insert_element_default<T: Element + 'static>(&self) {
// find the next unused identifier in the default sequence
let default_id = T::default_id();
let mut id_num = 1;
let mut id = format!("sphere{}", id_num);
let mut id = format!("{default_id}{id_num}");
while self.elements_by_id.with_untracked(
|elts_by_id| elts_by_id.contains_key(&id)
) {
id_num += 1;
id = format!("sphere{}", id_num);
id = format!("{default_id}{id_num}");
}
// create and insert a sphere
let _ = self.insert_sphere_unchecked(
Element::new(
id,
format!("Sphere {}", id_num),
[0.75_f32, 0.75_f32, 0.75_f32],
sphere(0.0, 0.0, 0.0, 1.0)
)
);
// create and insert the default example of `T`
let _ = self.insert_element_unchecked(T::default(id, id_num));
}
pub fn insert_regulator<T: Regulator + 'static>(&self, regulator: T) {
pub fn insert_regulator(&self, regulator: Rc<dyn Regulator>) {
// add the regulator to the assembly's regulator list
let regulator_rc = Rc::new(regulator);
let key = self.regulators.update(
|regs| regs.insert(regulator_rc.clone())
|regs| regs.insert(regulator.clone())
);
// add the regulator to each subject's regulator list
let subjects = regulator_rc.subjects();
let subjects = regulator.subjects();
let subject_regulators: Vec<_> = self.elements.with_untracked(
|elts| subjects.into_iter().map(
|subj| elts[subj].regulators
|subj| elts[subj].regulators()
).collect()
);
for regulators in subject_regulators {
@ -390,10 +527,10 @@ impl Assembly {
/* DEBUG */
// log the regulator update
console::log_1(&JsValue::from(
format!("Updated regulator with subjects {:?}", regulator_rc.subjects())
format!("Updated regulator with subjects {:?}", regulator.subjects())
));
if regulator_rc.try_activate(&self_for_effect) {
if regulator.try_activate(&self_for_effect) {
self_for_effect.realize();
}
});
@ -427,7 +564,7 @@ impl Assembly {
// index the elements
self.elements.update_silent(|elts| {
for (index, (_, elt)) in elts.into_iter().enumerate() {
elt.column_index = Some(index);
elt.set_column_index(index);
}
});
@ -482,8 +619,8 @@ impl Assembly {
if success {
// read out the solution
for (_, elt) in self.elements.get_clone_untracked() {
elt.representation.update(
|rep| rep.set_column(0, &config.column(elt.column_index.unwrap()))
elt.representation().update(
|rep| rep.set_column(0, &config.column(elt.column_index().unwrap()))
);
}
@ -521,8 +658,8 @@ impl Assembly {
let mut next_column_index = realized_dim;
for elt_motion in motion.iter() {
let moving_elt = &mut elts[elt_motion.key];
if moving_elt.column_index.is_none() {
moving_elt.column_index = Some(next_column_index);
if moving_elt.column_index().is_none() {
moving_elt.set_column_index(next_column_index);
next_column_index += 1;
}
}
@ -539,7 +676,7 @@ impl Assembly {
// we can unwrap the column index because we know that every moving
// element has one at this point
let column_index = self.elements.with_untracked(
|elts| elts[elt_motion.key].column_index.unwrap()
|elts| elts[elt_motion.key].column_index().unwrap()
);
if column_index < realized_dim {
@ -555,7 +692,7 @@ impl Assembly {
let mut target_column = motion_proj.column_mut(column_index);
let unif_to_std = self.elements.with_untracked(
|elts| {
elts[elt_motion.key].representation.with_untracked(
elts[elt_motion.key].representation().with_untracked(
|rep| local_unif_to_std(rep.as_view())
)
}
@ -567,26 +704,27 @@ impl Assembly {
// step the assembly along the deformation. this changes the elements'
// normalizations, so we restore those afterward
/* KLUDGE */
// since our test assemblies only include spheres, we assume that every
// element is on the 1 mass shell
// for now, we only restore the normalizations of spheres
for (_, elt) in self.elements.get_clone_untracked() {
elt.representation.update_silent(|rep| {
match elt.column_index {
elt.representation().update_silent(|rep| {
match elt.column_index() {
Some(column_index) => {
// step the assembly along the deformation
*rep += motion_proj.column(column_index);
// restore normalization by contracting toward the last
// coordinate axis
let q_sp = rep.fixed_rows::<3>(0).norm_squared();
let half_q_lt = -2.0 * rep[3] * rep[4];
let half_q_lt_sq = half_q_lt * half_q_lt;
let scaling = half_q_lt + (q_sp + half_q_lt_sq).sqrt();
rep.fixed_rows_mut::<4>(0).scale_mut(1.0 / scaling);
if elt.type_id() == TypeId::of::<Sphere>() {
// restore normalization by contracting toward the
// last coordinate axis
let q_sp = rep.fixed_rows::<3>(0).norm_squared();
let half_q_lt = -2.0 * rep[3] * rep[4];
let half_q_lt_sq = half_q_lt * half_q_lt;
let scaling = half_q_lt + (q_sp + half_q_lt_sq).sqrt();
rep.fixed_rows_mut::<4>(0).scale_mut(1.0 / scaling);
}
},
None => {
console::log_1(&JsValue::from(
format!("No velocity to unpack for fresh element \"{}\"", elt.id)
format!("No velocity to unpack for fresh element \"{}\"", elt.id())
))
}
};
@ -602,20 +740,14 @@ impl Assembly {
#[cfg(test)]
mod tests {
use crate::engine;
use super::*;
#[test]
#[should_panic(expected = "Element \"sphere\" should be indexed before writing problem data")]
#[should_panic(expected = "Sphere \"sphere\" should be indexed before writing problem data")]
fn unindexed_element_test() {
let _ = create_root(|| {
Element::new(
"sphere".to_string(),
"Sphere".to_string(),
[1.0_f32, 1.0_f32, 1.0_f32],
engine::sphere(0.0, 0.0, 0.0, 1.0)
).pose(&mut ConstraintProblem::new(1), &Slab::new());
let elt = Sphere::default("sphere".to_string(), 0);
elt.pose(&mut ConstraintProblem::new(1), &Slab::new());
});
}
@ -623,18 +755,13 @@ mod tests {
#[should_panic(expected = "Subjects should be indexed before inversive distance regulator writes problem data")]
fn unindexed_subject_test_inversive_distance() {
let _ = create_root(|| {
let mut elts = Slab::new();
let mut elts = Slab::<Rc<dyn Element>>::new();
let subjects = [0, 1].map(|k| {
elts.insert(
Element::new(
format!("sphere{k}"),
format!("Sphere {k}"),
[1.0_f32, 1.0_f32, 1.0_f32],
engine::sphere(0.0, 0.0, 0.0, 1.0)
)
Rc::new(Sphere::default(format!("sphere{k}"), k))
)
});
elts[subjects[0]].column_index = Some(0);
elts[subjects[0]].set_column_index(0);
InversiveDistanceRegulator {
subjects: subjects,
measurement: create_memo(|| 0.0),