diff --git a/app-proto/examples/irisawa-hexlet.rs b/app-proto/examples/irisawa-hexlet.rs index fc14f91..2bc94c0 100644 --- a/app-proto/examples/irisawa-hexlet.rs +++ b/app-proto/examples/irisawa-hexlet.rs @@ -2,7 +2,7 @@ use dyna3::engine::{Q, irisawa::realize_irisawa_hexlet}; fn main() { const SCALED_TOL: f64 = 1.0e-12; - let (config, success, history) = realize_irisawa_hexlet(SCALED_TOL); + let (config, _, success, history) = realize_irisawa_hexlet(SCALED_TOL); print!("\nCompleted Gram matrix:{}", config.tr_mul(&*Q) * &config); if success { println!("Target accuracy achieved!"); diff --git a/app-proto/examples/point-on-sphere.rs b/app-proto/examples/point-on-sphere.rs index 0e4765a..13040e5 100644 --- a/app-proto/examples/point-on-sphere.rs +++ b/app-proto/examples/point-on-sphere.rs @@ -18,7 +18,7 @@ fn main() { ]); let frozen = [(3, 0)]; println!(); - let (config, success, history) = realize_gram( + let (config, _, success, history) = realize_gram( &gram, guess, &frozen, 1.0e-12, 0.5, 0.9, 1.1, 200, 110 ); diff --git a/app-proto/examples/three-spheres.rs b/app-proto/examples/three-spheres.rs index d348b18..19acfd1 100644 --- a/app-proto/examples/three-spheres.rs +++ b/app-proto/examples/three-spheres.rs @@ -21,7 +21,7 @@ fn main() { ]) }; println!(); - let (config, success, history) = realize_gram( + let (config, _, success, history) = realize_gram( &gram, guess, &[], 1.0e-12, 0.5, 0.9, 1.1, 200, 110 ); diff --git a/app-proto/src/assembly.rs b/app-proto/src/assembly.rs index fb5bbf7..278a8f9 100644 --- a/app-proto/src/assembly.rs +++ b/app-proto/src/assembly.rs @@ -1,11 +1,11 @@ -use nalgebra::{DMatrix, DVector, Vector3}; +use nalgebra::{DMatrix, DVector, DVectorView, Vector3}; use rustc_hash::FxHashMap; use slab::Slab; use std::{collections::BTreeSet, sync::atomic::{AtomicU64, Ordering}}; use sycamore::prelude::*; use web_sys::{console, wasm_bindgen::JsValue}; /* DEBUG */ -use crate::engine::{realize_gram, PartialMatrix}; +use crate::engine::{realize_gram, ConfigSubspace, PartialMatrix, Q}; // the types of the keys we use to access an assembly's elements and constraints pub type ElementKey = usize; @@ -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 element has never + // been through a realization + 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 } } @@ -109,7 +110,6 @@ impl Element { } } } - #[derive(Clone)] pub struct Constraint { @@ -120,6 +120,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 { @@ -127,6 +134,18 @@ pub struct Assembly { pub elements: Signal>, pub constraints: Signal>, + // 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 pub elements_by_id: Signal> } @@ -136,6 +155,7 @@ impl Assembly { Assembly { elements: create_signal(Slab::new()), constraints: create_signal(Slab::new()), + tangent: create_signal(ConfigSubspace::zero(0)), elements_by_id: create_signal(FxHashMap::default()) } } @@ -199,7 +219,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); } }); @@ -211,8 +231,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()); } } @@ -222,7 +242,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()); } @@ -247,7 +267,7 @@ impl Assembly { } // look for a configuration with the given Gram matrix - let (config, success, history) = realize_gram( + let (config, tangent, success, history) = realize_gram( &gram, guess, &[], 1.0e-12, 0.5, 0.9, 1.1, 200, 110 ); @@ -263,14 +283,111 @@ impl Assembly { )); console::log_2(&JsValue::from("Steps:"), &JsValue::from(history.scaled_loss.len() - 1)); console::log_2(&JsValue::from("Loss:"), &JsValue::from(*history.scaled_loss.last().unwrap())); + console::log_2(&JsValue::from("Tangent dimension:"), &JsValue::from(tangent.dim())); 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)) + |rep| rep.set_column(0, &config.column(elt.column_index.unwrap())) ); } + + // save the tangent space + self.tangent.set_silent(tangent); } } + + // --- 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 + // 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")); + } + + // 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 + // velocity in the deformation found above + /* KLUDGE */ + // since our test assemblies only include spheres, we assume that every + // element is on the 1 mass shell + for (_, elt) in self.elements.get_clone_untracked() { + elt.representation.update_silent(|rep| { + 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) + )) + } + }; + }); + } + + // 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 diff --git a/app-proto/src/display.rs b/app-proto/src/display.rs index c39e575..d8ee23c 100644 --- a/app-proto/src/display.rs +++ b/app-proto/src/display.rs @@ -1,5 +1,5 @@ use core::array; -use nalgebra::{DMatrix, Rotation3, Vector3}; +use nalgebra::{DMatrix, DVector, Rotation3, Vector3}; use sycamore::{prelude::*, motion::create_raf}; use web_sys::{ console, @@ -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, @@ -123,6 +123,14 @@ pub fn Display() -> View { let zoom_out = create_signal(0.0); let turntable = create_signal(false); /* BENCHMARKING */ + // manipulation + let translate_neg_x = create_signal(0.0); + let translate_pos_x = create_signal(0.0); + let translate_neg_y = create_signal(0.0); + let translate_pos_y = create_signal(0.0); + let translate_neg_z = create_signal(0.0); + let translate_pos_z = create_signal(0.0); + // change listener let scene_changed = create_signal(true); create_effect(move || { @@ -141,6 +149,7 @@ pub fn Display() -> View { let mut frames_since_last_sample = 0; let mean_frame_interval = create_signal(0.0); + let assembly_for_raf = state.assembly.clone(); on_mount(move || { // timing let mut last_time = 0.0; @@ -153,6 +162,9 @@ pub fn Display() -> View { let mut rotation = DMatrix::::identity(5, 5); let mut location_z: f64 = 5.0; + // manipulation + const TRANSLATION_SPEED: f64 = 0.15; // in length units per second + // display parameters const OPACITY: f32 = 0.5; /* SCAFFOLDING */ const HIGHLIGHT: f32 = 0.2; /* SCAFFOLDING */ @@ -273,6 +285,14 @@ pub fn Display() -> View { let zoom_out_val = zoom_out.get(); let turntable_val = turntable.get(); /* BENCHMARKING */ + // get the manipulation state + let translate_neg_x_val = translate_neg_x.get(); + let translate_pos_x_val = translate_pos_x.get(); + let translate_neg_y_val = translate_neg_y.get(); + let translate_pos_y_val = translate_pos_y.get(); + let translate_neg_z_val = translate_neg_z.get(); + let translate_pos_z_val = translate_pos_z.get(); + // update the assembly's orientation let ang_vel = { let pitch = pitch_up_val - pitch_down_val; @@ -298,6 +318,41 @@ pub fn Display() -> View { let zoom = zoom_out_val - zoom_in_val; location_z *= (time_step * ZOOM_SPEED * zoom).exp(); + // manipulate the assembly + if state.selection.with(|sel| sel.len() == 1) { + let sel = state.selection.with( + |sel| *sel.into_iter().next().unwrap() + ); + let rep = state.assembly.elements.with_untracked( + |elts| elts[sel].representation.get_clone_untracked() + ); + let translate_x = translate_pos_x_val - translate_neg_x_val; + let translate_y = translate_pos_y_val - translate_neg_y_val; + let translate_z = translate_pos_z_val - translate_neg_z_val; + if translate_x != 0.0 || translate_y != 0.0 || translate_z != 0.0 { + let vel_field = { + let u = Vector3::new(translate_x, translate_y, translate_z).normalize(); + DMatrix::from_column_slice(5, 5, &[ + 0.0, 0.0, 0.0, 0.0, u[0], + 0.0, 0.0, 0.0, 0.0, u[1], + 0.0, 0.0, 0.0, 0.0, u[2], + 2.0*u[0], 2.0*u[1], 2.0*u[2], 0.0, 0.0, + 0.0, 0.0, 0.0, 0.0, 0.0 + ]) + }; + let elt_motion: DVector = time_step * TRANSLATION_SPEED * vel_field * rep; + assembly_for_raf.deform( + vec![ + ElementMotion { + key: sel, + velocity: elt_motion.as_view() + } + ] + ); + scene_changed.set(true); + } + } + if scene_changed.get() { /* INSTRUMENTS */ // measure mean frame interval @@ -416,7 +471,7 @@ pub fn Display() -> View { start_animation_loop(); }); - let set_nav_signal = move |event: KeyboardEvent, value: f64| { + let set_nav_signal = move |event: &KeyboardEvent, value: f64| { let mut navigating = true; let shift = event.shift_key(); match event.key().as_str() { @@ -436,6 +491,23 @@ pub fn Display() -> View { } }; + let set_manip_signal = move |event: &KeyboardEvent, value: f64| { + let mut manipulating = true; + let shift = event.shift_key(); + match event.key().as_str() { + "d" | "D" => translate_pos_x.set(value), + "a" | "A" => translate_neg_x.set(value), + "w" | "W" if shift => translate_neg_z.set(value), + "s" | "S" if shift => translate_pos_z.set(value), + "w" | "W" => translate_pos_y.set(value), + "s" | "S" => translate_neg_y.set(value), + _ => manipulating = false + }; + if manipulating { + event.prevent_default(); + } + }; + view! { /* TO DO */ // switch back to integer-valued parameters when that becomes possible @@ -447,6 +519,7 @@ pub fn Display() -> View { tabindex="0", on:keydown=move |event: KeyboardEvent| { if event.key() == "Shift" { + // swap navigation inputs roll_cw.set(yaw_right.get()); roll_ccw.set(yaw_left.get()); zoom_in.set(pitch_up.get()); @@ -455,16 +528,24 @@ pub fn Display() -> View { yaw_left.set(0.0); pitch_up.set(0.0); pitch_down.set(0.0); + + // swap manipulation inputs + translate_pos_z.set(translate_neg_y.get()); + translate_neg_z.set(translate_pos_y.get()); + translate_pos_y.set(0.0); + translate_neg_y.set(0.0); } else { if event.key() == "Enter" { /* BENCHMARKING */ turntable.set_fn(|turn| !turn); scene_changed.set(true); } - set_nav_signal(event, 1.0); + set_nav_signal(&event, 1.0); + set_manip_signal(&event, 1.0); } }, on:keyup=move |event: KeyboardEvent| { if event.key() == "Shift" { + // swap navigation inputs yaw_right.set(roll_cw.get()); yaw_left.set(roll_ccw.get()); pitch_up.set(zoom_in.get()); @@ -473,8 +554,15 @@ pub fn Display() -> View { roll_ccw.set(0.0); zoom_in.set(0.0); zoom_out.set(0.0); + + // swap manipulation inputs + translate_pos_y.set(translate_neg_z.get()); + translate_neg_y.set(translate_pos_z.get()); + translate_pos_z.set(0.0); + translate_neg_z.set(0.0); } else { - set_nav_signal(event, 0.0); + set_nav_signal(&event, 0.0); + set_manip_signal(&event, 0.0); } }, on:blur=move |_| { diff --git a/app-proto/src/engine.rs b/app-proto/src/engine.rs index 285a9c6..6a0202e 100644 --- a/app-proto/src/engine.rs +++ b/app-proto/src/engine.rs @@ -1,5 +1,5 @@ use lazy_static::lazy_static; -use nalgebra::{Const, DMatrix, DVector, Dyn}; +use nalgebra::{Const, DMatrix, DVector, DVectorView, Dyn, SymmetricEigen}; use web_sys::{console, wasm_bindgen::JsValue}; /* DEBUG */ // --- elements --- @@ -85,6 +85,75 @@ impl PartialMatrix { } } +// --- configuration subspaces --- + +#[derive(Clone)] +pub struct ConfigSubspace { + assembly_dim: usize, + basis: Vec> +} + +impl ConfigSubspace { + 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 + // space for `assembly_dim` elements. we consider an eigenvector to be part + // of the kernel if its eigenvalue is smaller than the constant `THRESHOLD` + fn symmetric_kernel(a: DMatrix, assembly_dim: usize) -> ConfigSubspace { + const ELEMENT_DIM: usize = 5; + const THRESHOLD: f64 = 1.0e-4; + let eig = SymmetricEigen::new(a); + let eig_vecs = eig.eigenvectors.column_iter(); + let eig_pairs = eig.eigenvalues.iter().zip(eig_vecs); + let basis = eig_pairs.filter_map( + |(λ, v)| (λ.abs() < THRESHOLD).then_some( + Into::>::into( + v.reshape_generic(Dyn(ELEMENT_DIM), Dyn(assembly_dim)) + ) + ) + ); + + /* DEBUG */ + // print the eigenvalues + #[cfg(all(target_family = "wasm", target_os = "unknown"))] + console::log_1(&JsValue::from( + format!("Eigenvalues used to find kernel:{}", eig.eigenvalues) + )); + + ConfigSubspace { + assembly_dim: assembly_dim, + basis: basis.collect() + } + } + + pub fn dim(&self) -> usize { + self.basis.len() + } + + pub fn assembly_dim(&self) -> usize { + self.assembly_dim + } + + // find the projection onto this subspace, with respect to the Euclidean + // inner product, of the motion where the element with the given column + // index has velocity `v` + pub fn proj(&self, v: &DVectorView, column_index: usize) -> DMatrix { + 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() + } + } +} + // --- descent history --- pub struct DescentHistory { @@ -181,7 +250,7 @@ pub fn realize_gram( reg_scale: f64, max_descent_steps: i32, max_backoff_steps: i32 -) -> (DMatrix, bool, DescentHistory) { +) -> (DMatrix, ConfigSubspace, bool, DescentHistory) { // start the descent history let mut history = DescentHistory::new(); @@ -201,12 +270,8 @@ pub fn realize_gram( // use Newton's method with backtracking and gradient descent backup let mut state = SearchState::from_config(gram, guess); + let mut hess = DMatrix::zeros(element_dim, assembly_dim); for _ in 0..max_descent_steps { - // stop if the loss is tolerably low - history.config.push(state.config.clone()); - history.scaled_loss.push(state.loss / scale_adjustment); - if state.loss < tol { break; } - // find the negative gradient of the loss function let neg_grad = 4.0 * &*Q * &state.config * &state.err_proj; let mut neg_grad_stacked = neg_grad.clone().reshape_generic(Dyn(total_dim), Const::<1>); @@ -229,7 +294,7 @@ pub fn realize_gram( hess_cols.push(deriv_grad.reshape_generic(Dyn(total_dim), Const::<1>)); } } - let mut hess = DMatrix::from_columns(hess_cols.as_slice()); + hess = DMatrix::from_columns(hess_cols.as_slice()); // regularize the Hessian let min_eigval = hess.symmetric_eigenvalues().min(); @@ -249,6 +314,11 @@ pub fn realize_gram( hess[(k, k)] = 1.0; } + // stop if the loss is tolerably low + history.config.push(state.config.clone()); + history.scaled_loss.push(state.loss / scale_adjustment); + if state.loss < tol { break; } + // compute the Newton step /* we need to either handle or eliminate the case where the minimum @@ -256,7 +326,7 @@ pub fn realize_gram( singular. right now, this causes the Cholesky decomposition to return `None`, leading to a panic when we unrap */ - let base_step_stacked = hess.cholesky().unwrap().solve(&neg_grad_stacked); + let base_step_stacked = hess.clone().cholesky().unwrap().solve(&neg_grad_stacked); let base_step = base_step_stacked.reshape_generic(Dyn(element_dim), Dyn(assembly_dim)); history.base_step.push(base_step.clone()); @@ -269,10 +339,16 @@ pub fn realize_gram( state = better_state; history.backoff_steps.push(backoff_steps); }, - None => return (state.config, false, history) + None => return (state.config, ConfigSubspace::zero(assembly_dim), false, history) }; } - (state.config, state.loss < tol, history) + let success = state.loss < tol; + let tangent = if success { + ConfigSubspace::symmetric_kernel(hess, assembly_dim) + } else { + ConfigSubspace::zero(assembly_dim) + }; + (state.config, tangent, success, history) } // --- tests --- @@ -291,7 +367,7 @@ pub mod irisawa { use super::*; - pub fn realize_irisawa_hexlet(scaled_tol: f64) -> (DMatrix, bool, DescentHistory) { + pub fn realize_irisawa_hexlet(scaled_tol: f64) -> (DMatrix, ConfigSubspace, bool, DescentHistory) { let gram = { let mut gram_to_be = PartialMatrix::new(); for s in 0..9 { @@ -399,7 +475,7 @@ mod tests { fn irisawa_hexlet_test() { // solve Irisawa's problem const SCALED_TOL: f64 = 1.0e-12; - let (config, _, _) = realize_irisawa_hexlet(SCALED_TOL); + let (config, _, _, _) = realize_irisawa_hexlet(SCALED_TOL); // check against Irisawa's solution let entry_tol = SCALED_TOL.sqrt(); @@ -409,6 +485,61 @@ mod tests { } } + #[test] + fn tangent_test() { + const SCALED_TOL: f64 = 1.0e-12; + const ELEMENT_DIM: usize = 5; + const ASSEMBLY_DIM: usize = 3; + let gram = { + let mut gram_to_be = PartialMatrix::new(); + for j in 0..3 { + for k in j..3 { + gram_to_be.push_sym(j, k, if j == k { 1.0 } else { -1.0 }); + } + } + gram_to_be + }; + let guess = DMatrix::from_columns(&[ + sphere(0.0, 0.0, 0.0, -2.0), + sphere(0.0, 0.0, 1.0, 1.0), + sphere(0.0, 0.0, -1.0, 1.0) + ]); + let frozen: [_; 5] = std::array::from_fn(|k| (k, 0)); + let (config, tangent, success, history) = realize_gram( + &gram, guess.clone(), &frozen, + SCALED_TOL, 0.5, 0.9, 1.1, 200, 110 + ); + assert_eq!(config, guess); + assert_eq!(success, true); + assert_eq!(history.scaled_loss.len(), 1); + + // confirm that the tangent space has dimension five or less + let ConfigSubspace(ref tangent_basis) = tangent; + assert_eq!(tangent_basis.len(), 5); + + // confirm that the tangent space contains all the motions we expect it + // to. since we've already bounded the dimension of the tangent space, + // this confirms that the tangent space is what we expect it to be + let tangent_motions = vec![ + basis_matrix((0, 1), ELEMENT_DIM, ASSEMBLY_DIM), + basis_matrix((1, 1), ELEMENT_DIM, ASSEMBLY_DIM), + basis_matrix((0, 2), ELEMENT_DIM, ASSEMBLY_DIM), + basis_matrix((1, 2), ELEMENT_DIM, ASSEMBLY_DIM), + DMatrix::::from_column_slice(ELEMENT_DIM, 3, &[ + 0.0, 0.0, 0.0, 0.0, 0.0, + 0.0, 0.0, -1.0, -0.25, -1.0, + 0.0, 0.0, -1.0, 0.25, 1.0 + ]) + ]; + let tol_sq = ((ELEMENT_DIM * ASSEMBLY_DIM) as f64) * SCALED_TOL * SCALED_TOL; + for motion in tangent_motions { + let motion_proj: DMatrix<_> = motion.column_iter().enumerate().map( + |(k, v)| tangent.proj(&v, k) + ).sum(); + assert!((motion - motion_proj).norm_squared() < tol_sq); + } + } + // at the frozen indices, the optimization steps should have exact zeros, // and the realized configuration should match the initial guess #[test] @@ -428,7 +559,7 @@ mod tests { ]); let frozen = [(3, 0), (3, 1)]; println!(); - let (config, success, history) = realize_gram( + let (config, _, success, history) = realize_gram( &gram, guess.clone(), &frozen, 1.0e-12, 0.5, 0.9, 1.1, 200, 110 ); diff --git a/app-proto/src/main.rs b/app-proto/src/main.rs index 8a012d3..f961504 100644 --- a/app-proto/src/main.rs +++ b/app-proto/src/main.rs @@ -46,6 +46,10 @@ impl AppState { } fn main() { + // set the console error panic hook + #[cfg(feature = "console_error_panic_hook")] + console_error_panic_hook::set_once(); + sycamore::render(|| { provide_context(AppState::new());