From 4fd79b9e47bce0a829aff4e3178374af217fa487 Mon Sep 17 00:00:00 2001 From: Aaron Fenyes Date: Tue, 17 Dec 2024 18:21:53 -0800 Subject: [PATCH 1/6] Add structures for element and assembly motions --- app-proto/src/assembly.rs | 15 +++++++++++---- app-proto/src/display.rs | 11 +++++++++-- 2 files changed, 20 insertions(+), 6 deletions(-) diff --git a/app-proto/src/assembly.rs b/app-proto/src/assembly.rs index 1684716..07b0aba 100644 --- a/app-proto/src/assembly.rs +++ b/app-proto/src/assembly.rs @@ -119,6 +119,13 @@ pub struct Constraint { pub active: Signal } +pub struct ElementMotion<'a> { + pub key: ElementKey, + pub velocity: DVectorView<'a, f64> +} + +type AssemblyMotion<'a> = Vec>; + // a complete, view-independent description of an assembly #[derive(Clone)] pub struct Assembly { @@ -290,7 +297,7 @@ impl Assembly { // --- deformation --- - pub fn deform(&self, element_motions: Vec<(ElementKey, DVectorView)>) { + pub fn deform(&self, motion: AssemblyMotion) { /* KLUDGE */ // when the tangent space is zero, we currently need to avoid calling // its `proj` method, because it will panic rather than returning zero. @@ -307,9 +314,9 @@ impl Assembly { // project the element motions onto the tangent space of the solution // variety, and sum them to get a deformation of the whole assembly - for (elt_key, elt_motion) in element_motions { - let column_index = self.elements.with(|elts| elts[elt_key].column_index); - motion_proj += self.tangent.with(|tan| tan.proj(&elt_motion, column_index)); + for elt_motion in motion { + let column_index = self.elements.with(|elts| elts[elt_motion.key].column_index); + motion_proj += self.tangent.with(|tan| tan.proj(&elt_motion.velocity, column_index)); } // step each element along the mass shell geodesic that matches its diff --git a/app-proto/src/display.rs b/app-proto/src/display.rs index 673985c..d8ee23c 100644 --- a/app-proto/src/display.rs +++ b/app-proto/src/display.rs @@ -14,7 +14,7 @@ use web_sys::{ wasm_bindgen::{JsCast, JsValue} }; -use crate::{AppState, assembly::ElementKey}; +use crate::{AppState, assembly::{ElementKey, ElementMotion}}; fn compile_shader( context: &WebGl2RenderingContext, @@ -341,7 +341,14 @@ pub fn Display() -> View { ]) }; let elt_motion: DVector = time_step * TRANSLATION_SPEED * vel_field * rep; - assembly_for_raf.deform(vec![(sel, elt_motion.as_view())]); + assembly_for_raf.deform( + vec![ + ElementMotion { + key: sel, + velocity: elt_motion.as_view() + } + ] + ); scene_changed.set(true); } } From 971a7ca7e279827e7a7b6129200899a3a20aac8a Mon Sep 17 00:00:00 2001 From: Aaron Fenyes Date: Tue, 17 Dec 2024 21:24:38 -0800 Subject: [PATCH 2/6] Check tangent space sync when deforming Only give elements column indices once they've actually been through a realization. Ignore motions of elements that haven't been through a realization. Get the dimensions of the projected motion matrix from the saved tangent space, not the current number of elements. --- app-proto/src/assembly.rs | 53 ++++++++++++++++++++++++++++----------- app-proto/src/engine.rs | 31 +++++++++++++++-------- 2 files changed, 59 insertions(+), 25 deletions(-) diff --git a/app-proto/src/assembly.rs b/app-proto/src/assembly.rs index 07b0aba..ec480ac 100644 --- a/app-proto/src/assembly.rs +++ b/app-proto/src/assembly.rs @@ -33,8 +33,9 @@ pub struct Element { pub serial: u64, // the configuration matrix column index that was assigned to this element - // last time the assembly was realized - column_index: usize + // last time the assembly was realized, or `None` if the assembly has never + // been realized with this element in it + column_index: Option } impl Element { @@ -62,7 +63,7 @@ impl Element { representation: create_signal(representation), constraints: create_signal(BTreeSet::default()), serial: serial, - column_index: 0 + column_index: None } } @@ -145,7 +146,7 @@ impl Assembly { Assembly { elements: create_signal(Slab::new()), constraints: create_signal(Slab::new()), - tangent: create_signal(ConfigSubspace::zero()), + tangent: create_signal(ConfigSubspace::zero(0)), elements_by_id: create_signal(FxHashMap::default()) } } @@ -216,7 +217,7 @@ impl Assembly { // index the elements self.elements.update_silent(|elts| { for (index, (_, elt)) in elts.into_iter().enumerate() { - elt.column_index = index; + elt.column_index = Some(index); } }); @@ -228,8 +229,8 @@ impl Assembly { for (_, cst) in csts { if cst.active.get_untracked() && cst.lorentz_prod_valid.get_untracked() { let subjects = cst.subjects; - let row = elts[subjects.0].column_index; - let col = elts[subjects.1].column_index; + let row = elts[subjects.0].column_index.unwrap(); + let col = elts[subjects.1].column_index.unwrap(); gram_to_be.push_sym(row, col, cst.lorentz_prod.get_untracked()); } } @@ -239,7 +240,7 @@ impl Assembly { // Gram matrix let mut guess_to_be = DMatrix::::zeros(5, elts.len()); for (_, elt) in elts { - let index = elt.column_index; + let index = elt.column_index.unwrap(); gram_to_be.push_sym(index, index, 1.0); guess_to_be.set_column(index, &elt.representation.get_clone_untracked()); } @@ -286,7 +287,7 @@ impl Assembly { // 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)) + |rep| rep.set_column(0, &config.column(elt.column_index.unwrap())) ); } @@ -309,14 +310,27 @@ impl Assembly { } const ELEMENT_DIM: usize = 5; - let assembly_dim = self.elements.with(|elts| elts.len()); + let assembly_dim = self.tangent.with(|tan| tan.assembly_dim()); let mut motion_proj = DMatrix::zeros(ELEMENT_DIM, assembly_dim); // project the element motions onto the tangent space of the solution // variety, and sum them to get a deformation of the whole assembly for elt_motion in motion { - let column_index = self.elements.with(|elts| elts[elt_motion.key].column_index); - motion_proj += self.tangent.with(|tan| tan.proj(&elt_motion.velocity, column_index)); + match self.elements.with_untracked(|elts| elts[elt_motion.key].column_index) { + Some(column_index) => { + motion_proj += self.tangent.with( + |tan| tan.proj(&elt_motion.velocity, column_index) + ) + }, + None => { + console::log_1(&JsValue::from( + format!( + "Ignoring motion of fresh element \"{}\"", + self.elements.with_untracked(|elts| elts[elt_motion.key].id.clone()) + ) + )) + } + } } // step each element along the mass shell geodesic that matches its @@ -326,9 +340,18 @@ impl Assembly { // element is on the 1 mass shell for (_, elt) in self.elements.get_clone_untracked() { elt.representation.update_silent(|rep| { - let rep_next = &*rep + motion_proj.column(elt.column_index); - let normalizer = rep_next.dot(&(&*Q * &rep_next)); - rep.set_column(0, &(rep_next / normalizer)); + match elt.column_index { + Some(column_index) => { + let rep_next = &*rep + motion_proj.column(column_index); + let normalizer = rep_next.dot(&(&*Q * &rep_next)); + rep.set_column(0, &(rep_next / normalizer)); + }, + None => { + console::log_1(&JsValue::from( + format!("No velocity to unpack for fresh element \"{}\"", elt.id) + )) + } + }; }); } diff --git a/app-proto/src/engine.rs b/app-proto/src/engine.rs index 2f06035..36998bd 100644 --- a/app-proto/src/engine.rs +++ b/app-proto/src/engine.rs @@ -88,11 +88,17 @@ impl PartialMatrix { // --- configuration subspaces --- #[derive(Clone)] -pub struct ConfigSubspace(Vec>); +pub struct ConfigSubspace { + assembly_dim: usize, + basis: Vec> +} impl ConfigSubspace { - pub fn zero() -> ConfigSubspace { - ConfigSubspace(Vec::new()) + pub fn zero(assembly_dim: usize) -> ConfigSubspace { + ConfigSubspace { + assembly_dim: assembly_dim, + basis: Vec::new() + } } // approximate the kernel of a symmetric endomorphism of the configuration @@ -119,12 +125,18 @@ impl ConfigSubspace { format!("Eigenvalues used to find kernel: {}", eig.eigenvalues) )); - ConfigSubspace(basis.collect()) + ConfigSubspace { + assembly_dim: assembly_dim, + basis: basis.collect() + } } pub fn dim(&self) -> usize { - let ConfigSubspace(basis) = self; - basis.len() + self.basis.len() + } + + pub fn assembly_dim(&self) -> usize { + self.assembly_dim } // find the projection onto this subspace of the motion where the element @@ -133,8 +145,7 @@ impl ConfigSubspace { // for the zero subspace, this method's behavior doesn't match its name: it // panics rather than returning zero pub fn proj(&self, v: &DVectorView, column_index: usize) -> DMatrix { - let ConfigSubspace(basis) = self; - basis.into_iter().map( + self.basis.iter().map( |b| b.column(column_index).dot(&v) * b ).sum() } @@ -325,14 +336,14 @@ pub fn realize_gram( state = better_state; history.backoff_steps.push(backoff_steps); }, - None => return (state.config, ConfigSubspace::zero(), false, history) + None => return (state.config, ConfigSubspace::zero(assembly_dim), false, history) }; } let success = state.loss < tol; let tangent = if success { ConfigSubspace::symmetric_kernel(hess, assembly_dim) } else { - ConfigSubspace::zero() + ConfigSubspace::zero(assembly_dim) }; (state.config, tangent, success, history) } From dc067976eb4efbf9fd5dbc663bd97850c2a30121 Mon Sep 17 00:00:00 2001 From: Aaron Fenyes Date: Wed, 18 Dec 2024 00:25:15 -0800 Subject: [PATCH 3/6] Implement projection onto the zero subspace --- app-proto/src/engine.rs | 14 ++++++++------ 1 file changed, 8 insertions(+), 6 deletions(-) diff --git a/app-proto/src/engine.rs b/app-proto/src/engine.rs index 36998bd..e96ea52 100644 --- a/app-proto/src/engine.rs +++ b/app-proto/src/engine.rs @@ -141,13 +141,15 @@ impl ConfigSubspace { // find the projection onto this subspace of the motion where the element // with the given column index has velocity `v` - /* TO DO */ - // for the zero subspace, this method's behavior doesn't match its name: it - // panics rather than returning zero pub fn proj(&self, v: &DVectorView, column_index: usize) -> DMatrix { - self.basis.iter().map( - |b| b.column(column_index).dot(&v) * b - ).sum() + if self.dim() == 0 { + const ELEMENT_DIM: usize = 5; + DMatrix::zeros(ELEMENT_DIM, self.assembly_dim) + } else { + self.basis.iter().map( + |b| b.column(column_index).dot(&v) * b + ).sum() + } } } From 967daa595dd3b3ffd7564d6848319e4fb13f6ded Mon Sep 17 00:00:00 2001 From: Aaron Fenyes Date: Wed, 18 Dec 2024 00:34:25 -0800 Subject: [PATCH 4/6] Deform fresh elements too Implement deformation of elements that haven't gone through realization. --- app-proto/src/assembly.rs | 51 +++++++++++++++++++++------------------ 1 file changed, 27 insertions(+), 24 deletions(-) diff --git a/app-proto/src/assembly.rs b/app-proto/src/assembly.rs index ec480ac..fa0a48c 100644 --- a/app-proto/src/assembly.rs +++ b/app-proto/src/assembly.rs @@ -33,8 +33,8 @@ pub struct 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 assembly has never - // been realized with this element in it + // last time the assembly was realized, or `None` if the element has never + // been through a realization column_index: Option } @@ -160,13 +160,6 @@ impl Assembly { let id = elt.id.clone(); let key = self.elements.update(|elts| elts.insert(elt)); self.elements_by_id.update(|elts_by_id| elts_by_id.insert(id, key)); - - // realize to update the tangent space - /* KLUDGE */ - // since the newly inserted element is unconstrained, we should be able - // to update the tangent space without recomputing the Hessian and its - // eigendecomposition - self.realize(); } pub fn try_insert_element(&self, elt: Element) -> bool { @@ -300,35 +293,43 @@ impl Assembly { pub fn deform(&self, motion: AssemblyMotion) { /* KLUDGE */ - // when the tangent space is zero, we currently need to avoid calling - // its `proj` method, because it will panic rather than returning zero. - // in the future, we'll want a more intentionally designed system for - // handling this case - if self.tangent.with(|tan| tan.dim() <= 0) { + // when the tangent space is zero, deformation won't do anything, but + // the attempt to deform should be registered in the UI. this console + // message will do for now + if self.tangent.with(|tan| tan.dim() <= 0 && tan.assembly_dim() > 0) { console::log_1(&JsValue::from("The assembly is rigid")); - return; } const ELEMENT_DIM: usize = 5; - let assembly_dim = self.tangent.with(|tan| tan.assembly_dim()); + let assembly_dim = self.elements.with(|elts| elts.len()); let mut motion_proj = DMatrix::zeros(ELEMENT_DIM, assembly_dim); + let tan_assembly_dim = self.tangent.with(|tan| tan.assembly_dim()); // project the element motions onto the tangent space of the solution // variety, and sum them to get a deformation of the whole assembly + let mut next_column_index = tan_assembly_dim; for elt_motion in motion { match self.elements.with_untracked(|elts| elts[elt_motion.key].column_index) { Some(column_index) => { - motion_proj += self.tangent.with( + let mut target_columns = motion_proj.columns_mut(0, tan_assembly_dim); + target_columns += self.tangent.with( |tan| tan.proj(&elt_motion.velocity, column_index) ) }, None => { - console::log_1(&JsValue::from( - format!( - "Ignoring motion of fresh element \"{}\"", - self.elements.with_untracked(|elts| elts[elt_motion.key].id.clone()) - ) - )) + // give the element a column index, even though it's never + // been through a realization, and add its motion to that + // column of the projected motion. this puts the element + // column indices out of sync with the saved tangent space, + // but that's okay: the realization we do to get back onto + // the solution variety will re-index the elements and + // recompute the tangent space + let mut target_column = motion_proj.column_mut(next_column_index); + target_column += elt_motion.velocity; + self.elements.update_silent( + |elts| elts[elt_motion.key].column_index = Some(next_column_index) + ); + next_column_index += 1; } } } @@ -355,7 +356,9 @@ impl Assembly { }); } - // bring the configuration back onto the solution variety + // bring the configuration back onto the solution variety. this also + // gets the elements' column indices and the saved tangent space back in + // sync self.realize(); } } \ No newline at end of file From e2c5ba0fc7520eb3da54a7f919ed65fb937785ed Mon Sep 17 00:00:00 2001 From: Aaron Fenyes Date: Wed, 18 Dec 2024 09:49:14 -0800 Subject: [PATCH 5/6] Set out invariants for column indices MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit This should make it safe to use the elements' column indices outside the realization method—for unpacking tangent vectors, at least. --- app-proto/src/assembly.rs | 27 +++++++++++++++++++++------ 1 file changed, 21 insertions(+), 6 deletions(-) diff --git a/app-proto/src/assembly.rs b/app-proto/src/assembly.rs index fa0a48c..b4b5329 100644 --- a/app-proto/src/assembly.rs +++ b/app-proto/src/assembly.rs @@ -134,7 +134,16 @@ pub struct Assembly { pub elements: Signal>, pub constraints: Signal>, - // solution variety tangent space + // solution variety tangent space. the basis vectors are stored in + // configuration matrix format, ordered according to the elements' column + // indices. when you realize the assembly, every element that's present + // during realization gets a column index and is reflected in the tangent + // space. since the methods in this module never assign column indices + // without later realizing the assembly, we get the following invariant: + // + // (1) if an element has a column index, its tangent motions can be found + // in that column of the tangent space basis matrices + // pub tangent: Signal, // indexing @@ -291,6 +300,14 @@ impl Assembly { // --- deformation --- + // project the given motion to the tangent space of the solution variety and + // move the assembly along it. the implementation is based on invariant (1) + // from above and the following additional invariant: + // + // (2) if an element is affected by a constraint, it has a column index + // + // we have this invariant because the assembly gets realized each time you + // add a constraint pub fn deform(&self, motion: AssemblyMotion) { /* KLUDGE */ // when the tangent space is zero, deformation won't do anything, but @@ -319,11 +336,9 @@ impl Assembly { None => { // give the element a column index, even though it's never // been through a realization, and add its motion to that - // column of the projected motion. this puts the element - // column indices out of sync with the saved tangent space, - // but that's okay: the realization we do to get back onto - // the solution variety will re-index the elements and - // recompute the tangent space + // column of the projected motion. this temporarily breaks + // invariant (1), but the invariant will be restored when we + // realize the assembly at the end of the deformation let mut target_column = motion_proj.column_mut(next_column_index); target_column += elt_motion.velocity; self.elements.update_silent( From 6df0e855cf3a7ca144c2388965215ad6888a8bcb Mon Sep 17 00:00:00 2001 From: Aaron Fenyes Date: Wed, 18 Dec 2024 11:43:54 -0800 Subject: [PATCH 6/6] Make the deformation matrix just the right size Also, correct the check for whether an element had a column index when we started. The previous revision would've gotten the wrong answer for an element without a column index that appeared more than once in the motion. --- app-proto/src/assembly.rs | 68 +++++++++++++++++++++++---------------- 1 file changed, 41 insertions(+), 27 deletions(-) diff --git a/app-proto/src/assembly.rs b/app-proto/src/assembly.rs index b4b5329..278a8f9 100644 --- a/app-proto/src/assembly.rs +++ b/app-proto/src/assembly.rs @@ -317,36 +317,50 @@ impl Assembly { console::log_1(&JsValue::from("The assembly is rigid")); } - const ELEMENT_DIM: usize = 5; - let assembly_dim = self.elements.with(|elts| elts.len()); - let mut motion_proj = DMatrix::zeros(ELEMENT_DIM, assembly_dim); - let tan_assembly_dim = self.tangent.with(|tan| tan.assembly_dim()); - - // project the element motions onto the tangent space of the solution - // variety, and sum them to get a deformation of the whole assembly - let mut next_column_index = tan_assembly_dim; - for elt_motion in motion { - match self.elements.with_untracked(|elts| elts[elt_motion.key].column_index) { - Some(column_index) => { - let mut target_columns = motion_proj.columns_mut(0, tan_assembly_dim); - target_columns += self.tangent.with( - |tan| tan.proj(&elt_motion.velocity, column_index) - ) - }, - None => { - // give the element a column index, even though it's never - // been through a realization, and add its motion to that - // column of the projected motion. this temporarily breaks - // invariant (1), but the invariant will be restored when we - // realize the assembly at the end of the deformation - let mut target_column = motion_proj.column_mut(next_column_index); - target_column += elt_motion.velocity; - self.elements.update_silent( - |elts| elts[elt_motion.key].column_index = Some(next_column_index) - ); + // give a column index to each moving element that doesn't have one yet. + // this temporarily breaks invariant (1), but the invariant will be + // restored when we realize the assembly at the end of the deformation. + // in the process, we find out how many matrix columns we'll need to + // hold the deformation + let realized_dim = self.tangent.with(|tan| tan.assembly_dim()); + let motion_dim = self.elements.update_silent(|elts| { + 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); next_column_index += 1; } } + next_column_index + }); + + // project the element motions onto the tangent space of the solution + // variety and sum them to get a deformation of the whole assembly. the + // matrix `motion_proj` that holds the deformation has extra columns for + // any moving elements that aren't reflected in the saved tangent space + const ELEMENT_DIM: usize = 5; + let mut motion_proj = DMatrix::zeros(ELEMENT_DIM, motion_dim); + for elt_motion in motion { + // 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() + ); + + if column_index < realized_dim { + // this element had a column index when we started, so by + // invariant (1), it's reflected in the tangent space + let mut target_columns = motion_proj.columns_mut(0, realized_dim); + target_columns += self.tangent.with( + |tan| tan.proj(&elt_motion.velocity, column_index) + ); + } else { + // this element didn't have a column index when we started, so + // by invariant (2), it's unconstrained + let mut target_column = motion_proj.column_mut(column_index); + target_column += elt_motion.velocity; + } } // step each element along the mass shell geodesic that matches its