StudioInfinity/dyna3
2
0
Fork 2
Code Issues 65 Pull requests 2 Projects 3 Releases 4 Wiki Activity Actions

Compare commits

base: StudioInfinity:19e5458e1b135a762ac292de07b66d1b331487a6
Branches Tags
StudioInfinity:main
StudioInfinity:readme-toplevel
StudioInfinity:kaleidocycle-example
StudioInfinity:outline-update-fix
StudioInfinity:lang-trials
StudioInfinity:seed-problem
StudioInfinity:0.3pre
StudioInfinity:0.2
StudioInfinity:v0.1.1
StudioInfinity:v0.1.0
...
compare: StudioInfinity:2f6c93ec4294752b328e621169def02d26a62c75
Branches Tags
StudioInfinity:main
StudioInfinity:readme-toplevel
StudioInfinity:kaleidocycle-example
StudioInfinity:outline-update-fix
StudioInfinity:lang-trials
StudioInfinity:seed-problem
StudioInfinity:0.3pre
StudioInfinity:0.2
StudioInfinity:v0.1.1
StudioInfinity:v0.1.0

2 commits
19e5458e1b ... 2f6c93ec42

Author SHA1 Message Date
Glen Whitney
2f6c93ec42 chore:wrap lines at 80 characters
All checks were successful
/ test (pull_request) Successful in 3m39s
Details
2025-10-10 22:34:18 -07:00
Glen Whitney
3405797d14 chore: remove trailing whitespace, add CR at end of file
All checks were successful
/ test (pull_request) Successful in 3m36s
Details
2025-10-10 18:08:31 -07:00
12 changed files with 629 additions and 468 deletions
Download patch file Download diff file Expand all files Collapse all files

218
app-proto/src/assembly.rs
View file

@ -45,7 +45,7 @@ static NEXT_SERIAL: AtomicU64 = AtomicU64::new(0);
pub trait Serial { pub trait Serial {
// a serial number that uniquely identifies this element // a serial number that uniquely identifies this element
fn serial(&self) -> u64; fn serial(&self) -> u64;
// take the next serial number, panicking if that was the last one left // take the next serial number, panicking if that was the last one left
fn next_serial() -> u64 where Self: Sized { fn next_serial() -> u64 where Self: Sized {
// the technique we use to panic on overflow is taken from _Rust Atomics // the technique we use to panic on overflow is taken from _Rust Atomics
@ -101,33 +101,33 @@ pub trait ProblemPoser {
pub trait Element: Serial + ProblemPoser + DisplayItem { pub trait Element: Serial + ProblemPoser + DisplayItem {
// the default identifier for an element of this type // the default identifier for an element of this type
fn default_id() -> String where Self: Sized; fn default_id() -> String where Self: Sized;
// the default example of an element of this type // the default example of an element of this type
fn default(id: String, id_num: u64) -> Self where Self: Sized; fn default(id: String, id_num: u64) -> Self where Self: Sized;
// the default regulators that come with this element // the default regulators that come with this element
fn default_regulators(self: Rc<Self>) -> Vec<Rc<dyn Regulator>> { fn default_regulators(self: Rc<Self>) -> Vec<Rc<dyn Regulator>> {
Vec::new() Vec::new()
} }
fn id(&self) -> &String; fn id(&self) -> &String;
fn label(&self) -> &String; fn label(&self) -> &String;
fn representation(&self) -> Signal<DVector<f64>>; fn representation(&self) -> Signal<DVector<f64>>;
fn ghost(&self) -> Signal<bool>; fn ghost(&self) -> Signal<bool>;
// the regulators the element is subject to. the assembly that owns the // the regulators the element is subject to. the assembly that owns the
// element is responsible for keeping this set up to date // element is responsible for keeping this set up to date
fn regulators(&self) -> Signal<BTreeSet<Rc<dyn Regulator>>>; fn regulators(&self) -> Signal<BTreeSet<Rc<dyn Regulator>>>;
// project a representation vector for this kind of element onto its // project a representation vector for this kind of element onto its
// normalization variety // normalization variety
fn project_to_normalized(&self, rep: &mut DVector<f64>); fn project_to_normalized(&self, rep: &mut DVector<f64>);
// the configuration matrix column index that was assigned to the element // the configuration matrix column index that was assigned to the element
// last time the assembly was realized, or `None` if the element has never // last time the assembly was realized, or `None` if the element has never
// been through a realization // been through a realization
fn column_index(&self) -> Option<usize>; fn column_index(&self) -> Option<usize>;
// assign the element a configuration matrix column index. this method must // assign the element a configuration matrix column index. this method must
// be used carefully to preserve invariant (1), described in the comment on // be used carefully to preserve invariant (1), described in the comment on
// the `tangent` field of the `Assembly` structure // the `tangent` field of the `Assembly` structure
@ -179,7 +179,7 @@ pub struct Sphere {
impl Sphere { impl Sphere {
const CURVATURE_COMPONENT: usize = 3; const CURVATURE_COMPONENT: usize = 3;
pub fn new( pub fn new(
id: String, id: String,
label: String, label: String,
@ -203,7 +203,7 @@ impl Element for Sphere {
fn default_id() -> String { fn default_id() -> String {
"sphere".to_string() "sphere".to_string()
} }
fn default(id: String, id_num: u64) -> Self { fn default(id: String, id_num: u64) -> Self {
Self::new( Self::new(
id, id,
@ -212,39 +212,39 @@ impl Element for Sphere {
sphere(0.0, 0.0, 0.0, 1.0), sphere(0.0, 0.0, 0.0, 1.0),
) )
} }
fn default_regulators(self: Rc<Self>) -> Vec<Rc<dyn Regulator>> { fn default_regulators(self: Rc<Self>) -> Vec<Rc<dyn Regulator>> {
vec![Rc::new(HalfCurvatureRegulator::new(self))] vec![Rc::new(HalfCurvatureRegulator::new(self))]
} }
fn id(&self) -> &String { fn id(&self) -> &String {
&self.id &self.id
} }
fn label(&self) -> &String { fn label(&self) -> &String {
&self.label &self.label
} }
fn representation(&self) -> Signal<DVector<f64>> { fn representation(&self) -> Signal<DVector<f64>> {
self.representation self.representation
} }
fn ghost(&self) -> Signal<bool> { fn ghost(&self) -> Signal<bool> {
self.ghost self.ghost
} }
fn regulators(&self) -> Signal<BTreeSet<Rc<dyn Regulator>>> { fn regulators(&self) -> Signal<BTreeSet<Rc<dyn Regulator>>> {
self.regulators self.regulators
} }
fn project_to_normalized(&self, rep: &mut DVector<f64>) { fn project_to_normalized(&self, rep: &mut DVector<f64>) {
project_sphere_to_normalized(rep); project_sphere_to_normalized(rep);
} }
fn column_index(&self) -> Option<usize> { fn column_index(&self) -> Option<usize> {
self.column_index.get() self.column_index.get()
} }
fn set_column_index(&self, index: usize) { fn set_column_index(&self, index: usize) {
self.column_index.set(Some(index)); self.column_index.set(Some(index));
} }
@ -261,7 +261,8 @@ impl ProblemPoser for Sphere {
let index = self.column_index().expect( let index = self.column_index().expect(
indexing_error("Sphere", &self.id, "it").as_str()); indexing_error("Sphere", &self.id, "it").as_str());
problem.gram.push_sym(index, index, 1.0); problem.gram.push_sym(index, index, 1.0);
problem.guess.set_column(index, &self.representation.get_clone_untracked()); problem.guess.set_column(
index, &self.representation.get_clone_untracked());
} }
} }
@ -279,7 +280,7 @@ pub struct Point {
impl Point { impl Point {
const WEIGHT_COMPONENT: usize = 3; const WEIGHT_COMPONENT: usize = 3;
const NORM_COMPONENT: usize = 4; const NORM_COMPONENT: usize = 4;
pub fn new( pub fn new(
id: String, id: String,
label: String, label: String,
@ -303,7 +304,7 @@ impl Element for Point {
fn default_id() -> String { fn default_id() -> String {
"point".to_string() "point".to_string()
} }
fn default(id: String, id_num: u64) -> Self { fn default(id: String, id_num: u64) -> Self {
Self::new( Self::new(
id, id,
@ -321,35 +322,35 @@ impl Element for Point {
}) })
.collect() .collect()
} }
fn id(&self) -> &String { fn id(&self) -> &String {
&self.id &self.id
} }
fn label(&self) -> &String { fn label(&self) -> &String {
&self.label &self.label
} }
fn representation(&self) -> Signal<DVector<f64>> { fn representation(&self) -> Signal<DVector<f64>> {
self.representation self.representation
} }
fn ghost(&self) -> Signal<bool> { fn ghost(&self) -> Signal<bool> {
self.ghost self.ghost
} }
fn regulators(&self) -> Signal<BTreeSet<Rc<dyn Regulator>>> { fn regulators(&self) -> Signal<BTreeSet<Rc<dyn Regulator>>> {
self.regulators self.regulators
} }
fn project_to_normalized(&self, rep: &mut DVector<f64>) { fn project_to_normalized(&self, rep: &mut DVector<f64>) {
project_point_to_normalized(rep); project_point_to_normalized(rep);
} }
fn column_index(&self) -> Option<usize> { fn column_index(&self) -> Option<usize> {
self.column_index.get() self.column_index.get()
} }
fn set_column_index(&self, index: usize) { fn set_column_index(&self, index: usize) {
self.column_index.set(Some(index)); self.column_index.set(Some(index));
} }
@ -367,7 +368,8 @@ impl ProblemPoser for Point {
indexing_error("Point", &self.id, "it").as_str()); indexing_error("Point", &self.id, "it").as_str());
problem.gram.push_sym(index, index, 0.0); problem.gram.push_sym(index, index, 0.0);
problem.frozen.push(Self::WEIGHT_COMPONENT, index, 0.5); problem.frozen.push(Self::WEIGHT_COMPONENT, index, 0.5);
problem.guess.set_column(index, &self.representation.get_clone_untracked()); problem.guess.set_column(
index, &self.representation.get_clone_untracked());
} }
} }
@ -412,7 +414,8 @@ pub struct InversiveDistanceRegulator {
impl InversiveDistanceRegulator { impl InversiveDistanceRegulator {
pub fn new(subjects: [Rc<dyn Element>; 2]) -> Self { pub fn new(subjects: [Rc<dyn Element>; 2]) -> Self {
let representations = subjects.each_ref().map(|subj| subj.representation()); let representations = subjects.each_ref().map(
|subj| subj.representation());
let measurement = create_memo(move || { let measurement = create_memo(move || {
representations[0].with(|rep_0| representations[0].with(|rep_0|
representations[1].with(|rep_1| representations[1].with(|rep_1|
@ -420,10 +423,10 @@ impl InversiveDistanceRegulator {
) )
) )
}); });
let set_point = create_signal(SpecifiedValue::from_empty_spec()); let set_point = create_signal(SpecifiedValue::from_empty_spec());
let serial = Self::next_serial(); let serial = Self::next_serial();
Self { subjects, measurement, set_point, serial } Self { subjects, measurement, set_point, serial }
} }
} }
@ -432,11 +435,11 @@ impl Regulator for InversiveDistanceRegulator {
fn subjects(&self) -> Vec<Rc<dyn Element>> { fn subjects(&self) -> Vec<Rc<dyn Element>> {
self.subjects.clone().into() self.subjects.clone().into()
} }
fn measurement(&self) -> ReadSignal<f64> { fn measurement(&self) -> ReadSignal<f64> {
self.measurement self.measurement
} }
fn set_point(&self) -> Signal<SpecifiedValue> { fn set_point(&self) -> Signal<SpecifiedValue> {
self.set_point self.set_point
} }
@ -475,10 +478,10 @@ impl HalfCurvatureRegulator {
let measurement = subject.representation().map( let measurement = subject.representation().map(
|rep| rep[Sphere::CURVATURE_COMPONENT] |rep| rep[Sphere::CURVATURE_COMPONENT]
); );
let set_point = create_signal(SpecifiedValue::from_empty_spec()); let set_point = create_signal(SpecifiedValue::from_empty_spec());
let serial = Self::next_serial(); let serial = Self::next_serial();
Self { subject, measurement, set_point, serial } Self { subject, measurement, set_point, serial }
} }
} }
@ -487,11 +490,11 @@ impl Regulator for HalfCurvatureRegulator {
fn subjects(&self) -> Vec<Rc<dyn Element>> { fn subjects(&self) -> Vec<Rc<dyn Element>> {
vec![self.subject.clone()] vec![self.subject.clone()]
} }
fn measurement(&self) -> ReadSignal<f64> { fn measurement(&self) -> ReadSignal<f64> {
self.measurement self.measurement
} }
fn set_point(&self) -> Signal<SpecifiedValue> { fn set_point(&self) -> Signal<SpecifiedValue> {
self.set_point self.set_point
} }
@ -545,7 +548,9 @@ impl PointCoordinateRegulator {
move |rep| rep[axis as usize] move |rep| rep[axis as usize]
); );
let set_point = create_signal(SpecifiedValue::from_empty_spec()); let set_point = create_signal(SpecifiedValue::from_empty_spec());
Self { subject, axis, measurement, set_point, serial: Self::next_serial() } Self {
subject, axis, measurement, set_point, serial: Self::next_serial()
}
} }
} }
@ -579,8 +584,8 @@ impl ProblemPoser for PointCoordinateRegulator {
} }
if nset == Axis::CARDINALITY { if nset == Axis::CARDINALITY {
let [x, y, z] = coords; let [x, y, z] = coords;
problem.frozen.push( problem.frozen.push(Point::NORM_COMPONENT,
Point::NORM_COMPONENT, col, point(x,y,z)[Point::NORM_COMPONENT]); col, point(x,y,z)[Point::NORM_COMPONENT]);
} }
} }
}); });
@ -601,7 +606,7 @@ pub struct Assembly {
// elements and regulators // elements and regulators
pub elements: Signal<BTreeSet<Rc<dyn Element>>>, pub elements: Signal<BTreeSet<Rc<dyn Element>>>,
pub regulators: Signal<BTreeSet<Rc<dyn Regulator>>>, pub regulators: Signal<BTreeSet<Rc<dyn Regulator>>>,
// solution variety tangent space. the basis vectors are stored in // solution variety tangent space. the basis vectors are stored in
// configuration matrix format, ordered according to the elements' column // configuration matrix format, ordered according to the elements' column
// indices. when you realize the assembly, every element that's present // indices. when you realize the assembly, every element that's present
@ -613,13 +618,13 @@ pub struct Assembly {
// in that column of the tangent space basis matrices // in that column of the tangent space basis matrices
// //
pub tangent: Signal<ConfigSubspace>, pub tangent: Signal<ConfigSubspace>,
// indexing // indexing
pub elements_by_id: Signal<BTreeMap<String, Rc<dyn Element>>>, pub elements_by_id: Signal<BTreeMap<String, Rc<dyn Element>>>,
// realization control // realization control
pub realization_trigger: Signal<()>, pub realization_trigger: Signal<()>,
// realization diagnostics // realization diagnostics
pub realization_status: Signal<Result<(), String>>, pub realization_status: Signal<Result<(), String>>,
pub descent_history: Signal<DescentHistory>, pub descent_history: Signal<DescentHistory>,
@ -639,7 +644,7 @@ impl Assembly {
descent_history: create_signal(DescentHistory::new()), descent_history: create_signal(DescentHistory::new()),
step: create_signal(SpecifiedValue::from_empty_spec()), step: create_signal(SpecifiedValue::from_empty_spec()),
}; };
// realize the assembly whenever the element list, the regulator list, // realize the assembly whenever the element list, the regulator list,
// a regulator's set point, or the realization trigger is updated // a regulator's set point, or the realization trigger is updated
let assembly_for_realization = assembly.clone(); let assembly_for_realization = assembly.clone();
@ -653,7 +658,7 @@ impl Assembly {
assembly_for_realization.realization_trigger.track(); assembly_for_realization.realization_trigger.track();
assembly_for_realization.realize(); assembly_for_realization.realize();
}); });
// load a configuration from the descent history whenever the active // load a configuration from the descent history whenever the active
// step is updated // step is updated
let assembly_for_step_selection = assembly.clone(); let assembly_for_step_selection = assembly.clone();
@ -665,12 +670,12 @@ impl Assembly {
assembly_for_step_selection.load_config(&config) assembly_for_step_selection.load_config(&config)
} }
}); });
assembly assembly
} }
// --- inserting elements and regulators --- // --- inserting elements and regulators ---
// insert an element 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 // 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 // same identifier will get kicked out of the `elements_by_id` index
@ -679,22 +684,24 @@ impl Assembly {
let id = elt.id().clone(); let id = elt.id().clone();
let elt_rc = Rc::new(elt); let elt_rc = Rc::new(elt);
self.elements.update(|elts| elts.insert(elt_rc.clone())); self.elements.update(|elts| elts.insert(elt_rc.clone()));
self.elements_by_id.update(|elts_by_id| elts_by_id.insert(id, elt_rc.clone())); self.elements_by_id.update(
|elts_by_id| elts_by_id.insert(id, elt_rc.clone()));
// create and insert the element's default regulators // create and insert the element's default regulators
for reg in elt_rc.default_regulators() { for reg in elt_rc.default_regulators() {
self.insert_regulator(reg); self.insert_regulator(reg);
} }
} }
pub fn try_insert_element(&self, elt: impl Element + 'static) -> bool { pub fn try_insert_element(&self, elt: impl Element + 'static) -> bool {
let can_insert = self.elements_by_id.with_untracked(|elts_by_id| !elts_by_id.contains_key(elt.id())); let can_insert = self.elements_by_id.with_untracked(
|elts_by_id| !elts_by_id.contains_key(elt.id()));
if can_insert { if can_insert {
self.insert_element_unchecked(elt); self.insert_element_unchecked(elt);
} }
can_insert can_insert
} }
pub fn insert_element_default<T: Element + 'static>(&self) { pub fn insert_element_default<T: Element + 'static>(&self) {
// find the next unused identifier in the default sequence // find the next unused identifier in the default sequence
let default_id = T::default_id(); let default_id = T::default_id();
@ -706,17 +713,17 @@ impl Assembly {
id_num += 1; id_num += 1;
id = format!("{default_id}{id_num}"); id = format!("{default_id}{id_num}");
} }
// create and insert the default example of `T` // create and insert the default example of `T`
let _ = self.insert_element_unchecked(T::default(id, id_num)); let _ = self.insert_element_unchecked(T::default(id, id_num));
} }
pub fn insert_regulator(&self, regulator: Rc<dyn Regulator>) { pub fn insert_regulator(&self, regulator: Rc<dyn Regulator>) {
// add the regulator to the assembly's regulator list // add the regulator to the assembly's regulator list
self.regulators.update( self.regulators.update(
|regs| regs.insert(regulator.clone()) |regs| regs.insert(regulator.clone())
); );
// add the regulator to each subject's regulator list // add the regulator to each subject's regulator list
let subject_regulators: Vec<_> = regulator.subjects().into_iter().map( let subject_regulators: Vec<_> = regulator.subjects().into_iter().map(
|subj| subj.regulators() |subj| subj.regulators()
@ -724,7 +731,7 @@ impl Assembly {
for regulators in subject_regulators { for regulators in subject_regulators {
regulators.update(|regs| regs.insert(regulator.clone())); regulators.update(|regs| regs.insert(regulator.clone()));
} }
/* DEBUG */ /* DEBUG */
// print an updated list of regulators // print an updated list of regulators
console_log!("Regulators:"); console_log!("Regulators:");
@ -747,19 +754,20 @@ impl Assembly {
} }
}); });
} }
// --- updating the configuration --- // --- updating the configuration ---
pub fn load_config(&self, config: &DMatrix<f64>) { pub fn load_config(&self, config: &DMatrix<f64>) {
for elt in self.elements.get_clone_untracked() { for elt in self.elements.get_clone_untracked() {
elt.representation().update( elt.representation().update(
|rep| rep.set_column(0, &config.column(elt.column_index().unwrap())) |rep| rep.set_column(
0, &config.column(elt.column_index().unwrap()))
); );
} }
} }
// --- realization --- // --- realization ---
pub fn realize(&self) { pub fn realize(&self) {
// index the elements // index the elements
self.elements.update_silent(|elts| { self.elements.update_silent(|elts| {
@ -767,7 +775,7 @@ impl Assembly {
elt.set_column_index(index); elt.set_column_index(index);
} }
}); });
// set up the constraint problem // set up the constraint problem
let problem = self.elements.with_untracked(|elts| { let problem = self.elements.with_untracked(|elts| {
let mut problem = ConstraintProblem::new(elts.len()); let mut problem = ConstraintProblem::new(elts.len());
@ -781,21 +789,21 @@ impl Assembly {
}); });
problem problem
}); });
/* DEBUG */ /* DEBUG */
// log the Gram matrix // log the Gram matrix
console_log!("Gram matrix:\n{}", problem.gram); console_log!("Gram matrix:\n{}", problem.gram);
console_log!("Frozen entries:\n{}", problem.frozen); console_log!("Frozen entries:\n{}", problem.frozen);
/* DEBUG */ /* DEBUG */
// log the initial configuration matrix // log the initial configuration matrix
console_log!("Old configuration:{:>8.3}", problem.guess); console_log!("Old configuration:{:>8.3}", problem.guess);
// look for a configuration with the given Gram matrix // look for a configuration with the given Gram matrix
let Realization { result, history } = realize_gram( let Realization { result, history } = realize_gram(
&problem, 1.0e-12, 0.5, 0.9, 1.1, 200, 110 &problem, 1.0e-12, 0.5, 0.9, 1.1, 200, 110
); );
/* DEBUG */ /* DEBUG */
// report the outcome of the search in the browser console // report the outcome of the search in the browser console
if let Err(ref message) = result { if let Err(ref message) = result {
@ -807,20 +815,20 @@ impl Assembly {
console_log!("Steps: {}", history.scaled_loss.len() - 1); console_log!("Steps: {}", history.scaled_loss.len() - 1);
console_log!("Loss: {}", history.scaled_loss.last().unwrap()); console_log!("Loss: {}", history.scaled_loss.last().unwrap());
} }
// report the descent history // report the descent history
let step_cnt = history.config.len(); let step_cnt = history.config.len();
self.descent_history.set(history); self.descent_history.set(history);
match result { match result {
Ok(ConfigNeighborhood { nbhd: tangent, .. }) => { Ok(ConfigNeighborhood { nbhd: tangent, .. }) => {
/* DEBUG */ /* DEBUG */
// report the tangent dimension // report the tangent dimension
console_log!("Tangent dimension: {}", tangent.dim()); console_log!("Tangent dimension: {}", tangent.dim());
// report the realization status // report the realization status
self.realization_status.set(Ok(())); self.realization_status.set(Ok(()));
// display the last realization step // display the last realization step
self.step.set( self.step.set(
if step_cnt > 0 { if step_cnt > 0 {
@ -830,7 +838,7 @@ impl Assembly {
SpecifiedValue::from_empty_spec() SpecifiedValue::from_empty_spec()
} }
); );
// save the tangent space // save the tangent space
self.tangent.set_silent(tangent); self.tangent.set_silent(tangent);
}, },
@ -840,15 +848,15 @@ impl Assembly {
// `Err(message)` we received from the match: we're changing the // `Err(message)` we received from the match: we're changing the
// `Ok` type from `Realization` to `()` // `Ok` type from `Realization` to `()`
self.realization_status.set(Err(message)); self.realization_status.set(Err(message));
// display the initial guess // display the initial guess
self.step.set(SpecifiedValue::from(Some(0.0))); self.step.set(SpecifiedValue::from(Some(0.0)));
}, },
} }
} }
// --- deformation --- // --- deformation ---
// project the given motion to the tangent space of the solution variety and // 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) // move the assembly along it. the implementation is based on invariant (1)
// from above and the following additional invariant: // from above and the following additional invariant:
@ -865,7 +873,7 @@ impl Assembly {
if self.tangent.with(|tan| tan.dim() <= 0 && tan.assembly_dim() > 0) { if self.tangent.with(|tan| tan.dim() <= 0 && tan.assembly_dim() > 0) {
console::log_1(&JsValue::from("The assembly is rigid")); console::log_1(&JsValue::from("The assembly is rigid"));
} }
// give a column index to each moving element that doesn't have one yet. // give a column index to each moving element that doesn't have one yet.
// this temporarily breaks invariant (1), but the invariant will be // this temporarily breaks invariant (1), but the invariant will be
// restored when we realize the assembly at the end of the deformation. // restored when we realize the assembly at the end of the deformation.
@ -883,7 +891,7 @@ impl Assembly {
} }
next_column_index next_column_index
}; };
// project the element motions onto the tangent space of the solution // project the element motions onto the tangent space of the solution
// variety and sum them to get a deformation of the whole assembly. the // variety and sum them to get a deformation of the whole assembly. the
// matrix `motion_proj` that holds the deformation has extra columns for // matrix `motion_proj` that holds the deformation has extra columns for
@ -894,11 +902,12 @@ impl Assembly {
// we can unwrap the column index because we know that every moving // we can unwrap the column index because we know that every moving
// element has one at this point // element has one at this point
let column_index = elt_motion.element.column_index().unwrap(); let column_index = elt_motion.element.column_index().unwrap();
if column_index < realized_dim { if column_index < realized_dim {
// this element had a column index when we started, so by // this element had a column index when we started, so by
// invariant (1), it's reflected in the tangent space // invariant (1), it's reflected in the tangent space
let mut target_columns = motion_proj.columns_mut(0, realized_dim); let mut target_columns =
motion_proj.columns_mut(0, realized_dim);
target_columns += self.tangent.with( target_columns += self.tangent.with(
|tan| tan.proj(&elt_motion.velocity, column_index) |tan| tan.proj(&elt_motion.velocity, column_index)
); );
@ -906,13 +915,14 @@ impl Assembly {
// this element didn't have a column index when we started, so // this element didn't have a column index when we started, so
// by invariant (2), it's unconstrained // by invariant (2), it's unconstrained
let mut target_column = motion_proj.column_mut(column_index); let mut target_column = motion_proj.column_mut(column_index);
let unif_to_std = elt_motion.element.representation().with_untracked( let unif_to_std =
|rep| local_unif_to_std(rep.as_view()) elt_motion.element.representation().with_untracked(
); |rep| local_unif_to_std(rep.as_view())
);
target_column += unif_to_std * elt_motion.velocity; target_column += unif_to_std * elt_motion.velocity;
} }
} }
// step the assembly along the deformation. this changes the elements' // step the assembly along the deformation. this changes the elements'
// normalizations, so we restore those afterward // normalizations, so we restore those afterward
for elt in self.elements.get_clone_untracked() { for elt in self.elements.get_clone_untracked() {
@ -925,12 +935,15 @@ impl Assembly {
elt.project_to_normalized(rep); elt.project_to_normalized(rep);
}, },
None => { None => {
console_log!("No velocity to unpack for fresh element \"{}\"", elt.id()) console_log!(
"No velocity to unpack for fresh element \"{}\"",
elt.id()
)
}, },
}; };
}); });
} }
// trigger a realization to bring the configuration back onto the // trigger a realization to bring the configuration back onto the
// solution variety. this also gets the elements' column indices and the // solution variety. this also gets the elements' column indices and the
// saved tangent space back in sync // saved tangent space back in sync
@ -941,9 +954,9 @@ impl Assembly {
#[cfg(test)] #[cfg(test)]
mod tests { mod tests {
use super::*; use super::*;
use crate::engine; use crate::engine;
#[test] #[test]
#[should_panic(expected = #[should_panic(expected =
"Sphere \"sphere\" must be indexed before it writes problem data")] "Sphere \"sphere\" must be indexed before it writes problem data")]
@ -953,25 +966,27 @@ mod tests {
elt.pose(&mut ConstraintProblem::new(1)); elt.pose(&mut ConstraintProblem::new(1));
}); });
} }
#[test] #[test]
#[should_panic(expected = "Subject \"sphere1\" must be indexed before \ #[should_panic(expected = "Subject \"sphere1\" must be indexed before \
inversive distance regulator writes problem data")] inversive distance regulator writes problem data")]
fn unindexed_subject_test_inversive_distance() { fn unindexed_subject_test_inversive_distance() {
let _ = create_root(|| { let _ = create_root(|| {
let subjects = [0, 1].map( let subjects = [0, 1].map(
|k| Rc::new(Sphere::default(format!("sphere{k}"), k)) as Rc<dyn Element> |k| Rc::new(
Sphere::default(format!("sphere{k}"), k)) as Rc<dyn Element>
); );
subjects[0].set_column_index(0); subjects[0].set_column_index(0);
InversiveDistanceRegulator { InversiveDistanceRegulator {
subjects: subjects, subjects: subjects,
measurement: create_memo(|| 0.0), measurement: create_memo(|| 0.0),
set_point: create_signal(SpecifiedValue::try_from("0.0".to_string()).unwrap()), set_point: create_signal(
SpecifiedValue::try_from("0.0".to_string()).unwrap()),
serial: InversiveDistanceRegulator::next_serial() serial: InversiveDistanceRegulator::next_serial()
}.pose(&mut ConstraintProblem::new(2)); }.pose(&mut ConstraintProblem::new(2));
}); });
} }
#[test] #[test]
fn curvature_drift_test() { fn curvature_drift_test() {
const INITIAL_RADIUS: f64 = 0.25; const INITIAL_RADIUS: f64 = 0.25;
@ -991,12 +1006,14 @@ inversive distance regulator writes problem data")]
engine::sphere(0.0, 0.0, 0.0, INITIAL_RADIUS), engine::sphere(0.0, 0.0, 0.0, INITIAL_RADIUS),
) )
); );
// nudge the sphere repeatedly along the `z` axis // nudge the sphere repeatedly along the `z` axis
const STEP_SIZE: f64 = 0.0025; const STEP_SIZE: f64 = 0.0025;
const STEP_CNT: usize = 400; const STEP_CNT: usize = 400;
let sphere = assembly.elements_by_id.with(|elts_by_id| elts_by_id[sphere_id].clone()); let sphere = assembly.elements_by_id.with(
let velocity = DVector::from_column_slice(&[0.0, 0.0, STEP_SIZE, 0.0]); |elts_by_id| elts_by_id[sphere_id].clone());
let velocity =
DVector::from_column_slice(&[0.0, 0.0, STEP_SIZE, 0.0]);
for _ in 0..STEP_CNT { for _ in 0..STEP_CNT {
assembly.deform( assembly.deform(
vec![ vec![
@ -1007,14 +1024,15 @@ inversive distance regulator writes problem data")]
] ]
); );
} }
// check how much the sphere's curvature has drifted // check how much the sphere's curvature has drifted
const INITIAL_HALF_CURV: f64 = 0.5 / INITIAL_RADIUS; const INITIAL_HALF_CURV: f64 = 0.5 / INITIAL_RADIUS;
const DRIFT_TOL: f64 = 0.015; const DRIFT_TOL: f64 = 0.015;
let final_half_curv = sphere.representation().with_untracked( let final_half_curv = sphere.representation().with_untracked(
|rep| rep[Sphere::CURVATURE_COMPONENT] |rep| rep[Sphere::CURVATURE_COMPONENT]
); );
assert!((final_half_curv / INITIAL_HALF_CURV - 1.0).abs() < DRIFT_TOL); assert!((final_half_curv / INITIAL_HALF_CURV - 1.0).abs()
< DRIFT_TOL);
}); });
} }
} }

4
app-proto/src/components/add_remove.rs
View file

@ -39,7 +39,9 @@ pub fn AddRemove() -> View {
} }
) { "Add point" } ) { "Add point" }
button( button(
class = "emoji", /* KLUDGE */ // for convenience, we're using an emoji as a temporary icon for this button /* KLUDGE */ // for convenience, we're using an emoji as an
// icon for this button
class = "emoji",
disabled = { disabled = {
let state = use_context::<AppState>(); let state = use_context::<AppState>();
state.selection.with(|sel| sel.len() != 2) state.selection.with(|sel| sel.len() != 2)

30
app-proto/src/components/diagnostics.rs
View file

@ -54,7 +54,7 @@ fn StepInput() -> View {
// get the assembly // get the assembly
let state = use_context::<AppState>(); let state = use_context::<AppState>();
let assembly = state.assembly; let assembly = state.assembly;
// the `last_step` signal holds the index of the last step // the `last_step` signal holds the index of the last step
let last_step = assembly.descent_history.map( let last_step = assembly.descent_history.map(
|history| match history.config.len() { |history| match history.config.len() {
@ -63,15 +63,15 @@ fn StepInput() -> View {
} }
); );
let input_max = last_step.map(|last| last.unwrap_or(0)); let input_max = last_step.map(|last| last.unwrap_or(0));
// these signals hold the entered step number // these signals hold the entered step number
let value = create_signal(String::new()); let value = create_signal(String::new());
let value_as_number = create_signal(0.0); let value_as_number = create_signal(0.0);
create_effect(move || { create_effect(move || {
value.set(assembly.step.with(|n| n.spec.clone())); value.set(assembly.step.with(|n| n.spec.clone()));
}); });
view! { view! {
div(id = "step-input") { div(id = "step-input") {
label { "Step" } label { "Step" }
@ -98,7 +98,7 @@ fn StepInput() -> View {
|val| val.clamp(0.0, input_max.get() as f64) |val| val.clamp(0.0, input_max.get() as f64)
) )
); );
// set the input string and the assembly's active step // set the input string and the assembly's active step
value.set(step.spec.clone()); value.set(step.spec.clone());
assembly.step.set(step); assembly.step.set(step);
@ -124,7 +124,7 @@ fn LossHistory() -> View {
const CONTAINER_ID: &str = "loss-history"; const CONTAINER_ID: &str = "loss-history";
let state = use_context::<AppState>(); let state = use_context::<AppState>();
let renderer = WasmRenderer::new_opt(None, Some(178)); let renderer = WasmRenderer::new_opt(None, Some(178));
on_mount(move || { on_mount(move || {
create_effect(move || { create_effect(move || {
// get the loss history // get the loss history
@ -136,13 +136,13 @@ fn LossHistory() -> View {
.map(into_log10_time_point) .map(into_log10_time_point)
.collect() .collect()
); );
// initialize the chart axes // initialize the chart axes
let step_axis = Axis::new() let step_axis = Axis::new()
.type_(AxisType::Category) .type_(AxisType::Category)
.boundary_gap(false); .boundary_gap(false);
let scaled_loss_axis = Axis::new(); let scaled_loss_axis = Axis::new();
// load the chart data. when there's no history, we load the data // load the chart data. when there's no history, we load the data
// point (0, None) to clear the chart. it would feel more natural to // point (0, None) to clear the chart. it would feel more natural to
// load empty data vectors, but that turns out not to clear the // load empty data vectors, but that turns out not to clear the
@ -164,7 +164,7 @@ fn LossHistory() -> View {
renderer.render(CONTAINER_ID, &chart).unwrap(); renderer.render(CONTAINER_ID, &chart).unwrap();
}); });
}); });
view! { view! {
div(id = CONTAINER_ID, class = "diagnostics-chart") div(id = CONTAINER_ID, class = "diagnostics-chart")
} }
@ -176,7 +176,7 @@ fn SpectrumHistory() -> View {
const CONTAINER_ID: &str = "spectrum-history"; const CONTAINER_ID: &str = "spectrum-history";
let state = use_context::<AppState>(); let state = use_context::<AppState>();
let renderer = WasmRenderer::new(478, 178); let renderer = WasmRenderer::new(478, 178);
on_mount(move || { on_mount(move || {
create_effect(move || { create_effect(move || {
// get the spectrum of the Hessian at each step, split into its // get the spectrum of the Hessian at each step, split into its
@ -208,13 +208,13 @@ fn SpectrumHistory() -> View {
): (Vec<_>, Vec<_>) = hess_eigvals_nonzero ): (Vec<_>, Vec<_>) = hess_eigvals_nonzero
.into_iter() .into_iter()
.partition(|&(_, val)| val > 0.0); .partition(|&(_, val)| val > 0.0);
// initialize the chart axes // initialize the chart axes
let step_axis = Axis::new() let step_axis = Axis::new()
.type_(AxisType::Category) .type_(AxisType::Category)
.boundary_gap(false); .boundary_gap(false);
let eigval_axis = Axis::new(); let eigval_axis = Axis::new();
// load the chart data. when there's no history, we load the data // load the chart data. when there's no history, we load the data
// point (0, None) to clear the chart. it would feel more natural to // point (0, None) to clear the chart. it would feel more natural to
// load empty data vectors, but that turns out not to clear the // load empty data vectors, but that turns out not to clear the
@ -270,7 +270,7 @@ fn SpectrumHistory() -> View {
renderer.render(CONTAINER_ID, &chart).unwrap(); renderer.render(CONTAINER_ID, &chart).unwrap();
}); });
}); });
view! { view! {
div(id = CONTAINER_ID, class = "diagnostics-chart") div(id = CONTAINER_ID, class = "diagnostics-chart")
} }
@ -302,7 +302,7 @@ pub fn Diagnostics() -> View {
let diagnostics_state = DiagnosticsState::new("loss".to_string()); let diagnostics_state = DiagnosticsState::new("loss".to_string());
let active_tab = diagnostics_state.active_tab.clone(); let active_tab = diagnostics_state.active_tab.clone();
provide_context(diagnostics_state); provide_context(diagnostics_state);
view! { view! {
div(id = "diagnostics") { div(id = "diagnostics") {
div(id = "diagnostics-bar") { div(id = "diagnostics-bar") {
@ -317,4 +317,4 @@ pub fn Diagnostics() -> View {
DiagnosticsPanel(name = "spectrum") { SpectrumHistory {} } DiagnosticsPanel(name = "spectrum") { SpectrumHistory {} }
} }
} }
} }

278
app-proto/src/components/display.rs
View file

@ -48,11 +48,12 @@ impl SceneSpheres {
highlights: Vec::new(), highlights: Vec::new(),
} }
} }
fn len_i32(&self) -> i32 { fn len_i32(&self) -> i32 {
self.representations.len().try_into().expect("Number of spheres must fit in a 32-bit integer") self.representations.len().try_into().expect(
"Number of spheres must fit in a 32-bit integer")
} }
fn push( fn push(
&mut self, representation: DVector<f64>, &mut self, representation: DVector<f64>,
color: ElementColor, opacity: f32, highlight: f32, color: ElementColor, opacity: f32, highlight: f32,
@ -79,7 +80,7 @@ impl ScenePoints {
selections: Vec::new(), selections: Vec::new(),
} }
} }
fn push( fn push(
&mut self, representation: DVector<f64>, &mut self, representation: DVector<f64>,
color: ElementColor, opacity: f32, highlight: f32, selected: bool, color: ElementColor, opacity: f32, highlight: f32, selected: bool,
@ -107,7 +108,7 @@ impl Scene {
pub trait DisplayItem { pub trait DisplayItem {
fn show(&self, scene: &mut Scene, selected: bool); fn show(&self, scene: &mut Scene, selected: bool);
// the smallest positive depth, represented as a multiple of `dir`, where // the smallest positive depth, represented as a multiple of `dir`, where
// the line generated by `dir` hits the element. returns `None` if the line // the line generated by `dir` hits the element. returns `None` if the line
// misses the element // misses the element
@ -125,14 +126,18 @@ impl DisplayItem for Sphere {
const DEFAULT_OPACITY: f32 = 0.5; const DEFAULT_OPACITY: f32 = 0.5;
const GHOST_OPACITY: f32 = 0.2; const GHOST_OPACITY: f32 = 0.2;
const HIGHLIGHT: f32 = 0.2; const HIGHLIGHT: f32 = 0.2;
let representation = self.representation.get_clone_untracked(); let representation = self.representation.get_clone_untracked();
let color = if selected { self.color.map(|channel| 0.2 + 0.8*channel) } else { self.color }; let color =
let opacity = if self.ghost.get() { GHOST_OPACITY } else { DEFAULT_OPACITY }; if selected { self.color.map(|channel| 0.2 + 0.8*channel) }
else { self.color };
let opacity =
if self.ghost.get() { GHOST_OPACITY }
else { DEFAULT_OPACITY };
let highlight = if selected { 1.0 } else { HIGHLIGHT }; let highlight = if selected { 1.0 } else { HIGHLIGHT };
scene.spheres.push(representation, color, opacity, highlight); scene.spheres.push(representation, color, opacity, highlight);
} }
// this method should be kept synchronized with `sphere_cast` in // this method should be kept synchronized with `sphere_cast` in
// `spheres.frag`, which does essentially the same thing on the GPU side // `spheres.frag`, which does essentially the same thing on the GPU side
fn cast( fn cast(
@ -144,12 +149,13 @@ impl DisplayItem for Sphere {
// if `a/b` is less than this threshold, we approximate // if `a/b` is less than this threshold, we approximate
// `a*u^2 + b*u + c` by the linear function `b*u + c` // `a*u^2 + b*u + c` by the linear function `b*u + c`
const DEG_THRESHOLD: f64 = 1e-9; const DEG_THRESHOLD: f64 = 1e-9;
let rep = self.representation.with_untracked(|rep| assembly_to_world * rep); let rep = self.representation.with_untracked(
|rep| assembly_to_world * rep);
let a = -rep[3] * dir.norm_squared(); let a = -rep[3] * dir.norm_squared();
let b = rep.rows_range(..3).dot(&dir); let b = rep.rows_range(..3).dot(&dir);
let c = -rep[4]; let c = -rep[4];
let adjust = 4.0*a*c/(b*b); let adjust = 4.0*a*c/(b*b);
if adjust < 1.0 { if adjust < 1.0 {
// as long as `b` is non-zero, the linear approximation of // as long as `b` is non-zero, the linear approximation of
@ -184,14 +190,16 @@ impl DisplayItem for Point {
/* SCAFFOLDING */ /* SCAFFOLDING */
const GHOST_OPACITY: f32 = 0.4; const GHOST_OPACITY: f32 = 0.4;
const HIGHLIGHT: f32 = 0.5; const HIGHLIGHT: f32 = 0.5;
let representation = self.representation.get_clone_untracked(); let representation = self.representation.get_clone_untracked();
let color = if selected { self.color.map(|channel| 0.2 + 0.8*channel) } else { self.color }; let color =
if selected { self.color.map(|channel| 0.2 + 0.8*channel) }
else { self.color };
let opacity = if self.ghost.get() { GHOST_OPACITY } else { 1.0 }; let opacity = if self.ghost.get() { GHOST_OPACITY } else { 1.0 };
let highlight = if selected { 1.0 } else { HIGHLIGHT }; let highlight = if selected { 1.0 } else { HIGHLIGHT };
scene.points.push(representation, color, opacity, highlight, selected); scene.points.push(representation, color, opacity, highlight, selected);
} }
/* SCAFFOLDING */ /* SCAFFOLDING */
fn cast( fn cast(
&self, &self,
@ -199,20 +207,21 @@ impl DisplayItem for Point {
assembly_to_world: &DMatrix<f64>, assembly_to_world: &DMatrix<f64>,
pixel_size: f64, pixel_size: f64,
) -> Option<f64> { ) -> Option<f64> {
let rep = self.representation.with_untracked(|rep| assembly_to_world * rep); let rep = self.representation.with_untracked(
|rep| assembly_to_world * rep);
if rep[2] < 0.0 { if rep[2] < 0.0 {
// this constant should be kept synchronized with `point.frag` // this constant should be kept synchronized with `point.frag`
const POINT_RADIUS_PX: f64 = 4.0; const POINT_RADIUS_PX: f64 = 4.0;
// find the radius of the point in screen projection units // find the radius of the point in screen projection units
let point_radius_proj = POINT_RADIUS_PX * pixel_size; let point_radius_proj = POINT_RADIUS_PX * pixel_size;
// find the squared distance between the screen projections of the // find the squared distance between the screen projections of the
// ray and the point // ray and the point
let dir_proj = -dir.fixed_rows::<2>(0) / dir[2]; let dir_proj = -dir.fixed_rows::<2>(0) / dir[2];
let rep_proj = -rep.fixed_rows::<2>(0) / rep[2]; let rep_proj = -rep.fixed_rows::<2>(0) / rep[2];
let dist_sq = (dir_proj - rep_proj).norm_squared(); let dist_sq = (dir_proj - rep_proj).norm_squared();
// if the ray hits the point, return its depth // if the ray hits the point, return its depth
if dist_sq < point_radius_proj * point_radius_proj { if dist_sq < point_radius_proj * point_radius_proj {
Some(rep[2] / dir[2]) Some(rep[2] / dir[2])
@ -254,13 +263,13 @@ fn set_up_program(
WebGl2RenderingContext::FRAGMENT_SHADER, WebGl2RenderingContext::FRAGMENT_SHADER,
fragment_shader_source, fragment_shader_source,
); );
// create the program and attach the shaders // create the program and attach the shaders
let program = context.create_program().unwrap(); let program = context.create_program().unwrap();
context.attach_shader(&program, &vertex_shader); context.attach_shader(&program, &vertex_shader);
context.attach_shader(&program, &fragment_shader); context.attach_shader(&program, &fragment_shader);
context.link_program(&program); context.link_program(&program);
/* DEBUG */ /* DEBUG */
// report whether linking succeeded // report whether linking succeeded
let link_status = context let link_status = context
@ -273,7 +282,7 @@ fn set_up_program(
"Linking failed" "Linking failed"
}; };
console::log_1(&JsValue::from(link_msg)); console::log_1(&JsValue::from(link_msg));
program program
} }
@ -318,7 +327,7 @@ fn load_new_buffer(
// create a buffer and bind it to ARRAY_BUFFER // create a buffer and bind it to ARRAY_BUFFER
let buffer = context.create_buffer(); let buffer = context.create_buffer();
context.bind_buffer(WebGl2RenderingContext::ARRAY_BUFFER, buffer.as_ref()); context.bind_buffer(WebGl2RenderingContext::ARRAY_BUFFER, buffer.as_ref());
// load the given data into the buffer. this block is unsafe because // load the given data into the buffer. this block is unsafe because
// `Float32Array::view` creates a raw view into our module's // `Float32Array::view` creates a raw view into our module's
// `WebAssembly.Memory` buffer. allocating more memory will change the // `WebAssembly.Memory` buffer. allocating more memory will change the
@ -332,7 +341,7 @@ fn load_new_buffer(
WebGl2RenderingContext::STATIC_DRAW, WebGl2RenderingContext::STATIC_DRAW,
); );
} }
buffer buffer
} }
@ -353,15 +362,16 @@ fn event_dir(event: &MouseEvent) -> (Vector3<f64>, f64) {
let width = rect.width(); let width = rect.width();
let height = rect.height(); let height = rect.height();
let shortdim = width.min(height); let shortdim = width.min(height);
// this constant should be kept synchronized with `spheres.frag` and // this constant should be kept synchronized with `spheres.frag` and
// `point.vert` // `point.vert`
const FOCAL_SLOPE: f64 = 0.3; const FOCAL_SLOPE: f64 = 0.3;
let horizontal = f64::from(event.client_x()) - rect.left();
let vertical = rect.bottom() - f64::from(event.client_y());
( (
Vector3::new( Vector3::new(
FOCAL_SLOPE * (2.0*(f64::from(event.client_x()) - rect.left()) - width) / shortdim, FOCAL_SLOPE * (2.0*horizontal - width) / shortdim,
FOCAL_SLOPE * (2.0*(rect.bottom() - f64::from(event.client_y())) - height) / shortdim, FOCAL_SLOPE * (2.0*vertical - height) / shortdim,
-1.0, -1.0,
), ),
FOCAL_SLOPE * 2.0 / shortdim, FOCAL_SLOPE * 2.0 / shortdim,
@ -373,13 +383,13 @@ fn event_dir(event: &MouseEvent) -> (Vector3<f64>, f64) {
#[component] #[component]
pub fn Display() -> View { pub fn Display() -> View {
let state = use_context::<AppState>(); let state = use_context::<AppState>();
// canvas // canvas
let display = create_node_ref(); let display = create_node_ref();
// viewpoint // viewpoint
let assembly_to_world = create_signal(DMatrix::<f64>::identity(5, 5)); let assembly_to_world = create_signal(DMatrix::<f64>::identity(5, 5));
// navigation // navigation
let pitch_up = create_signal(0.0); let pitch_up = create_signal(0.0);
let pitch_down = create_signal(0.0); let pitch_down = create_signal(0.0);
@ -390,7 +400,7 @@ pub fn Display() -> View {
let zoom_in = create_signal(0.0); let zoom_in = create_signal(0.0);
let zoom_out = create_signal(0.0); let zoom_out = create_signal(0.0);
let turntable = create_signal(false); /* BENCHMARKING */ let turntable = create_signal(false); /* BENCHMARKING */
// manipulation // manipulation
let translate_neg_x = create_signal(0.0); let translate_neg_x = create_signal(0.0);
let translate_pos_x = create_signal(0.0); let translate_pos_x = create_signal(0.0);
@ -400,7 +410,7 @@ pub fn Display() -> View {
let translate_pos_z = create_signal(0.0); let translate_pos_z = create_signal(0.0);
let shrink_neg = create_signal(0.0); let shrink_neg = create_signal(0.0);
let shrink_pos = create_signal(0.0); let shrink_pos = create_signal(0.0);
// change listener // change listener
let scene_changed = create_signal(true); let scene_changed = create_signal(true);
create_effect(move || { create_effect(move || {
@ -413,18 +423,18 @@ pub fn Display() -> View {
state.selection.track(); state.selection.track();
scene_changed.set(true); scene_changed.set(true);
}); });
/* INSTRUMENTS */ /* INSTRUMENTS */
const SAMPLE_PERIOD: i32 = 60; const SAMPLE_PERIOD: i32 = 60;
let mut last_sample_time = 0.0; let mut last_sample_time = 0.0;
let mut frames_since_last_sample = 0; let mut frames_since_last_sample = 0;
let mean_frame_interval = create_signal(0.0); let mean_frame_interval = create_signal(0.0);
let assembly_for_raf = state.assembly.clone(); let assembly_for_raf = state.assembly.clone();
on_mount(move || { on_mount(move || {
// timing // timing
let mut last_time = 0.0; let mut last_time = 0.0;
// viewpoint // viewpoint
const ROT_SPEED: f64 = 0.4; // in radians per second const ROT_SPEED: f64 = 0.4; // in radians per second
const ZOOM_SPEED: f64 = 0.15; // multiplicative rate per second const ZOOM_SPEED: f64 = 0.15; // multiplicative rate per second
@ -432,48 +442,50 @@ pub fn Display() -> View {
let mut orientation = DMatrix::<f64>::identity(5, 5); let mut orientation = DMatrix::<f64>::identity(5, 5);
let mut rotation = DMatrix::<f64>::identity(5, 5); let mut rotation = DMatrix::<f64>::identity(5, 5);
let mut location_z: f64 = 5.0; let mut location_z: f64 = 5.0;
// manipulation // manipulation
const TRANSLATION_SPEED: f64 = 0.15; // in length units per second const TRANSLATION_SPEED: f64 = 0.15; // in length units per second
const SHRINKING_SPEED: f64 = 0.15; // in length units per second const SHRINKING_SPEED: f64 = 0.15; // in length units per second
// display parameters // display parameters
const LAYER_THRESHOLD: i32 = 0; /* DEBUG */ const LAYER_THRESHOLD: i32 = 0; /* DEBUG */
const DEBUG_MODE: i32 = 0; /* DEBUG */ const DEBUG_MODE: i32 = 0; /* DEBUG */
/* INSTRUMENTS */ /* INSTRUMENTS */
let performance = window().unwrap().performance().unwrap(); let performance = window().unwrap().performance().unwrap();
// get the display canvas // get the display canvas
let canvas = display.get().unchecked_into::<web_sys::HtmlCanvasElement>(); let canvas =
display.get().unchecked_into::<web_sys::HtmlCanvasElement>();
let ctx = canvas let ctx = canvas
.get_context("webgl2") .get_context("webgl2")
.unwrap() .unwrap()
.unwrap() .unwrap()
.dyn_into::<WebGl2RenderingContext>() .dyn_into::<WebGl2RenderingContext>()
.unwrap(); .unwrap();
// disable depth testing // disable depth testing
ctx.disable(WebGl2RenderingContext::DEPTH_TEST); ctx.disable(WebGl2RenderingContext::DEPTH_TEST);
// set blend mode // set blend mode
ctx.enable(WebGl2RenderingContext::BLEND); ctx.enable(WebGl2RenderingContext::BLEND);
ctx.blend_func(WebGl2RenderingContext::SRC_ALPHA, WebGl2RenderingContext::ONE_MINUS_SRC_ALPHA); ctx.blend_func(WebGl2RenderingContext::SRC_ALPHA,
WebGl2RenderingContext::ONE_MINUS_SRC_ALPHA);
// set up the sphere rendering program // set up the sphere rendering program
let sphere_program = set_up_program( let sphere_program = set_up_program(
&ctx, &ctx,
include_str!("identity.vert"), include_str!("identity.vert"),
include_str!("spheres.frag"), include_str!("spheres.frag"),
); );
// set up the point rendering program // set up the point rendering program
let point_program = set_up_program( let point_program = set_up_program(
&ctx, &ctx,
include_str!("point.vert"), include_str!("point.vert"),
include_str!("point.frag"), include_str!("point.frag"),
); );
/* DEBUG */ /* DEBUG */
// print the maximum number of vectors that can be passed as // print the maximum number of vectors that can be passed as
// uniforms to a fragment shader. the OpenGL ES 3.0 standard // uniforms to a fragment shader. the OpenGL ES 3.0 standard
@ -487,16 +499,20 @@ pub fn Display() -> View {
// machine, the the length of a float or genType array seems to be // machine, the the length of a float or genType array seems to be
// capped at 1024 elements // capped at 1024 elements
console::log_2( console::log_2(
&ctx.get_parameter(WebGl2RenderingContext::MAX_FRAGMENT_UNIFORM_VECTORS).unwrap(), &ctx.get_parameter(
WebGl2RenderingContext::MAX_FRAGMENT_UNIFORM_VECTORS).unwrap(),
&JsValue::from("uniform vectors available"), &JsValue::from("uniform vectors available"),
); );
// find the sphere program's vertex attribute // find the sphere program's vertex attribute
let viewport_position_attr = ctx.get_attrib_location(&sphere_program, "position") as u32; let viewport_position_attr =
ctx.get_attrib_location(&sphere_program, "position") as u32;
// find the sphere program's uniforms // find the sphere program's uniforms
const SPHERE_MAX: usize = 200; const SPHERE_MAX: usize = 200;
let sphere_cnt_loc = ctx.get_uniform_location(&sphere_program, "sphere_cnt"); let sphere_cnt_loc = ctx.get_uniform_location(
&sphere_program, "sphere_cnt"
);
let sphere_sp_locs = get_uniform_array_locations::<SPHERE_MAX>( let sphere_sp_locs = get_uniform_array_locations::<SPHERE_MAX>(
&ctx, &sphere_program, "sphere_list", Some("sp") &ctx, &sphere_program, "sphere_list", Some("sp")
); );
@ -509,11 +525,19 @@ pub fn Display() -> View {
let sphere_highlight_locs = get_uniform_array_locations::<SPHERE_MAX>( let sphere_highlight_locs = get_uniform_array_locations::<SPHERE_MAX>(
&ctx, &sphere_program, "highlight_list", None &ctx, &sphere_program, "highlight_list", None
); );
let resolution_loc = ctx.get_uniform_location(&sphere_program, "resolution"); let resolution_loc = ctx.get_uniform_location(
let shortdim_loc = ctx.get_uniform_location(&sphere_program, "shortdim"); &sphere_program, "resolution"
let layer_threshold_loc = ctx.get_uniform_location(&sphere_program, "layer_threshold"); );
let debug_mode_loc = ctx.get_uniform_location(&sphere_program, "debug_mode"); let shortdim_loc = ctx.get_uniform_location(
&sphere_program, "shortdim"
);
let layer_threshold_loc = ctx.get_uniform_location(
&sphere_program, "layer_threshold"
);
let debug_mode_loc = ctx.get_uniform_location(
&sphere_program, "debug_mode"
);
// load the viewport vertex positions into a new vertex buffer object // load the viewport vertex positions into a new vertex buffer object
const VERTEX_CNT: usize = 6; const VERTEX_CNT: usize = 6;
let viewport_positions: [f32; 3*VERTEX_CNT] = [ let viewport_positions: [f32; 3*VERTEX_CNT] = [
@ -526,21 +550,26 @@ pub fn Display() -> View {
1.0, 1.0, 0.0, 1.0, 1.0, 0.0,
1.0, -1.0, 0.0, 1.0, -1.0, 0.0,
]; ];
let viewport_position_buffer = load_new_buffer(&ctx, &viewport_positions); let viewport_position_buffer =
load_new_buffer(&ctx, &viewport_positions);
// find the point program's vertex attributes // find the point program's vertex attributes
let point_position_attr = ctx.get_attrib_location(&point_program, "position") as u32; let point_position_attr =
let point_color_attr = ctx.get_attrib_location(&point_program, "color") as u32; ctx.get_attrib_location(&point_program, "position") as u32;
let point_highlight_attr = ctx.get_attrib_location(&point_program, "highlight") as u32; let point_color_attr =
let point_selection_attr = ctx.get_attrib_location(&point_program, "selected") as u32; ctx.get_attrib_location(&point_program, "color") as u32;
let point_highlight_attr =
ctx.get_attrib_location(&point_program, "highlight") as u32;
let point_selection_attr =
ctx.get_attrib_location(&point_program, "selected") as u32;
// set up a repainting routine // set up a repainting routine
let (_, start_animation_loop, _) = create_raf(move || { let (_, start_animation_loop, _) = create_raf(move || {
// get the time step // get the time step
let time = performance.now(); let time = performance.now();
let time_step = 0.001*(time - last_time); let time_step = 0.001*(time - last_time);
last_time = time; last_time = time;
// get the navigation state // get the navigation state
let pitch_up_val = pitch_up.get(); let pitch_up_val = pitch_up.get();
let pitch_down_val = pitch_down.get(); let pitch_down_val = pitch_down.get();
@ -551,7 +580,7 @@ pub fn Display() -> View {
let zoom_in_val = zoom_in.get(); let zoom_in_val = zoom_in.get();
let zoom_out_val = zoom_out.get(); let zoom_out_val = zoom_out.get();
let turntable_val = turntable.get(); /* BENCHMARKING */ let turntable_val = turntable.get(); /* BENCHMARKING */
// get the manipulation state // get the manipulation state
let translate_neg_x_val = translate_neg_x.get(); let translate_neg_x_val = translate_neg_x.get();
let translate_pos_x_val = translate_pos_x.get(); let translate_pos_x_val = translate_pos_x.get();
@ -561,7 +590,7 @@ pub fn Display() -> View {
let translate_pos_z_val = translate_pos_z.get(); let translate_pos_z_val = translate_pos_z.get();
let shrink_neg_val = shrink_neg.get(); let shrink_neg_val = shrink_neg.get();
let shrink_pos_val = shrink_pos.get(); let shrink_pos_val = shrink_pos.get();
// update the assembly's orientation // update the assembly's orientation
let ang_vel = { let ang_vel = {
let pitch = pitch_up_val - pitch_down_val; let pitch = pitch_up_val - pitch_down_val;
@ -582,11 +611,11 @@ pub fn Display() -> View {
Rotation3::from_scaled_axis(time_step * ang_vel).matrix() Rotation3::from_scaled_axis(time_step * ang_vel).matrix()
); );
orientation = &rotation * &orientation; orientation = &rotation * &orientation;
// update the assembly's location // update the assembly's location
let zoom = zoom_out_val - zoom_in_val; let zoom = zoom_out_val - zoom_in_val;
location_z *= (time_step * ZOOM_SPEED * zoom).exp(); location_z *= (time_step * ZOOM_SPEED * zoom).exp();
// manipulate the assembly // manipulate the assembly
/* KLUDGE */ /* KLUDGE */
// to avoid the complexity of making tangent space projection // to avoid the complexity of making tangent space projection
@ -596,7 +625,8 @@ pub fn Display() -> View {
let realization_successful = state.assembly.realization_status.with( let realization_successful = state.assembly.realization_status.with(
|status| status.is_ok() |status| status.is_ok()
); );
let step_val = state.assembly.step.with_untracked(|step| step.value); let step_val =
state.assembly.step.with_untracked(|step| step.value);
let on_init_step = step_val.is_some_and(|n| n == 0.0); let on_init_step = step_val.is_some_and(|n| n == 0.0);
let on_last_step = step_val.is_some_and( let on_last_step = step_val.is_some_and(
|n| state.assembly.descent_history.with_untracked( |n| state.assembly.descent_history.with_untracked(
@ -606,7 +636,8 @@ pub fn Display() -> View {
let on_manipulable_step = let on_manipulable_step =
!realization_successful && on_init_step !realization_successful && on_init_step
|| realization_successful && on_last_step; || realization_successful && on_last_step;
if on_manipulable_step && state.selection.with(|sel| sel.len() == 1) { if on_manipulable_step
&& state.selection.with(|sel| sel.len() == 1) {
let sel = state.selection.with( let sel = state.selection.with(
|sel| sel.into_iter().next().unwrap().clone() |sel| sel.into_iter().next().unwrap().clone()
); );
@ -642,24 +673,25 @@ pub fn Display() -> View {
scene_changed.set(true); scene_changed.set(true);
} }
} }
if scene_changed.get() { if scene_changed.get() {
const SPACE_DIM: usize = 3; const SPACE_DIM: usize = 3;
const COLOR_SIZE: usize = 3; const COLOR_SIZE: usize = 3;
/* INSTRUMENTS */ /* INSTRUMENTS */
// measure mean frame interval // measure mean frame interval
frames_since_last_sample += 1; frames_since_last_sample += 1;
if frames_since_last_sample >= SAMPLE_PERIOD { if frames_since_last_sample >= SAMPLE_PERIOD {
mean_frame_interval.set((time - last_sample_time) / (SAMPLE_PERIOD as f64)); mean_frame_interval.set(
(time - last_sample_time) / (SAMPLE_PERIOD as f64));
last_sample_time = time; last_sample_time = time;
frames_since_last_sample = 0; frames_since_last_sample = 0;
} }
// --- get the assembly --- // --- get the assembly ---
let mut scene = Scene::new(); let mut scene = Scene::new();
// find the map from assembly space to world space // find the map from assembly space to world space
let location = { let location = {
let u = -location_z; let u = -location_z;
@ -672,35 +704,37 @@ pub fn Display() -> View {
]) ])
}; };
let asm_to_world = &location * &orientation; let asm_to_world = &location * &orientation;
// set up the scene // set up the scene
state.assembly.elements.with_untracked( state.assembly.elements.with_untracked(
|elts| for elt in elts { |elts| for elt in elts {
let selected = state.selection.with(|sel| sel.contains(elt)); let selected =
state.selection.with(|sel| sel.contains(elt));
elt.show(&mut scene, selected); elt.show(&mut scene, selected);
} }
); );
let sphere_cnt = scene.spheres.len_i32(); let sphere_cnt = scene.spheres.len_i32();
// --- draw the spheres --- // --- draw the spheres ---
// use the sphere rendering program // use the sphere rendering program
ctx.use_program(Some(&sphere_program)); ctx.use_program(Some(&sphere_program));
// enable the sphere program's vertex attribute // enable the sphere program's vertex attribute
ctx.enable_vertex_attrib_array(viewport_position_attr); ctx.enable_vertex_attrib_array(viewport_position_attr);
// write the spheres in world coordinates // write the spheres in world coordinates
let sphere_reps_world: Vec<_> = scene.spheres.representations.into_iter().map( let sphere_reps_world: Vec<_> =
|rep| (&asm_to_world * rep).cast::<f32>() scene.spheres.representations.into_iter().map(
).collect(); |rep| (&asm_to_world * rep).cast::<f32>()
).collect();
// set the resolution // set the resolution
let width = canvas.width() as f32; let width = canvas.width() as f32;
let height = canvas.height() as f32; let height = canvas.height() as f32;
ctx.uniform2f(resolution_loc.as_ref(), width, height); ctx.uniform2f(resolution_loc.as_ref(), width, height);
ctx.uniform1f(shortdim_loc.as_ref(), width.min(height)); ctx.uniform1f(shortdim_loc.as_ref(), width.min(height));
// pass the scene data // pass the scene data
ctx.uniform1i(sphere_cnt_loc.as_ref(), sphere_cnt); ctx.uniform1i(sphere_cnt_loc.as_ref(), sphere_cnt);
for n in 0..sphere_reps_world.len() { for n in 0..sphere_reps_world.len() {
@ -722,33 +756,35 @@ pub fn Display() -> View {
scene.spheres.highlights[n], scene.spheres.highlights[n],
); );
} }
// pass the display parameters // pass the display parameters
ctx.uniform1i(layer_threshold_loc.as_ref(), LAYER_THRESHOLD); ctx.uniform1i(layer_threshold_loc.as_ref(), LAYER_THRESHOLD);
ctx.uniform1i(debug_mode_loc.as_ref(), DEBUG_MODE); ctx.uniform1i(debug_mode_loc.as_ref(), DEBUG_MODE);
// bind the viewport vertex position buffer to the position // bind the viewport vertex position buffer to the position
// attribute in the vertex shader // attribute in the vertex shader
bind_to_attribute(&ctx, viewport_position_attr, SPACE_DIM as i32, &viewport_position_buffer); bind_to_attribute(&ctx, viewport_position_attr,
SPACE_DIM as i32, &viewport_position_buffer);
// draw the scene // draw the scene
ctx.draw_arrays(WebGl2RenderingContext::TRIANGLES, 0, VERTEX_CNT as i32); ctx.draw_arrays(WebGl2RenderingContext::TRIANGLES, 0,
VERTEX_CNT as i32);
// disable the sphere program's vertex attribute // disable the sphere program's vertex attribute
ctx.disable_vertex_attrib_array(viewport_position_attr); ctx.disable_vertex_attrib_array(viewport_position_attr);
// --- draw the points --- // --- draw the points ---
if !scene.points.representations.is_empty() { if !scene.points.representations.is_empty() {
// use the point rendering program // use the point rendering program
ctx.use_program(Some(&point_program)); ctx.use_program(Some(&point_program));
// enable the point program's vertex attributes // enable the point program's vertex attributes
ctx.enable_vertex_attrib_array(point_position_attr); ctx.enable_vertex_attrib_array(point_position_attr);
ctx.enable_vertex_attrib_array(point_color_attr); ctx.enable_vertex_attrib_array(point_color_attr);
ctx.enable_vertex_attrib_array(point_highlight_attr); ctx.enable_vertex_attrib_array(point_highlight_attr);
ctx.enable_vertex_attrib_array(point_selection_attr); ctx.enable_vertex_attrib_array(point_selection_attr);
// write the points in world coordinates // write the points in world coordinates
let asm_to_world_sp = asm_to_world.rows(0, SPACE_DIM); let asm_to_world_sp = asm_to_world.rows(0, SPACE_DIM);
let point_positions = DMatrix::from_columns( let point_positions = DMatrix::from_columns(
@ -756,30 +792,36 @@ pub fn Display() -> View {
|rep| &asm_to_world_sp * rep |rep| &asm_to_world_sp * rep
).collect::<Vec<_>>().as_slice() ).collect::<Vec<_>>().as_slice()
).cast::<f32>(); ).cast::<f32>();
// load the point positions and colors into new buffers and // load the point positions and colors into new buffers and
// bind them to the corresponding attributes in the vertex // bind them to the corresponding attributes in the vertex
// shader // shader
bind_new_buffer_to_attribute(&ctx, point_position_attr, SPACE_DIM as i32, point_positions.as_slice()); bind_new_buffer_to_attribute(&ctx, point_position_attr,
bind_new_buffer_to_attribute(&ctx, point_color_attr, (COLOR_SIZE + 1) as i32, scene.points.colors_with_opacity.concat().as_slice()); SPACE_DIM as i32, point_positions.as_slice());
bind_new_buffer_to_attribute(&ctx, point_highlight_attr, 1 as i32, scene.points.highlights.as_slice()); bind_new_buffer_to_attribute(&ctx, point_color_attr,
bind_new_buffer_to_attribute(&ctx, point_selection_attr, 1 as i32, scene.points.selections.as_slice()); (COLOR_SIZE + 1) as i32,
scene.points.colors_with_opacity.concat().as_slice());
bind_new_buffer_to_attribute(&ctx, point_highlight_attr,
1i32, scene.points.highlights.as_slice());
bind_new_buffer_to_attribute(&ctx, point_selection_attr,
1i32, scene.points.selections.as_slice());
// draw the scene // draw the scene
ctx.draw_arrays(WebGl2RenderingContext::POINTS, 0, point_positions.ncols() as i32); ctx.draw_arrays(WebGl2RenderingContext::POINTS, 0,
point_positions.ncols() as i32);
// disable the point program's vertex attributes // disable the point program's vertex attributes
ctx.disable_vertex_attrib_array(point_position_attr); ctx.disable_vertex_attrib_array(point_position_attr);
ctx.disable_vertex_attrib_array(point_color_attr); ctx.disable_vertex_attrib_array(point_color_attr);
ctx.disable_vertex_attrib_array(point_highlight_attr); ctx.disable_vertex_attrib_array(point_highlight_attr);
ctx.disable_vertex_attrib_array(point_selection_attr); ctx.disable_vertex_attrib_array(point_selection_attr);
} }
// --- update the display state --- // --- update the display state ---
// update the viewpoint // update the viewpoint
assembly_to_world.set(asm_to_world); assembly_to_world.set(asm_to_world);
// clear the scene change flag // clear the scene change flag
scene_changed.set( scene_changed.set(
pitch_up_val != 0.0 pitch_up_val != 0.0
@ -799,7 +841,7 @@ pub fn Display() -> View {
}); });
start_animation_loop(); 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 mut navigating = true;
let shift = event.shift_key(); let shift = event.shift_key();
@ -819,7 +861,7 @@ pub fn Display() -> View {
event.prevent_default(); event.prevent_default();
} }
}; };
let set_manip_signal = move |event: &KeyboardEvent, value: f64| { let set_manip_signal = move |event: &KeyboardEvent, value: f64| {
let mut manipulating = true; let mut manipulating = true;
let shift = event.shift_key(); let shift = event.shift_key();
@ -838,7 +880,7 @@ pub fn Display() -> View {
event.prevent_default(); event.prevent_default();
} }
}; };
view! { view! {
/* TO DO */ /* TO DO */
// switch back to integer-valued parameters when that becomes possible // switch back to integer-valued parameters when that becomes possible
@ -860,7 +902,7 @@ pub fn Display() -> View {
yaw_left.set(0.0); yaw_left.set(0.0);
pitch_up.set(0.0); pitch_up.set(0.0);
pitch_down.set(0.0); pitch_down.set(0.0);
// swap manipulation inputs // swap manipulation inputs
translate_pos_z.set(translate_neg_y.get()); translate_pos_z.set(translate_neg_y.get());
translate_neg_z.set(translate_pos_y.get()); translate_neg_z.set(translate_pos_y.get());
@ -886,7 +928,7 @@ pub fn Display() -> View {
roll_ccw.set(0.0); roll_ccw.set(0.0);
zoom_in.set(0.0); zoom_in.set(0.0);
zoom_out.set(0.0); zoom_out.set(0.0);
// swap manipulation inputs // swap manipulation inputs
translate_pos_y.set(translate_neg_z.get()); translate_pos_y.set(translate_neg_z.get());
translate_neg_y.set(translate_pos_z.get()); translate_neg_y.set(translate_pos_z.get());
@ -915,7 +957,9 @@ pub fn Display() -> View {
.into_iter() .into_iter()
.filter(|elt| !elt.ghost().get()); .filter(|elt| !elt.ghost().get());
for elt in tangible_elts { for elt in tangible_elts {
match assembly_to_world.with(|asm_to_world| elt.cast(dir, asm_to_world, pixel_size)) { let target = assembly_to_world.with(
|asm_to_world| elt.cast(dir, asm_to_world, pixel_size));
match target {
Some(depth) => match clicked { Some(depth) => match clicked {
Some((_, best_depth)) => { Some((_, best_depth)) => {
if depth < best_depth { if depth < best_depth {
@ -927,7 +971,7 @@ pub fn Display() -> View {
None => (), None => (),
}; };
} }
// if we clicked something, select it // if we clicked something, select it
match clicked { match clicked {
Some((elt, _)) => state.select(&elt, event.shift_key()), Some((elt, _)) => state.select(&elt, event.shift_key()),
@ -936,4 +980,4 @@ pub fn Display() -> View {
}, },
) )
} }
} }

37
app-proto/src/components/outline.rs
View file

@ -21,16 +21,16 @@ fn RegulatorInput(regulator: Rc<dyn Regulator>) -> View {
// get the regulator's measurement and set point signals // get the regulator's measurement and set point signals
let measurement = regulator.measurement(); let measurement = regulator.measurement();
let set_point = regulator.set_point(); let set_point = regulator.set_point();
// the `valid` signal tracks whether the last entered value is a valid set // the `valid` signal tracks whether the last entered value is a valid set
// point specification // point specification
let valid = create_signal(true); let valid = create_signal(true);
// the `value` signal holds the current set point specification // the `value` signal holds the current set point specification
let value = create_signal( let value = create_signal(
set_point.with_untracked(|set_pt| set_pt.spec.clone()) set_point.with_untracked(|set_pt| set_pt.spec.clone())
); );
// this `reset_value` closure resets the input value to the regulator's set // this `reset_value` closure resets the input value to the regulator's set
// point specification // point specification
let reset_value = move || { let reset_value = move || {
@ -39,11 +39,11 @@ fn RegulatorInput(regulator: Rc<dyn Regulator>) -> View {
value.set(set_point.with(|set_pt| set_pt.spec.clone())); value.set(set_point.with(|set_pt| set_pt.spec.clone()));
}) })
}; };
// reset the input value whenever the regulator's set point specification // reset the input value whenever the regulator's set point specification
// is updated // is updated
create_effect(reset_value); create_effect(reset_value);
view! { view! {
input( input(
r#type = "text", r#type = "text",
@ -63,8 +63,10 @@ fn RegulatorInput(regulator: Rc<dyn Regulator>) -> View {
placeholder = measurement.with(|result| result.to_string()), placeholder = measurement.with(|result| result.to_string()),
bind:value = value, bind:value = value,
on:change = move |_| { on:change = move |_| {
let specification =
SpecifiedValue::try_from(value.get_clone_untracked());
valid.set( valid.set(
match SpecifiedValue::try_from(value.get_clone_untracked()) { match specification {
Ok(set_pt) => { Ok(set_pt) => {
set_point.set(set_pt); set_point.set(set_pt);
true true
@ -141,7 +143,9 @@ fn ElementOutlineItem(element: Rc<dyn Element>) -> View {
let class = { let class = {
let element_for_class = element.clone(); let element_for_class = element.clone();
state.selection.map( state.selection.map(
move |sel| if sel.contains(&element_for_class) { "selected" } else { "" } move |sel|
if sel.contains(&element_for_class) { "selected" }
else { "" }
) )
}; };
let label = element.label().clone(); let label = element.label().clone();
@ -175,7 +179,8 @@ fn ElementOutlineItem(element: Rc<dyn Element>) -> View {
move |event: KeyboardEvent| { move |event: KeyboardEvent| {
match event.key().as_str() { match event.key().as_str() {
"Enter" => { "Enter" => {
state.select(&element_for_handler, event.shift_key()); state.select(&element_for_handler,
event.shift_key());
event.prevent_default(); event.prevent_default();
}, },
"ArrowRight" if regulated.get() => { "ArrowRight" if regulated.get() => {
@ -205,18 +210,22 @@ fn ElementOutlineItem(element: Rc<dyn Element>) -> View {
let state_for_handler = state.clone(); let state_for_handler = state.clone();
let element_for_handler = element.clone(); let element_for_handler = element.clone();
move |event: MouseEvent| { move |event: MouseEvent| {
state_for_handler.select(&element_for_handler, event.shift_key()); state_for_handler.select(&element_for_handler,
event.shift_key());
event.stop_propagation(); event.stop_propagation();
event.prevent_default(); event.prevent_default();
} }
} }
) { ) {
div(class = "element-label") { (label) } div(class = "element-label") { (label) }
div(class = "element-representation") { (rep_components) } div(class = "element-representation") {
(rep_components)
}
input( input(
r#type = "checkbox", r#type = "checkbox",
bind:checked = element.ghost(), bind:checked = element.ghost(),
on:click = |event: MouseEvent| event.stop_propagation() on:click =
|event: MouseEvent| event.stop_propagation()
) )
} }
} }
@ -241,7 +250,7 @@ fn ElementOutlineItem(element: Rc<dyn Element>) -> View {
#[component] #[component]
pub fn Outline() -> View { pub fn Outline() -> View {
let state = use_context::<AppState>(); let state = use_context::<AppState>();
// list the elements alphabetically by ID // list the elements alphabetically by ID
/* TO DO */ /* TO DO */
// this code is designed to generalize easily to other sort keys. if we only // this code is designed to generalize easily to other sort keys. if we only
@ -254,7 +263,7 @@ pub fn Outline() -> View {
.sorted_by_key(|elt| elt.id().clone()) .sorted_by_key(|elt| elt.id().clone())
.collect::<Vec<_>>() .collect::<Vec<_>>()
); );
view! { view! {
ul( ul(
id = "outline", id = "outline",
@ -272,4 +281,4 @@ pub fn Outline() -> View {
) )
} }
} }
} }

4
app-proto/src/components/point.frag
View file

@ -10,10 +10,10 @@ out vec4 outColor;
void main() { void main() {
float r = total_radius * length(2.*gl_PointCoord - vec2(1.)); float r = total_radius * length(2.*gl_PointCoord - vec2(1.));
const float POINT_RADIUS = 4.; const float POINT_RADIUS = 4.;
float border = smoothstep(POINT_RADIUS - 1., POINT_RADIUS, r); float border = smoothstep(POINT_RADIUS - 1., POINT_RADIUS, r);
float disk = 1. - smoothstep(total_radius - 1., total_radius, r); float disk = 1. - smoothstep(total_radius - 1., total_radius, r);
vec4 color = mix(point_color, vec4(1.), border * point_highlight); vec4 color = mix(point_color, vec4(1.), border * point_highlight);
outColor = vec4(vec3(1.), disk) * color; outColor = vec4(vec3(1.), disk) * color;
} }

6
app-proto/src/components/point.vert
View file

@ -14,11 +14,11 @@ const float focal_slope = 0.3;
void main() { void main() {
total_radius = 5. + 0.5*selected; total_radius = 5. + 0.5*selected;
float depth = -focal_slope * position.z; float depth = -focal_slope * position.z;
gl_Position = vec4(position.xy / depth, 0., 1.); gl_Position = vec4(position.xy / depth, 0., 1.);
gl_PointSize = 2.*total_radius; gl_PointSize = 2.*total_radius;
point_color = color; point_color = color;
point_highlight = highlight; point_highlight = highlight;
} }

36
app-proto/src/components/spheres.frag
View file

@ -75,7 +75,7 @@ Fragment sphere_shading(vecInv v, vec3 pt, vec4 base_color) {
// point. i calculated it in my head and decided that the result looked good // point. i calculated it in my head and decided that the result looked good
// enough for now // enough for now
vec3 normal = normalize(-v.sp + 2.*v.lt.s*pt); vec3 normal = normalize(-v.sp + 2.*v.lt.s*pt);
float incidence = dot(normal, light_dir); float incidence = dot(normal, light_dir);
float illum = mix(0.4, 1.0, max(incidence, 0.0)); float illum = mix(0.4, 1.0, max(incidence, 0.0));
return Fragment(pt, normal, vec4(illum * base_color.rgb, base_color.a)); return Fragment(pt, normal, vec4(illum * base_color.rgb, base_color.a));
@ -110,7 +110,7 @@ vec2 sphere_cast(vecInv v, vec3 dir) {
float a = -v.lt.s * dot(dir, dir); float a = -v.lt.s * dot(dir, dir);
float b = dot(v.sp, dir); float b = dot(v.sp, dir);
float c = -v.lt.t; float c = -v.lt.t;
float adjust = 4.*a*c/(b*b); float adjust = 4.*a*c/(b*b);
if (adjust < 1.) { if (adjust < 1.) {
// as long as `b` is non-zero, the linear approximation of // as long as `b` is non-zero, the linear approximation of
@ -136,7 +136,7 @@ vec2 sphere_cast(vecInv v, vec3 dir) {
void main() { void main() {
vec2 scr = (2.*gl_FragCoord.xy - resolution) / shortdim; vec2 scr = (2.*gl_FragCoord.xy - resolution) / shortdim;
vec3 dir = vec3(focal_slope * scr, -1.); vec3 dir = vec3(focal_slope * scr, -1.);
// cast rays through the spheres // cast rays through the spheres
const int LAYER_MAX = 12; const int LAYER_MAX = 12;
TaggedDepth top_hits [LAYER_MAX]; TaggedDepth top_hits [LAYER_MAX];
@ -144,7 +144,7 @@ void main() {
for (int id = 0; id < sphere_cnt; ++id) { for (int id = 0; id < sphere_cnt; ++id) {
// find out where the ray hits the sphere // find out where the ray hits the sphere
vec2 hit_depths = sphere_cast(sphere_list[id], dir); vec2 hit_depths = sphere_cast(sphere_list[id], dir);
// insertion-sort the points we hit into the hit list // insertion-sort the points we hit into the hit list
float dimming = 1.; float dimming = 1.;
for (int side = 0; side < 2; ++side) { for (int side = 0; side < 2; ++side) {
@ -169,14 +169,15 @@ void main() {
} }
} }
} }
/* DEBUG */ /* DEBUG */
// in debug mode, show the layer count instead of the shaded image // in debug mode, show the layer count instead of the shaded image
if (debug_mode) { if (debug_mode) {
// at the bottom of the screen, show the color scale instead of the // at the bottom of the screen, show the color scale instead of the
// layer count // layer count
if (gl_FragCoord.y < 10.) layer_cnt = int(16. * gl_FragCoord.x / resolution.x); if (gl_FragCoord.y < 10.) {
layer_cnt = int(16. * gl_FragCoord.x / resolution.x);
}
// convert number to color // convert number to color
ivec3 bits = layer_cnt / ivec3(1, 2, 4); ivec3 bits = layer_cnt / ivec3(1, 2, 4);
vec3 color = mod(vec3(bits), 2.); vec3 color = mod(vec3(bits), 2.);
@ -186,7 +187,7 @@ void main() {
outColor = vec4(color, 1.); outColor = vec4(color, 1.);
return; return;
} }
// composite the sphere fragments // composite the sphere fragments
vec3 color = vec3(0.); vec3 color = vec3(0.);
int layer = layer_cnt - 1; int layer = layer_cnt - 1;
@ -203,7 +204,7 @@ void main() {
// load the current fragment // load the current fragment
Fragment frag = frag_next; Fragment frag = frag_next;
float highlight = highlight_next; float highlight = highlight_next;
// shade the next fragment // shade the next fragment
hit = top_hits[layer]; hit = top_hits[layer];
sphere_color = color_list[hit.id]; sphere_color = color_list[hit.id];
@ -213,23 +214,26 @@ void main() {
vec4(hit.dimming * sphere_color.rgb, sphere_color.a) vec4(hit.dimming * sphere_color.rgb, sphere_color.a)
); );
highlight_next = highlight_list[hit.id]; highlight_next = highlight_list[hit.id];
// highlight intersections // highlight intersections
float ixn_dist = intersection_dist(frag, frag_next); float ixn_dist = intersection_dist(frag, frag_next);
float max_highlight = max(highlight, highlight_next); float max_highlight = max(highlight, highlight_next);
float ixn_highlight = 0.5 * max_highlight * (1. - smoothstep(2./3.*ixn_threshold, 1.5*ixn_threshold, ixn_dist)); float ixn_highlight = 0.5 * max_highlight * (1. - smoothstep(
2./3.*ixn_threshold, 1.5*ixn_threshold, ixn_dist));
frag.color = mix(frag.color, vec4(1.), ixn_highlight); frag.color = mix(frag.color, vec4(1.), ixn_highlight);
frag_next.color = mix(frag_next.color, vec4(1.), ixn_highlight); frag_next.color = mix(frag_next.color, vec4(1.), ixn_highlight);
// highlight cusps // highlight cusps
float cusp_cos = abs(dot(dir, frag.normal)); float cusp_cos = abs(dot(dir, frag.normal));
float cusp_threshold = 2.*sqrt(ixn_threshold * sphere_list[hit.id].lt.s); float cusp_threshold = 2.*sqrt(
float cusp_highlight = highlight * (1. - smoothstep(2./3.*cusp_threshold, 1.5*cusp_threshold, cusp_cos)); ixn_threshold * sphere_list[hit.id].lt.s);
float cusp_highlight = highlight * (1. - smoothstep(
2./3.*cusp_threshold, 1.5*cusp_threshold, cusp_cos));
frag.color = mix(frag.color, vec4(1.), cusp_highlight); frag.color = mix(frag.color, vec4(1.), cusp_highlight);
// composite the current fragment // composite the current fragment
color = mix(color, frag.color.rgb, frag.color.a); color = mix(color, frag.color.rgb, frag.color.a);
} }
color = mix(color, frag_next.color.rgb, frag_next.color.a); color = mix(color, frag_next.color.rgb, frag_next.color.a);
outColor = vec4(sRGB(color), 1.); outColor = vec4(sRGB(color), 1.);
} }

215
app-proto/src/components/test_assembly_chooser.rs
View file

@ -144,7 +144,7 @@ fn load_low_curvature(assembly: &Assembly) {
engine::sphere(2.0/3.0, 4.0/3.0 * a, 0.0, 1.0/3.0), engine::sphere(2.0/3.0, 4.0/3.0 * a, 0.0, 1.0/3.0),
) )
); );
// impose the desired tangencies and make the sides planar // impose the desired tangencies and make the sides planar
let index_range = 1..=3; let index_range = 1..=3;
let [central, assemb_plane] = ["central", "assemb_plane"].map( let [central, assemb_plane] = ["central", "assemb_plane"].map(
@ -167,29 +167,36 @@ fn load_low_curvature(assembly: &Assembly) {
let curvature = plane.regulators().with_untracked( let curvature = plane.regulators().with_untracked(
|regs| regs.first().unwrap().clone() |regs| regs.first().unwrap().clone()
); );
curvature.set_point().set(SpecifiedValue::try_from("0".to_string()).unwrap()); curvature.set_point().set(
SpecifiedValue::try_from("0".to_string()).unwrap());
} }
let all_perpendicular = [central.clone()].into_iter() let all_perpendicular = [central.clone()].into_iter()
.chain(sides.clone()) .chain(sides.clone())
.chain(corners.clone()); .chain(corners.clone());
for sphere in all_perpendicular { for sphere in all_perpendicular {
// make each side and packed sphere perpendicular to the assembly plane // make each side and packed sphere perpendicular to the assembly plane
let right_angle = InversiveDistanceRegulator::new([sphere, assemb_plane.clone()]); let right_angle = InversiveDistanceRegulator::new(
right_angle.set_point.set(SpecifiedValue::try_from("0".to_string()).unwrap()); [sphere, assemb_plane.clone()]);
right_angle.set_point.set(
SpecifiedValue::try_from("0".to_string()).unwrap());
assembly.insert_regulator(Rc::new(right_angle)); assembly.insert_regulator(Rc::new(right_angle));
} }
for sphere in sides.clone().chain(corners.clone()) { for sphere in sides.clone().chain(corners.clone()) {
// make each side and corner sphere tangent to the central sphere // make each side and corner sphere tangent to the central sphere
let tangency = InversiveDistanceRegulator::new([sphere.clone(), central.clone()]); let tangency = InversiveDistanceRegulator::new(
tangency.set_point.set(SpecifiedValue::try_from("-1".to_string()).unwrap()); [sphere.clone(), central.clone()]);
tangency.set_point.set(
SpecifiedValue::try_from("-1".to_string()).unwrap());
assembly.insert_regulator(Rc::new(tangency)); assembly.insert_regulator(Rc::new(tangency));
} }
for (side_index, side) in sides.enumerate() { for (side_index, side) in sides.enumerate() {
// make each side tangent to the two adjacent corner spheres // make each side tangent to the two adjacent corner spheres
for (corner_index, corner) in corners.clone().enumerate() { for (corner_index, corner) in corners.clone().enumerate() {
if side_index != corner_index { if side_index != corner_index {
let tangency = InversiveDistanceRegulator::new([side.clone(), corner]); let tangency = InversiveDistanceRegulator::new(
tangency.set_point.set(SpecifiedValue::try_from("-1".to_string()).unwrap()); [side.clone(), corner]);
tangency.set_point.set(
SpecifiedValue::try_from("-1".to_string()).unwrap());
assembly.insert_regulator(Rc::new(tangency)); assembly.insert_regulator(Rc::new(tangency));
} }
} }
@ -217,21 +224,24 @@ fn load_pointed(assembly: &Assembly) {
for index_y in 0..=1 { for index_y in 0..=1 {
let x = index_x as f64 - 0.5; let x = index_x as f64 - 0.5;
let y = index_y as f64 - 0.5; let y = index_y as f64 - 0.5;
let x32 = x as f32;
let y32 = y as f32;
let coords =
[0.5*(1.0 + x32), 0.5*(1.0 + y32), 0.5*(1.0 - x32*y32)];
let _ = assembly.try_insert_element( let _ = assembly.try_insert_element(
Sphere::new( Sphere::new(
format!("sphere{index_x}{index_y}"), format!("sphere{index_x}{index_y}"),
format!("Sphere {index_x}{index_y}"), format!("Sphere {index_x}{index_y}"),
[0.5*(1.0 + x) as f32, 0.5*(1.0 + y) as f32, 0.5*(1.0 - x*y) as f32], coords,
engine::sphere(x, y, 0.0, 1.0), engine::sphere(x, y, 0.0, 1.0),
) )
); );
let _ = assembly.try_insert_element( let _ = assembly.try_insert_element(
Point::new( Point::new(
format!("point{index_x}{index_y}"), format!("point{index_x}{index_y}"),
format!("Point {index_x}{index_y}"), format!("Point {index_x}{index_y}"),
[0.5*(1.0 + x) as f32, 0.5*(1.0 + y) as f32, 0.5*(1.0 - x*y) as f32], coords,
engine::point(x, y, 0.0), engine::point(x, y, 0.0),
) )
); );
@ -310,7 +320,7 @@ fn load_tridiminished_icosahedron(assembly: &Assembly) {
for vertex in vertices { for vertex in vertices {
let _ = assembly.try_insert_element(vertex); let _ = assembly.try_insert_element(vertex);
} }
// create the faces // create the faces
const COLOR_FACE: ElementColor = [0.75_f32, 0.75_f32, 0.75_f32]; const COLOR_FACE: ElementColor = [0.75_f32, 0.75_f32, 0.75_f32];
let frac_1_sqrt_6 = 1.0 / 6.0_f64.sqrt(); let frac_1_sqrt_6 = 1.0 / 6.0_f64.sqrt();
@ -320,26 +330,32 @@ fn load_tridiminished_icosahedron(assembly: &Assembly) {
"face1".to_string(), "face1".to_string(),
"Face 1".to_string(), "Face 1".to_string(),
COLOR_FACE, COLOR_FACE,
engine::sphere_with_offset(frac_2_sqrt_6, -frac_1_sqrt_6, -frac_1_sqrt_6, -frac_1_sqrt_6, 0.0), engine::sphere_with_offset(
frac_2_sqrt_6, -frac_1_sqrt_6, -frac_1_sqrt_6,
-frac_1_sqrt_6, 0.0),
), ),
Sphere::new( Sphere::new(
"face2".to_string(), "face2".to_string(),
"Face 2".to_string(), "Face 2".to_string(),
COLOR_FACE, COLOR_FACE,
engine::sphere_with_offset(-frac_1_sqrt_6, frac_2_sqrt_6, -frac_1_sqrt_6, -frac_1_sqrt_6, 0.0), engine::sphere_with_offset(
-frac_1_sqrt_6, frac_2_sqrt_6, -frac_1_sqrt_6,
-frac_1_sqrt_6, 0.0),
), ),
Sphere::new( Sphere::new(
"face3".to_string(), "face3".to_string(),
"Face 3".to_string(), "Face 3".to_string(),
COLOR_FACE, COLOR_FACE,
engine::sphere_with_offset(-frac_1_sqrt_6, -frac_1_sqrt_6, frac_2_sqrt_6, -frac_1_sqrt_6, 0.0), engine::sphere_with_offset(
-frac_1_sqrt_6, -frac_1_sqrt_6, frac_2_sqrt_6,
-frac_1_sqrt_6, 0.0),
), ),
]; ];
for face in faces { for face in faces {
face.ghost().set(true); face.ghost().set(true);
let _ = assembly.try_insert_element(face); let _ = assembly.try_insert_element(face);
} }
let index_range = 1..=3; let index_range = 1..=3;
for j in index_range.clone() { for j in index_range.clone() {
// make each face planar // make each face planar
@ -352,15 +368,17 @@ fn load_tridiminished_icosahedron(assembly: &Assembly) {
curvature_regulator.set_point().set( curvature_regulator.set_point().set(
SpecifiedValue::try_from("0".to_string()).unwrap() SpecifiedValue::try_from("0".to_string()).unwrap()
); );
// put each A vertex on the face it belongs to // put each A vertex on the face it belongs to
let vertex_a = assembly.elements_by_id.with_untracked( let vertex_a = assembly.elements_by_id.with_untracked(
|elts_by_id| elts_by_id[&format!("a{j}")].clone() |elts_by_id| elts_by_id[&format!("a{j}")].clone()
); );
let incidence_a = InversiveDistanceRegulator::new([face.clone(), vertex_a.clone()]); let incidence_a = InversiveDistanceRegulator::new(
incidence_a.set_point.set(SpecifiedValue::try_from("0".to_string()).unwrap()); [face.clone(), vertex_a.clone()]);
incidence_a.set_point.set(
SpecifiedValue::try_from("0".to_string()).unwrap());
assembly.insert_regulator(Rc::new(incidence_a)); assembly.insert_regulator(Rc::new(incidence_a));
// regulate the B-C vertex distances // regulate the B-C vertex distances
let vertices_bc = ["b", "c"].map( let vertices_bc = ["b", "c"].map(
|series| assembly.elements_by_id.with_untracked( |series| assembly.elements_by_id.with_untracked(
@ -370,27 +388,30 @@ fn load_tridiminished_icosahedron(assembly: &Assembly) {
assembly.insert_regulator( assembly.insert_regulator(
Rc::new(InversiveDistanceRegulator::new(vertices_bc)) Rc::new(InversiveDistanceRegulator::new(vertices_bc))
); );
// get the pair of indices adjacent to `j` // get the pair of indices adjacent to `j`
let adjacent_indices = [j % 3 + 1, (j + 1) % 3 + 1]; let adjacent_indices = [j % 3 + 1, (j + 1) % 3 + 1];
for k in adjacent_indices.clone() { for k in adjacent_indices.clone() {
for series in ["b", "c"] { for series in ["b", "c"] {
// put each B and C vertex on the faces it belongs to // put each B and C vertex on the faces it belongs to
let vertex = assembly.elements_by_id.with_untracked( let vertex = assembly.elements_by_id.with_untracked(
|elts_by_id| elts_by_id[&format!("{series}{k}")].clone() |elts_by_id| elts_by_id[&format!("{series}{k}")].clone()
); );
let incidence = InversiveDistanceRegulator::new([face.clone(), vertex.clone()]); let incidence = InversiveDistanceRegulator::new(
incidence.set_point.set(SpecifiedValue::try_from("0".to_string()).unwrap()); [face.clone(), vertex.clone()]);
incidence.set_point.set(
SpecifiedValue::try_from("0".to_string()).unwrap());
assembly.insert_regulator(Rc::new(incidence)); assembly.insert_regulator(Rc::new(incidence));
// regulate the A-B and A-C vertex distances // regulate the A-B and A-C vertex distances
assembly.insert_regulator( assembly.insert_regulator(
Rc::new(InversiveDistanceRegulator::new([vertex_a.clone(), vertex])) Rc::new(InversiveDistanceRegulator::new(
[vertex_a.clone(), vertex]))
); );
} }
} }
// regulate the A-A and C-C vertex distances // regulate the A-A and C-C vertex distances
let adjacent_pairs = ["a", "c"].map( let adjacent_pairs = ["a", "c"].map(
|series| adjacent_indices.map( |series| adjacent_indices.map(
@ -422,19 +443,20 @@ fn load_dodecahedral_packing(assembly: &Assembly) {
let substrate = assembly.elements_by_id.with_untracked( let substrate = assembly.elements_by_id.with_untracked(
|elts_by_id| elts_by_id["substrate"].clone() |elts_by_id| elts_by_id["substrate"].clone()
); );
// fix the substrate's curvature // fix the substrate's curvature
substrate.regulators().with_untracked( substrate.regulators().with_untracked(
|regs| regs.first().unwrap().clone() |regs| regs.first().unwrap().clone()
).set_point().set( ).set_point().set(
SpecifiedValue::try_from("0.5".to_string()).unwrap() SpecifiedValue::try_from("0.5".to_string()).unwrap()
); );
// add the circles to be packed // add the circles to be packed
const COLOR_A: ElementColor = [1.00_f32, 0.25_f32, 0.00_f32]; const COLOR_A: ElementColor = [1.00_f32, 0.25_f32, 0.00_f32];
const COLOR_B: ElementColor = [1.00_f32, 0.00_f32, 0.25_f32]; const COLOR_B: ElementColor = [1.00_f32, 0.00_f32, 0.25_f32];
const COLOR_C: ElementColor = [0.25_f32, 0.00_f32, 1.00_f32]; const COLOR_C: ElementColor = [0.25_f32, 0.00_f32, 1.00_f32];
let phi = 0.5 + 1.25_f64.sqrt(); /* TO DO */ // replace with std::f64::consts::PHI when that gets stabilized /* TO DO */ // replace with std::f64::consts::PHI when that gets stabilized
let phi = 0.5 + 1.25_f64.sqrt();
let phi_inv = 1.0 / phi; let phi_inv = 1.0 / phi;
let coord_scale = (phi + 2.0).sqrt(); let coord_scale = (phi + 2.0).sqrt();
let face_scales = [phi_inv, (13.0 / 12.0) / coord_scale]; let face_scales = [phi_inv, (13.0 / 12.0) / coord_scale];
@ -445,10 +467,10 @@ fn load_dodecahedral_packing(assembly: &Assembly) {
for k in 0..2 { for k in 0..2 {
let small_coord = face_scales[k] * (2.0*(j as f64) - 1.0); let small_coord = face_scales[k] * (2.0*(j as f64) - 1.0);
let big_coord = face_scales[k] * (2.0*(k as f64) - 1.0) * phi; let big_coord = face_scales[k] * (2.0*(k as f64) - 1.0) * phi;
let id_num = format!("{j}{k}"); let id_num = format!("{j}{k}");
let label_sub = format!("{}{}", subscripts[j], subscripts[k]); let label_sub = format!("{}{}", subscripts[j], subscripts[k]);
// add the A face // add the A face
let id_a = format!("a{id_num}"); let id_a = format!("a{id_num}");
let _ = assembly.try_insert_element( let _ = assembly.try_insert_element(
@ -464,7 +486,7 @@ fn load_dodecahedral_packing(assembly: &Assembly) {
|elts_by_id| elts_by_id[&id_a].clone() |elts_by_id| elts_by_id[&id_a].clone()
) )
); );
// add the B face // add the B face
let id_b = format!("b{id_num}"); let id_b = format!("b{id_num}");
let _ = assembly.try_insert_element( let _ = assembly.try_insert_element(
@ -480,7 +502,7 @@ fn load_dodecahedral_packing(assembly: &Assembly) {
|elts_by_id| elts_by_id[&id_b].clone() |elts_by_id| elts_by_id[&id_b].clone()
) )
); );
// add the C face // add the C face
let id_c = format!("c{id_num}"); let id_c = format!("c{id_num}");
let _ = assembly.try_insert_element( let _ = assembly.try_insert_element(
@ -498,16 +520,19 @@ fn load_dodecahedral_packing(assembly: &Assembly) {
); );
} }
} }
// make each face sphere perpendicular to the substrate // make each face sphere perpendicular to the substrate
for face in faces { for face in faces {
let right_angle = InversiveDistanceRegulator::new([face, substrate.clone()]); let right_angle = InversiveDistanceRegulator::new(
right_angle.set_point.set(SpecifiedValue::try_from("0".to_string()).unwrap()); [face, substrate.clone()]);
right_angle.set_point.set(
SpecifiedValue::try_from("0".to_string()).unwrap());
assembly.insert_regulator(Rc::new(right_angle)); assembly.insert_regulator(Rc::new(right_angle));
} }
// set up the tangencies that define the packing // set up the tangencies that define the packing
for [long_edge_plane, short_edge_plane] in [["a", "b"], ["b", "c"], ["c", "a"]] { for [long_edge_plane, short_edge_plane]
in [["a", "b"], ["b", "c"], ["c", "a"]] {
for k in 0..2 { for k in 0..2 {
let long_edge_ids = [ let long_edge_ids = [
format!("{long_edge_plane}{k}0"), format!("{long_edge_plane}{k}0"),
@ -524,14 +549,16 @@ fn load_dodecahedral_packing(assembly: &Assembly) {
) )
) )
); );
// set up the short-edge tangency // set up the short-edge tangency
let short_tangency = InversiveDistanceRegulator::new(short_edge.clone()); let short_tangency = InversiveDistanceRegulator::new(
short_edge.clone());
if k == 0 { if k == 0 {
short_tangency.set_point.set(SpecifiedValue::try_from("-1".to_string()).unwrap()); short_tangency.set_point.set(
SpecifiedValue::try_from("-1".to_string()).unwrap());
} }
assembly.insert_regulator(Rc::new(short_tangency)); assembly.insert_regulator(Rc::new(short_tangency));
// set up the side tangencies // set up the side tangencies
for i in 0..2 { for i in 0..2 {
for j in 0..2 { for j in 0..2 {
@ -539,7 +566,9 @@ fn load_dodecahedral_packing(assembly: &Assembly) {
[long_edge[i].clone(), short_edge[j].clone()] [long_edge[i].clone(), short_edge[j].clone()]
); );
if i == 0 && k == 0 { if i == 0 && k == 0 {
side_tangency.set_point.set(SpecifiedValue::try_from("-1".to_string()).unwrap()); side_tangency.set_point.set(
SpecifiedValue::try_from("-1".to_string()).unwrap()
);
} }
assembly.insert_regulator(Rc::new(side_tangency)); assembly.insert_regulator(Rc::new(side_tangency));
} }
@ -577,14 +606,14 @@ fn load_balanced(assembly: &Assembly) {
for sphere in spheres { for sphere in spheres {
let _ = assembly.try_insert_element(sphere); let _ = assembly.try_insert_element(sphere);
} }
// get references to the spheres // get references to the spheres
let [outer, a, b] = ["outer", "a", "b"].map( let [outer, a, b] = ["outer", "a", "b"].map(
|id| assembly.elements_by_id.with_untracked( |id| assembly.elements_by_id.with_untracked(
|elts_by_id| elts_by_id[id].clone() |elts_by_id| elts_by_id[id].clone()
) )
); );
// fix the diameters of the outer, sun, and moon spheres // fix the diameters of the outer, sun, and moon spheres
for (sphere, radius) in [ for (sphere, radius) in [
(outer.clone(), R_OUTER), (outer.clone(), R_OUTER),
@ -599,12 +628,13 @@ fn load_balanced(assembly: &Assembly) {
SpecifiedValue::try_from(curvature.to_string()).unwrap() SpecifiedValue::try_from(curvature.to_string()).unwrap()
); );
} }
// set the inversive distances between the spheres. as described above, the // set the inversive distances between the spheres. as described above, the
// initial configuration deliberately violates these constraints // initial configuration deliberately violates these constraints
for inner in [a, b] { for inner in [a, b] {
let tangency = InversiveDistanceRegulator::new([outer.clone(), inner]); let tangency = InversiveDistanceRegulator::new([outer.clone(), inner]);
tangency.set_point.set(SpecifiedValue::try_from("1".to_string()).unwrap()); tangency.set_point.set(
SpecifiedValue::try_from("1".to_string()).unwrap());
assembly.insert_regulator(Rc::new(tangency)); assembly.insert_regulator(Rc::new(tangency));
} }
} }
@ -629,14 +659,14 @@ fn load_off_center(assembly: &Assembly) {
engine::sphere(0.0, 0.0, 0.0, 1.0), engine::sphere(0.0, 0.0, 0.0, 1.0),
), ),
); );
// get references to the elements // get references to the elements
let point_and_sphere = ["point", "sphere"].map( let point_and_sphere = ["point", "sphere"].map(
|id| assembly.elements_by_id.with_untracked( |id| assembly.elements_by_id.with_untracked(
|elts_by_id| elts_by_id[id].clone() |elts_by_id| elts_by_id[id].clone()
) )
); );
// put the point on the sphere // put the point on the sphere
let incidence = InversiveDistanceRegulator::new(point_and_sphere); let incidence = InversiveDistanceRegulator::new(point_and_sphere);
incidence.set_point.set(SpecifiedValue::try_from("0".to_string()).unwrap()); incidence.set_point.set(SpecifiedValue::try_from("0".to_string()).unwrap());
@ -650,7 +680,7 @@ fn load_off_center(assembly: &Assembly) {
// inversive distance of 0 between the circumsphere and each vertex // inversive distance of 0 between the circumsphere and each vertex
fn load_radius_ratio(assembly: &Assembly) { fn load_radius_ratio(assembly: &Assembly) {
let index_range = 1..=4; let index_range = 1..=4;
// create the spheres // create the spheres
const GRAY: ElementColor = [0.75_f32, 0.75_f32, 0.75_f32]; const GRAY: ElementColor = [0.75_f32, 0.75_f32, 0.75_f32];
let spheres = [ let spheres = [
@ -670,7 +700,7 @@ fn load_radius_ratio(assembly: &Assembly) {
for sphere in spheres { for sphere in spheres {
let _ = assembly.try_insert_element(sphere); let _ = assembly.try_insert_element(sphere);
} }
// create the vertices // create the vertices
let vertices = izip!( let vertices = izip!(
index_range.clone(), index_range.clone(),
@ -699,7 +729,7 @@ fn load_radius_ratio(assembly: &Assembly) {
for vertex in vertices { for vertex in vertices {
let _ = assembly.try_insert_element(vertex); let _ = assembly.try_insert_element(vertex);
} }
// create the faces // create the faces
let base_dir = Vector3::new(1.0, 0.75, 1.0).normalize(); let base_dir = Vector3::new(1.0, 0.75, 1.0).normalize();
let offset = base_dir.dot(&Vector3::new(-0.6, 0.8, 0.6)); let offset = base_dir.dot(&Vector3::new(-0.6, 0.8, 0.6));
@ -712,10 +742,14 @@ fn load_radius_ratio(assembly: &Assembly) {
[0.25_f32, 0.00_f32, 1.00_f32], [0.25_f32, 0.00_f32, 1.00_f32],
].into_iter(), ].into_iter(),
[ [
engine::sphere_with_offset(base_dir[0], base_dir[1], base_dir[2], offset, 0.0), engine::sphere_with_offset(
engine::sphere_with_offset(base_dir[0], -base_dir[1], -base_dir[2], offset, 0.0), base_dir[0], base_dir[1], base_dir[2], offset, 0.0),
engine::sphere_with_offset(-base_dir[0], base_dir[1], -base_dir[2], offset, 0.0), engine::sphere_with_offset(
engine::sphere_with_offset(-base_dir[0], -base_dir[1], base_dir[2], offset, 0.0), base_dir[0], -base_dir[1], -base_dir[2], offset, 0.0),
engine::sphere_with_offset(
-base_dir[0], base_dir[1], -base_dir[2], offset, 0.0),
engine::sphere_with_offset(
-base_dir[0], -base_dir[1], base_dir[2], offset, 0.0),
].into_iter() ].into_iter()
).map( ).map(
|(k, color, representation)| { |(k, color, representation)| {
@ -731,7 +765,7 @@ fn load_radius_ratio(assembly: &Assembly) {
face.ghost().set(true); face.ghost().set(true);
let _ = assembly.try_insert_element(face); let _ = assembly.try_insert_element(face);
} }
// impose the constraints // impose the constraints
for j in index_range.clone() { for j in index_range.clone() {
let [face_j, vertex_j] = [ let [face_j, vertex_j] = [
@ -742,7 +776,7 @@ fn load_radius_ratio(assembly: &Assembly) {
|elts_by_id| elts_by_id[&id].clone() |elts_by_id| elts_by_id[&id].clone()
) )
); );
// make the faces planar // make the faces planar
let curvature_regulator = face_j.regulators().with_untracked( let curvature_regulator = face_j.regulators().with_untracked(
|regs| regs.first().unwrap().clone() |regs| regs.first().unwrap().clone()
@ -750,12 +784,12 @@ fn load_radius_ratio(assembly: &Assembly) {
curvature_regulator.set_point().set( curvature_regulator.set_point().set(
SpecifiedValue::try_from("0".to_string()).unwrap() SpecifiedValue::try_from("0".to_string()).unwrap()
); );
for k in index_range.clone().filter(|&index| index != j) { for k in index_range.clone().filter(|&index| index != j) {
let vertex_k = assembly.elements_by_id.with_untracked( let vertex_k = assembly.elements_by_id.with_untracked(
|elts_by_id| elts_by_id[&format!("v{k}")].clone() |elts_by_id| elts_by_id[&format!("v{k}")].clone()
); );
// fix the distances between the vertices // fix the distances between the vertices
if j < k { if j < k {
let distance_regulator = InversiveDistanceRegulator::new( let distance_regulator = InversiveDistanceRegulator::new(
@ -763,10 +797,12 @@ fn load_radius_ratio(assembly: &Assembly) {
); );
assembly.insert_regulator(Rc::new(distance_regulator)); assembly.insert_regulator(Rc::new(distance_regulator));
} }
// put the vertices on the faces // put the vertices on the faces
let incidence_regulator = InversiveDistanceRegulator::new([face_j.clone(), vertex_k.clone()]); let incidence_regulator = InversiveDistanceRegulator::new(
incidence_regulator.set_point.set(SpecifiedValue::try_from("0".to_string()).unwrap()); [face_j.clone(), vertex_k.clone()]);
incidence_regulator.set_point.set(
SpecifiedValue::try_from("0".to_string()).unwrap());
assembly.insert_regulator(Rc::new(incidence_regulator)); assembly.insert_regulator(Rc::new(incidence_regulator));
} }
} }
@ -799,7 +835,7 @@ fn load_irisawa_hexlet(assembly: &Assembly) {
[0.00_f32, 0.25_f32, 1.00_f32], [0.00_f32, 0.25_f32, 1.00_f32],
[0.25_f32, 0.00_f32, 1.00_f32], [0.25_f32, 0.00_f32, 1.00_f32],
].into_iter(); ].into_iter();
// create the spheres // create the spheres
let spheres = [ let spheres = [
Sphere::new( Sphere::new(
@ -836,7 +872,7 @@ fn load_irisawa_hexlet(assembly: &Assembly) {
for sphere in spheres { for sphere in spheres {
let _ = assembly.try_insert_element(sphere); let _ = assembly.try_insert_element(sphere);
} }
// put the outer sphere in ghost mode and fix its curvature // put the outer sphere in ghost mode and fix its curvature
let outer = assembly.elements_by_id.with_untracked( let outer = assembly.elements_by_id.with_untracked(
|elts_by_id| elts_by_id["outer"].clone() |elts_by_id| elts_by_id["outer"].clone()
@ -848,7 +884,7 @@ fn load_irisawa_hexlet(assembly: &Assembly) {
outer_curvature_regulator.set_point().set( outer_curvature_regulator.set_point().set(
SpecifiedValue::try_from((1.0 / 3.0).to_string()).unwrap() SpecifiedValue::try_from((1.0 / 3.0).to_string()).unwrap()
); );
// impose the desired tangencies // impose the desired tangencies
let [outer, sun, moon] = ["outer", "sun", "moon"].map( let [outer, sun, moon] = ["outer", "sun", "moon"].map(
|id| assembly.elements_by_id.with_untracked( |id| assembly.elements_by_id.with_untracked(
@ -860,25 +896,33 @@ fn load_irisawa_hexlet(assembly: &Assembly) {
|elts_by_id| elts_by_id[&format!("chain{k}")].clone() |elts_by_id| elts_by_id[&format!("chain{k}")].clone()
) )
); );
for (chain_sphere, chain_sphere_next) in chain.clone().zip(chain.cycle().skip(1)) { for (chain_sphere, chain_sphere_next)
in chain.clone().zip(chain.cycle().skip(1)) {
for (other_sphere, inversive_distance) in [ for (other_sphere, inversive_distance) in [
(outer.clone(), "1"), (outer.clone(), "1"),
(sun.clone(), "-1"), (sun.clone(), "-1"),
(moon.clone(), "-1"), (moon.clone(), "-1"),
(chain_sphere_next.clone(), "-1"), (chain_sphere_next.clone(), "-1"),
] { ] {
let tangency = InversiveDistanceRegulator::new([chain_sphere.clone(), other_sphere]); let tangency = InversiveDistanceRegulator::new(
tangency.set_point.set(SpecifiedValue::try_from(inversive_distance.to_string()).unwrap()); [chain_sphere.clone(), other_sphere]);
tangency.set_point.set(
SpecifiedValue::try_from(
inversive_distance.to_string()).unwrap());
assembly.insert_regulator(Rc::new(tangency)); assembly.insert_regulator(Rc::new(tangency));
} }
} }
let outer_sun_tangency = InversiveDistanceRegulator::new([outer.clone(), sun]); let outer_sun_tangency = InversiveDistanceRegulator::new(
outer_sun_tangency.set_point.set(SpecifiedValue::try_from("1".to_string()).unwrap()); [outer.clone(), sun]);
outer_sun_tangency.set_point.set(
SpecifiedValue::try_from("1".to_string()).unwrap());
assembly.insert_regulator(Rc::new(outer_sun_tangency)); assembly.insert_regulator(Rc::new(outer_sun_tangency));
let outer_moon_tangency = InversiveDistanceRegulator::new([outer.clone(), moon]); let outer_moon_tangency = InversiveDistanceRegulator::new(
outer_moon_tangency.set_point.set(SpecifiedValue::try_from("1".to_string()).unwrap()); [outer.clone(), moon]);
outer_moon_tangency.set_point.set(
SpecifiedValue::try_from("1".to_string()).unwrap());
assembly.insert_regulator(Rc::new(outer_moon_tangency)); assembly.insert_regulator(Rc::new(outer_moon_tangency));
} }
@ -895,24 +939,25 @@ pub fn TestAssemblyChooser() -> View {
console::log_1( console::log_1(
&JsValue::from(format!("Showing assembly \"{}\"", name.clone())) &JsValue::from(format!("Showing assembly \"{}\"", name.clone()))
); );
batch(|| { batch(|| {
let state = use_context::<AppState>(); let state = use_context::<AppState>();
let assembly = &state.assembly; let assembly = &state.assembly;
// clear state // clear state
assembly.regulators.update(|regs| regs.clear()); assembly.regulators.update(|regs| regs.clear());
assembly.elements.update(|elts| elts.clear()); assembly.elements.update(|elts| elts.clear());
assembly.elements_by_id.update(|elts_by_id| elts_by_id.clear()); assembly.elements_by_id.update(|elts_by_id| elts_by_id.clear());
assembly.descent_history.set(DescentHistory::new()); assembly.descent_history.set(DescentHistory::new());
state.selection.update(|sel| sel.clear()); state.selection.update(|sel| sel.clear());
// load assembly // load assembly
match name.as_str() { match name.as_str() {
"general" => load_general(assembly), "general" => load_general(assembly),
"low-curvature" => load_low_curvature(assembly), "low-curvature" => load_low_curvature(assembly),
"pointed" => load_pointed(assembly), "pointed" => load_pointed(assembly),
"tridiminished-icosahedron" => load_tridiminished_icosahedron(assembly), "tridiminished-icosahedron" =>
load_tridiminished_icosahedron(assembly),
"dodecahedral-packing" => load_dodecahedral_packing(assembly), "dodecahedral-packing" => load_dodecahedral_packing(assembly),
"balanced" => load_balanced(assembly), "balanced" => load_balanced(assembly),
"off-center" => load_off_center(assembly), "off-center" => load_off_center(assembly),
@ -922,14 +967,16 @@ pub fn TestAssemblyChooser() -> View {
}; };
}); });
}); });
// build the chooser // build the chooser
view! { view! {
select(bind:value = assembly_name) { select(bind:value = assembly_name) {
option(value = "general") { "General" } option(value = "general") { "General" }
option(value = "low-curvature") { "Low-curvature" } option(value = "low-curvature") { "Low-curvature" }
option(value = "pointed") { "Pointed" } option(value = "pointed") { "Pointed" }
option(value = "tridiminished-icosahedron") { "Tridiminished icosahedron" } option(value = "tridiminished-icosahedron") {
"Tridiminished icosahedron"
}
option(value = "dodecahedral-packing") { "Dodecahedral packing" } option(value = "dodecahedral-packing") { "Dodecahedral packing" }
option(value = "balanced") { "Balanced" } option(value = "balanced") { "Balanced" }
option(value = "off-center") { "Off-center" } option(value = "off-center") { "Off-center" }
@ -938,4 +985,4 @@ pub fn TestAssemblyChooser() -> View {
option(value = "empty") { "Empty" } option(value = "empty") { "Empty" }
} }
} }
} }

255
app-proto/src/engine.rs
View file

@ -9,8 +9,11 @@ pub fn point(x: f64, y: f64, z: f64) -> DVector<f64> {
} }
// the sphere with the given center and radius, with inward-pointing normals // the sphere with the given center and radius, with inward-pointing normals
pub fn sphere(center_x: f64, center_y: f64, center_z: f64, radius: f64) -> DVector<f64> { pub fn sphere(center_x: f64, center_y: f64, center_z: f64, radius: f64)
let center_norm_sq = center_x * center_x + center_y * center_y + center_z * center_z; -> DVector<f64>
{
let center_norm_sq =
center_x * center_x + center_y * center_y + center_z * center_z;
DVector::from_column_slice(&[ DVector::from_column_slice(&[
center_x / radius, center_x / radius,
center_y / radius, center_y / radius,
@ -23,7 +26,9 @@ pub fn sphere(center_x: f64, center_y: f64, center_z: f64, radius: f64) -> DVect
// the sphere of curvature `curv` whose closest point to the origin has position // the sphere of curvature `curv` whose closest point to the origin has position
// `off * dir` and normal `dir`, where `dir` is a unit vector. setting the // `off * dir` and normal `dir`, where `dir` is a unit vector. setting the
// curvature to zero gives a plane // curvature to zero gives a plane
pub fn sphere_with_offset(dir_x: f64, dir_y: f64, dir_z: f64, off: f64, curv: f64) -> DVector<f64> { pub fn sphere_with_offset(
dir_x: f64, dir_y: f64, dir_z: f64, off: f64, curv: f64) -> DVector<f64>
{
let norm_sp = 1.0 + off * curv; let norm_sp = 1.0 + off * curv;
DVector::from_column_slice(&[ DVector::from_column_slice(&[
norm_sp * dir_x, norm_sp * dir_x,
@ -62,19 +67,19 @@ impl PartialMatrix {
pub fn new() -> Self { pub fn new() -> Self {
Self(Vec::<MatrixEntry>::new()) Self(Vec::<MatrixEntry>::new())
} }
pub fn push(&mut self, row: usize, col: usize, value: f64) { pub fn push(&mut self, row: usize, col: usize, value: f64) {
let Self(entries) = self; let Self(entries) = self;
entries.push(MatrixEntry { index: (row, col), value }); entries.push(MatrixEntry { index: (row, col), value });
} }
pub fn push_sym(&mut self, row: usize, col: usize, value: f64) { pub fn push_sym(&mut self, row: usize, col: usize, value: f64) {
self.push(row, col, value); self.push(row, col, value);
if row != col { if row != col {
self.push(col, row, value); self.push(col, row, value);
} }
} }
fn freeze(&self, a: &DMatrix<f64>) -> DMatrix<f64> { fn freeze(&self, a: &DMatrix<f64>) -> DMatrix<f64> {
let mut result = a.clone(); let mut result = a.clone();
for &MatrixEntry { index, value } in self { for &MatrixEntry { index, value } in self {
@ -82,7 +87,7 @@ impl PartialMatrix {
} }
result result
} }
fn proj(&self, a: &DMatrix<f64>) -> DMatrix<f64> { fn proj(&self, a: &DMatrix<f64>) -> DMatrix<f64> {
let mut result = DMatrix::<f64>::zeros(a.nrows(), a.ncols()); let mut result = DMatrix::<f64>::zeros(a.nrows(), a.ncols());
for &MatrixEntry { index, .. } in self { for &MatrixEntry { index, .. } in self {
@ -90,7 +95,7 @@ impl PartialMatrix {
} }
result result
} }
fn sub_proj(&self, rhs: &DMatrix<f64>) -> DMatrix<f64> { fn sub_proj(&self, rhs: &DMatrix<f64>) -> DMatrix<f64> {
let mut result = DMatrix::<f64>::zeros(rhs.nrows(), rhs.ncols()); let mut result = DMatrix::<f64>::zeros(rhs.nrows(), rhs.ncols());
for &MatrixEntry { index, value } in self { for &MatrixEntry { index, value } in self {
@ -112,7 +117,7 @@ impl Display for PartialMatrix {
impl IntoIterator for PartialMatrix { impl IntoIterator for PartialMatrix {
type Item = MatrixEntry; type Item = MatrixEntry;
type IntoIter = std::vec::IntoIter<Self::Item>; type IntoIter = std::vec::IntoIter<Self::Item>;
fn into_iter(self) -> Self::IntoIter { fn into_iter(self) -> Self::IntoIter {
let Self(entries) = self; let Self(entries) = self;
entries.into_iter() entries.into_iter()
@ -122,7 +127,7 @@ impl IntoIterator for PartialMatrix {
impl<'a> IntoIterator for &'a PartialMatrix { impl<'a> IntoIterator for &'a PartialMatrix {
type Item = &'a MatrixEntry; type Item = &'a MatrixEntry;
type IntoIter = std::slice::Iter<'a, MatrixEntry>; type IntoIter = std::slice::Iter<'a, MatrixEntry>;
fn into_iter(self) -> Self::IntoIter { fn into_iter(self) -> Self::IntoIter {
let PartialMatrix(entries) = self; let PartialMatrix(entries) = self;
entries.into_iter() entries.into_iter()
@ -146,7 +151,7 @@ impl ConfigSubspace {
basis_std: Vec::new(), basis_std: Vec::new(),
} }
} }
// approximate the kernel of a symmetric endomorphism of the configuration // approximate the kernel of a symmetric endomorphism of the configuration
// space for `assembly_dim` elements. we consider an eigenvector to be part // space for `assembly_dim` elements. we consider an eigenvector to be part
// of the kernel if its eigenvalue is smaller than the constant `THRESHOLD` // of the kernel if its eigenvalue is smaller than the constant `THRESHOLD`
@ -167,10 +172,10 @@ impl ConfigSubspace {
|(λ, v)| (λ.abs() < THRESHOLD).then_some(v) |(λ, v)| (λ.abs() < THRESHOLD).then_some(v)
).collect::<Vec<_>>().as_slice() ).collect::<Vec<_>>().as_slice()
); );
// express the basis in the standard coordinates // express the basis in the standard coordinates
let basis_std = proj_to_std * &basis_proj; let basis_std = proj_to_std * &basis_proj;
const ELEMENT_DIM: usize = 5; const ELEMENT_DIM: usize = 5;
const UNIFORM_DIM: usize = 4; const UNIFORM_DIM: usize = 4;
Self { Self {
@ -187,20 +192,22 @@ impl ConfigSubspace {
).collect(), ).collect(),
} }
} }
pub fn dim(&self) -> usize { pub fn dim(&self) -> usize {
self.basis_std.len() self.basis_std.len()
} }
pub fn assembly_dim(&self) -> usize { pub fn assembly_dim(&self) -> usize {
self.assembly_dim self.assembly_dim
} }
// find the projection onto this subspace of the motion where the element // find the projection onto this subspace of the motion where the element
// with the given column index has velocity `v`. the velocity is given in // with the given column index has velocity `v`. the velocity is given in
// projection coordinates, and the projection is done with respect to the // projection coordinates, and the projection is done with respect to the
// projection inner product // projection inner product
pub fn proj(&self, v: &DVectorView<f64>, column_index: usize) -> DMatrix<f64> { pub fn proj(&self, v: &DVectorView<f64>, column_index: usize)
-> DMatrix<f64>
{
if self.dim() == 0 { if self.dim() == 0 {
const ELEMENT_DIM: usize = 5; const ELEMENT_DIM: usize = 5;
DMatrix::zeros(ELEMENT_DIM, self.assembly_dim) DMatrix::zeros(ELEMENT_DIM, self.assembly_dim)
@ -253,7 +260,7 @@ impl ConstraintProblem {
guess: DMatrix::<f64>::zeros(ELEMENT_DIM, element_count), guess: DMatrix::<f64>::zeros(ELEMENT_DIM, element_count),
} }
} }
#[cfg(feature = "dev")] #[cfg(feature = "dev")]
pub fn from_guess(guess_columns: &[DVector<f64>]) -> Self { pub fn from_guess(guess_columns: &[DVector<f64>]) -> Self {
Self { Self {
@ -291,7 +298,9 @@ impl SearchState {
} }
} }
fn basis_matrix(index: (usize, usize), nrows: usize, ncols: usize) -> DMatrix<f64> { fn basis_matrix(index: (usize, usize), nrows: usize, ncols: usize)
-> DMatrix<f64>
{
let mut result = DMatrix::<f64>::zeros(nrows, ncols); let mut result = DMatrix::<f64>::zeros(nrows, ncols);
result[index] = 1.0; result[index] = 1.0;
result result
@ -377,10 +386,10 @@ pub fn realize_gram(
) -> Realization { ) -> Realization {
// destructure the problem data // destructure the problem data
let ConstraintProblem { gram, guess, frozen } = problem; let ConstraintProblem { gram, guess, frozen } = problem;
// start the descent history // start the descent history
let mut history = DescentHistory::new(); let mut history = DescentHistory::new();
// handle the case where the assembly is empty. our general realization // handle the case where the assembly is empty. our general realization
// routine can't handle this case because it builds the Hessian using // routine can't handle this case because it builds the Hessian using
// `DMatrix::from_columns`, which panics when the list of columns is empty // `DMatrix::from_columns`, which panics when the list of columns is empty
@ -394,29 +403,30 @@ pub fn realize_gram(
); );
return Realization { result, history }; return Realization { result, history };
} }
// find the dimension of the search space // find the dimension of the search space
let element_dim = guess.nrows(); let element_dim = guess.nrows();
let total_dim = element_dim * assembly_dim; let total_dim = element_dim * assembly_dim;
// scale the tolerance // scale the tolerance
let scale_adjustment = (gram.0.len() as f64).sqrt(); let scale_adjustment = (gram.0.len() as f64).sqrt();
let tol = scale_adjustment * scaled_tol; let tol = scale_adjustment * scaled_tol;
// convert the frozen indices to stacked format // convert the frozen indices to stacked format
let frozen_stacked: Vec<usize> = frozen.into_iter().map( let frozen_stacked: Vec<usize> = frozen.into_iter().map(
|MatrixEntry { index: (row, col), .. }| col*element_dim + row |MatrixEntry { index: (row, col), .. }| col*element_dim + row
).collect(); ).collect();
// use a regularized Newton's method with backtracking // use a regularized Newton's method with backtracking
let mut state = SearchState::from_config(gram, frozen.freeze(guess)); let mut state = SearchState::from_config(gram, frozen.freeze(guess));
let mut hess = DMatrix::zeros(element_dim, assembly_dim); let mut hess = DMatrix::zeros(element_dim, assembly_dim);
for _ in 0..max_descent_steps { for _ in 0..max_descent_steps {
// find the negative gradient of the loss function // find the negative gradient of the loss function
let neg_grad = 4.0 * &*Q * &state.config * &state.err_proj; 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>); let mut neg_grad_stacked =
neg_grad.clone().reshape_generic(Dyn(total_dim), Const::<1>);
history.neg_grad.push(neg_grad.clone()); history.neg_grad.push(neg_grad.clone());
// find the negative Hessian of the loss function // find the negative Hessian of the loss function
let mut hess_cols = Vec::<DVector<f64>>::with_capacity(total_dim); let mut hess_cols = Vec::<DVector<f64>>::with_capacity(total_dim);
for col in 0..assembly_dim { for col in 0..assembly_dim {
@ -431,19 +441,21 @@ pub fn realize_gram(
-&basis_mat * &state.err_proj -&basis_mat * &state.err_proj
+ &state.config * &neg_d_err_proj + &state.config * &neg_d_err_proj
); );
hess_cols.push(deriv_grad.reshape_generic(Dyn(total_dim), Const::<1>)); hess_cols.push(
deriv_grad.reshape_generic(Dyn(total_dim), Const::<1>));
} }
} }
hess = DMatrix::from_columns(hess_cols.as_slice()); hess = DMatrix::from_columns(hess_cols.as_slice());
// regularize the Hessian // regularize the Hessian
let hess_eigvals = hess.symmetric_eigenvalues(); let hess_eigvals = hess.symmetric_eigenvalues();
let min_eigval = hess_eigvals.min(); let min_eigval = hess_eigvals.min();
if min_eigval <= 0.0 { if min_eigval <= 0.0 {
hess -= reg_scale * min_eigval * DMatrix::identity(total_dim, total_dim); hess -= reg_scale * min_eigval
* DMatrix::identity(total_dim, total_dim);
} }
history.hess_eigvals.push(hess_eigvals); history.hess_eigvals.push(hess_eigvals);
// project the negative gradient and negative Hessian onto the // project the negative gradient and negative Hessian onto the
// orthogonal complement of the frozen subspace // orthogonal complement of the frozen subspace
let zero_col = DVector::zeros(total_dim); let zero_col = DVector::zeros(total_dim);
@ -454,12 +466,12 @@ pub fn realize_gram(
hess.set_column(k, &zero_col); hess.set_column(k, &zero_col);
hess[(k, k)] = 1.0; hess[(k, k)] = 1.0;
} }
// stop if the loss is tolerably low // stop if the loss is tolerably low
history.config.push(state.config.clone()); history.config.push(state.config.clone());
history.scaled_loss.push(state.loss / scale_adjustment); history.scaled_loss.push(state.loss / scale_adjustment);
if state.loss < tol { break; } if state.loss < tol { break; }
// compute the Newton step // compute the Newton step
/* TO DO */ /* TO DO */
/* /*
@ -477,9 +489,10 @@ pub fn realize_gram(
}, },
}; };
let base_step_stacked = hess_cholesky.solve(&neg_grad_stacked); let base_step_stacked = hess_cholesky.solve(&neg_grad_stacked);
let base_step = base_step_stacked.reshape_generic(Dyn(element_dim), Dyn(assembly_dim)); let base_step = base_step_stacked.reshape_generic(
Dyn(element_dim), Dyn(assembly_dim));
history.base_step.push(base_step.clone()); history.base_step.push(base_step.clone());
// use backtracking line search to find a better configuration // use backtracking line search to find a better configuration
if let Some((better_state, backoff_steps)) = seek_better_config( if let Some((better_state, backoff_steps)) = seek_better_config(
gram, &state, &base_step, neg_grad.dot(&base_step), gram, &state, &base_step, neg_grad.dot(&base_step),
@ -505,11 +518,14 @@ pub fn realize_gram(
.view_mut(block_start, (element_dim, UNIFORM_DIM)) .view_mut(block_start, (element_dim, UNIFORM_DIM))
.copy_from(&local_unif_to_std(state.config.column(n))); .copy_from(&local_unif_to_std(state.config.column(n)));
} }
// find the kernel of the Hessian. give it the uniform inner product // find the kernel of the Hessian. give it the uniform inner product
let tangent = ConfigSubspace::symmetric_kernel(hess, unif_to_std, assembly_dim); let tangent =
ConfigSubspace::symmetric_kernel(hess, unif_to_std, assembly_dim);
Ok(ConfigNeighborhood { #[cfg(feature = "dev")] config: state.config, nbhd: tangent })
Ok(ConfigNeighborhood {
#[cfg(feature = "dev")] config: state.config, nbhd: tangent
})
} else { } else {
Err("Failed to reach target accuracy".to_string()) Err("Failed to reach target accuracy".to_string())
}; };
@ -521,9 +537,9 @@ pub fn realize_gram(
#[cfg(feature = "dev")] #[cfg(feature = "dev")]
pub mod examples { pub mod examples {
use std::f64::consts::PI; use std::f64::consts::PI;
use super::*; use super::*;
// this problem is from a sangaku by Irisawa Shintarō Hiroatsu. the article // this problem is from a sangaku by Irisawa Shintarō Hiroatsu. the article
// below includes a nice translation of the problem statement, which was // below includes a nice translation of the problem statement, which was
// recorded in Uchida Itsumi's book _Kokon sankan_ (_Mathematics, Past and // recorded in Uchida Itsumi's book _Kokon sankan_ (_Mathematics, Past and
@ -547,40 +563,40 @@ pub mod examples {
) )
).collect::<Vec<_>>().as_slice() ).collect::<Vec<_>>().as_slice()
); );
for s in 0..9 { for s in 0..9 {
// each sphere is represented by a spacelike vector // each sphere is represented by a spacelike vector
problem.gram.push_sym(s, s, 1.0); problem.gram.push_sym(s, s, 1.0);
// the circumscribing sphere is tangent to all of the other // the circumscribing sphere is tangent to all of the other
// spheres, with matching orientation // spheres, with matching orientation
if s > 0 { if s > 0 {
problem.gram.push_sym(0, s, 1.0); problem.gram.push_sym(0, s, 1.0);
} }
if s > 2 { if s > 2 {
// each chain sphere is tangent to the "sun" and "moon" // each chain sphere is tangent to the "sun" and "moon"
// spheres, with opposing orientation // spheres, with opposing orientation
for n in 1..3 { for n in 1..3 {
problem.gram.push_sym(s, n, -1.0); problem.gram.push_sym(s, n, -1.0);
} }
// each chain sphere is tangent to the next chain sphere, // each chain sphere is tangent to the next chain sphere,
// with opposing orientation // with opposing orientation
let s_next = 3 + (s-2) % 6; let s_next = 3 + (s-2) % 6;
problem.gram.push_sym(s, s_next, -1.0); problem.gram.push_sym(s, s_next, -1.0);
} }
} }
// the frozen entries fix the radii of the circumscribing sphere, the // the frozen entries fix the radii of the circumscribing sphere, the
// "sun" and "moon" spheres, and one of the chain spheres // "sun" and "moon" spheres, and one of the chain spheres
for k in 0..4 { for k in 0..4 {
problem.frozen.push(3, k, problem.guess[(3, k)]); problem.frozen.push(3, k, problem.guess[(3, k)]);
} }
realize_gram(&problem, scaled_tol, 0.5, 0.9, 1.1, 200, 110) realize_gram(&problem, scaled_tol, 0.5, 0.9, 1.1, 200, 110)
} }
// set up a kaleidocycle, made of points with fixed distances between them, // set up a kaleidocycle, made of points with fixed distances between them,
// and find its tangent space // and find its tangent space
pub fn realize_kaleidocycle(scaled_tol: f64) -> Realization { pub fn realize_kaleidocycle(scaled_tol: f64) -> Realization {
@ -601,27 +617,28 @@ pub mod examples {
} }
).collect::<Vec<_>>().as_slice() ).collect::<Vec<_>>().as_slice()
); );
const N_POINTS: usize = 2 * N_HINGES; const N_POINTS: usize = 2 * N_HINGES;
for block in (0..N_POINTS).step_by(2) { for block in (0..N_POINTS).step_by(2) {
let block_next = (block + 2) % N_POINTS; let block_next = (block + 2) % N_POINTS;
for j in 0..2 { for j in 0..2 {
// diagonal and hinge edges // diagonal and hinge edges
for k in j..2 { for k in j..2 {
problem.gram.push_sym(block + j, block + k, if j == k { 0.0 } else { -0.5 }); problem.gram.push_sym(
block + j, block + k, if j == k { 0.0 } else { -0.5 });
} }
// non-hinge edges // non-hinge edges
for k in 0..2 { for k in 0..2 {
problem.gram.push_sym(block + j, block_next + k, -0.625); problem.gram.push_sym(block + j, block_next + k, -0.625);
} }
} }
} }
for k in 0..N_POINTS { for k in 0..N_POINTS {
problem.frozen.push(3, k, problem.guess[(3, k)]) problem.frozen.push(3, k, problem.guess[(3, k)])
} }
realize_gram(&problem, scaled_tol, 0.5, 0.9, 1.1, 200, 110) realize_gram(&problem, scaled_tol, 0.5, 0.9, 1.1, 200, 110)
} }
} }
@ -630,9 +647,9 @@ pub mod examples {
mod tests { mod tests {
use nalgebra::Vector3; use nalgebra::Vector3;
use std::{f64::consts::{FRAC_1_SQRT_2, PI}, iter}; use std::{f64::consts::{FRAC_1_SQRT_2, PI}, iter};
use super::{*, examples::*}; use super::{*, examples::*};
#[test] #[test]
fn freeze_test() { fn freeze_test() {
let frozen = PartialMatrix(vec![ let frozen = PartialMatrix(vec![
@ -651,7 +668,7 @@ mod tests {
]); ]);
assert_eq!(frozen.freeze(&config), expected_result); assert_eq!(frozen.freeze(&config), expected_result);
} }
#[test] #[test]
fn sub_proj_test() { fn sub_proj_test() {
let target = PartialMatrix(vec![ let target = PartialMatrix(vec![
@ -670,7 +687,7 @@ mod tests {
]); ]);
assert_eq!(target.sub_proj(&attempt), expected_result); assert_eq!(target.sub_proj(&attempt), expected_result);
} }
#[test] #[test]
fn zero_loss_test() { fn zero_loss_test() {
let mut gram = PartialMatrix::new(); let mut gram = PartialMatrix::new();
@ -690,7 +707,7 @@ mod tests {
let state = SearchState::from_config(&gram, config); let state = SearchState::from_config(&gram, config);
assert!(state.loss.abs() < f64::EPSILON); assert!(state.loss.abs() < f64::EPSILON);
} }
/* TO DO */ /* TO DO */
// at the frozen indices, the optimization steps should have exact zeros, // at the frozen indices, the optimization steps should have exact zeros,
// and the realized configuration should have the desired values // and the realized configuration should have the desired values
@ -702,7 +719,8 @@ mod tests {
]); ]);
for j in 0..2 { for j in 0..2 {
for k in j..2 { for k in j..2 {
problem.gram.push_sym(j, k, if (j, k) == (1, 1) { 1.0 } else { 0.0 }); problem.gram.push_sym(
j, k, if (j, k) == (1, 1) { 1.0 } else { 0.0 });
} }
} }
problem.frozen.push(3, 0, problem.guess[(3, 0)]); problem.frozen.push(3, 0, problem.guess[(3, 0)]);
@ -720,21 +738,22 @@ mod tests {
assert_eq!(config[index], value); assert_eq!(config[index], value);
} }
} }
#[test] #[test]
fn irisawa_hexlet_test() { fn irisawa_hexlet_test() {
// solve Irisawa's problem // solve Irisawa's problem
const SCALED_TOL: f64 = 1.0e-12; const SCALED_TOL: f64 = 1.0e-12;
let config = realize_irisawa_hexlet(SCALED_TOL).result.unwrap().config; let config = realize_irisawa_hexlet(SCALED_TOL).result.unwrap().config;
// check against Irisawa's solution // check against Irisawa's solution
let entry_tol = SCALED_TOL.sqrt(); let entry_tol = SCALED_TOL.sqrt();
let solution_diams = [30.0, 10.0, 6.0, 5.0, 15.0, 10.0, 3.75, 2.5, 2.0 + 8.0/11.0]; let solution_diams =
[30.0, 10.0, 6.0, 5.0, 15.0, 10.0, 3.75, 2.5, 2.0 + 8.0/11.0];
for (k, diam) in solution_diams.into_iter().enumerate() { for (k, diam) in solution_diams.into_iter().enumerate() {
assert!((config[(3, k)] - 1.0 / diam).abs() < entry_tol); assert!((config[(3, k)] - 1.0 / diam).abs() < entry_tol);
} }
} }
#[test] #[test]
fn tangent_test_three_spheres() { fn tangent_test_three_spheres() {
const SCALED_TOL: f64 = 1.0e-12; const SCALED_TOL: f64 = 1.0e-12;
@ -758,7 +777,7 @@ mod tests {
let ConfigNeighborhood { config, nbhd: tangent } = result.unwrap(); let ConfigNeighborhood { config, nbhd: tangent } = result.unwrap();
assert_eq!(config, problem.guess); assert_eq!(config, problem.guess);
assert_eq!(history.scaled_loss.len(), 1); assert_eq!(history.scaled_loss.len(), 1);
// list some motions that should form a basis for the tangent space of // list some motions that should form a basis for the tangent space of
// the solution variety // the solution variety
const UNIFORM_DIM: usize = 4; const UNIFORM_DIM: usize = 4;
@ -786,30 +805,37 @@ mod tests {
0.0, 0.0, -1.0, 0.25, 1.0, 0.0, 0.0, -1.0, 0.25, 1.0,
]), ]),
]; ];
// confirm that the dimension of the tangent space is no greater than // confirm that the dimension of the tangent space is no greater than
// expected // expected
assert_eq!(tangent.basis_std.len(), tangent_motions_std.len()); assert_eq!(tangent.basis_std.len(), tangent_motions_std.len());
// confirm that the tangent space contains all the motions we expect it // confirm that the tangent space contains all the motions we expect it
// to. since we've already bounded the dimension of the tangent space, // 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 // this confirms that the tangent space is what we expect it to be
let tol_sq = ((element_dim * assembly_dim) as f64) * SCALED_TOL * SCALED_TOL; let tol_sq = ((element_dim * assembly_dim) as f64)
for (motion_unif, motion_std) in tangent_motions_unif.into_iter().zip(tangent_motions_std) { * SCALED_TOL * SCALED_TOL;
let motion_proj: DMatrix<_> = motion_unif.column_iter().enumerate().map( for (motion_unif, motion_std)
|(k, v)| tangent.proj(&v, k) in tangent_motions_unif.into_iter().zip(tangent_motions_std) {
).sum(); let motion_proj: DMatrix<_> =
motion_unif.column_iter().enumerate().map(
|(k, v)| tangent.proj(&v, k)
).sum();
assert!((motion_std - motion_proj).norm_squared() < tol_sq); assert!((motion_std - motion_proj).norm_squared() < tol_sq);
} }
} }
fn translation_motion_unif(vel: &Vector3<f64>, assembly_dim: usize) -> Vec<DVector<f64>> { fn translation_motion_unif(vel: &Vector3<f64>, assembly_dim: usize)
-> Vec<DVector<f64>>
{
let mut elt_motion = DVector::zeros(4); let mut elt_motion = DVector::zeros(4);
elt_motion.fixed_rows_mut::<3>(0).copy_from(vel); elt_motion.fixed_rows_mut::<3>(0).copy_from(vel);
iter::repeat(elt_motion).take(assembly_dim).collect() iter::repeat(elt_motion).take(assembly_dim).collect()
} }
fn rotation_motion_unif(ang_vel: &Vector3<f64>, points: Vec<DVectorView<f64>>) -> Vec<DVector<f64>> { fn rotation_motion_unif(
ang_vel: &Vector3<f64>, points: Vec<DVectorView<f64>>
) -> Vec<DVector<f64>> {
points.into_iter().map( points.into_iter().map(
|pt| { |pt| {
let vel = ang_vel.cross(&pt.fixed_rows::<3>(0)); let vel = ang_vel.cross(&pt.fixed_rows::<3>(0));
@ -819,7 +845,7 @@ mod tests {
} }
).collect() ).collect()
} }
#[test] #[test]
fn tangent_test_kaleidocycle() { fn tangent_test_kaleidocycle() {
// set up a kaleidocycle and find its tangent space // set up a kaleidocycle and find its tangent space
@ -827,7 +853,7 @@ mod tests {
let Realization { result, history } = realize_kaleidocycle(SCALED_TOL); let Realization { result, history } = realize_kaleidocycle(SCALED_TOL);
let ConfigNeighborhood { config, nbhd: tangent } = result.unwrap(); let ConfigNeighborhood { config, nbhd: tangent } = result.unwrap();
assert_eq!(history.scaled_loss.len(), 1); assert_eq!(history.scaled_loss.len(), 1);
// list some motions that should form a basis for the tangent space of // list some motions that should form a basis for the tangent space of
// the solution variety // the solution variety
const N_HINGES: usize = 6; const N_HINGES: usize = 6;
@ -838,12 +864,15 @@ mod tests {
translation_motion_unif(&Vector3::new(1.0, 0.0, 0.0), assembly_dim), translation_motion_unif(&Vector3::new(1.0, 0.0, 0.0), assembly_dim),
translation_motion_unif(&Vector3::new(0.0, 1.0, 0.0), assembly_dim), translation_motion_unif(&Vector3::new(0.0, 1.0, 0.0), assembly_dim),
translation_motion_unif(&Vector3::new(0.0, 0.0, 1.0), assembly_dim), translation_motion_unif(&Vector3::new(0.0, 0.0, 1.0), assembly_dim),
// the rotations about the coordinate axes // the rotations about the coordinate axes
rotation_motion_unif(&Vector3::new(1.0, 0.0, 0.0), config.column_iter().collect()), rotation_motion_unif(
rotation_motion_unif(&Vector3::new(0.0, 1.0, 0.0), config.column_iter().collect()), &Vector3::new(1.0, 0.0, 0.0), config.column_iter().collect()),
rotation_motion_unif(&Vector3::new(0.0, 0.0, 1.0), config.column_iter().collect()), rotation_motion_unif(
&Vector3::new(0.0, 1.0, 0.0), config.column_iter().collect()),
rotation_motion_unif(
&Vector3::new(0.0, 0.0, 1.0), config.column_iter().collect()),
// the twist motion. more precisely: a motion that keeps the center // the twist motion. more precisely: a motion that keeps the center
// of mass stationary and preserves the distances between the // of mass stationary and preserves the distances between the
// vertices to first order. this has to be the twist as long as: // vertices to first order. this has to be the twist as long as:
@ -859,8 +888,10 @@ mod tests {
[ [
DVector::from_column_slice(&[0.0, 0.0, 5.0, 0.0]), DVector::from_column_slice(&[0.0, 0.0, 5.0, 0.0]),
DVector::from_column_slice(&[0.0, 0.0, 1.0, 0.0]), DVector::from_column_slice(&[0.0, 0.0, 1.0, 0.0]),
DVector::from_column_slice(&[-vel_vert_x, -vel_vert_y, -3.0, 0.0]), DVector::from_column_slice(
DVector::from_column_slice(&[vel_vert_x, vel_vert_y, -3.0, 0.0]), &[-vel_vert_x, -vel_vert_y, -3.0, 0.0]),
DVector::from_column_slice(
&[vel_vert_x, vel_vert_y, -3.0, 0.0]),
] ]
} }
).collect::<Vec<_>>(), ).collect::<Vec<_>>(),
@ -872,23 +903,26 @@ mod tests {
).collect::<Vec<_>>() ).collect::<Vec<_>>()
) )
).collect::<Vec<_>>(); ).collect::<Vec<_>>();
// confirm that the dimension of the tangent space is no greater than // confirm that the dimension of the tangent space is no greater than
// expected // expected
assert_eq!(tangent.basis_std.len(), tangent_motions_unif.len()); assert_eq!(tangent.basis_std.len(), tangent_motions_unif.len());
// confirm that the tangent space contains all the motions we expect it // confirm that the tangent space contains all the motions we expect it
// to. since we've already bounded the dimension of the tangent space, // 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 // this confirms that the tangent space is what we expect it to be
let tol_sq = ((element_dim * assembly_dim) as f64) * SCALED_TOL * SCALED_TOL; let tol_sq = ((element_dim * assembly_dim) as f64)
for (motion_unif, motion_std) in tangent_motions_unif.into_iter().zip(tangent_motions_std) { * SCALED_TOL * SCALED_TOL;
let motion_proj: DMatrix<_> = motion_unif.into_iter().enumerate().map( for (motion_unif, motion_std)
|(k, v)| tangent.proj(&v.as_view(), k) in tangent_motions_unif.into_iter().zip(tangent_motions_std) {
).sum(); let motion_proj: DMatrix<_> =
motion_unif.into_iter().enumerate().map(
|(k, v)| tangent.proj(&v.as_view(), k)
).sum();
assert!((motion_std - motion_proj).norm_squared() < tol_sq); assert!((motion_std - motion_proj).norm_squared() < tol_sq);
} }
} }
fn translation(dis: Vector3<f64>) -> DMatrix<f64> { fn translation(dis: Vector3<f64>) -> DMatrix<f64> {
const ELEMENT_DIM: usize = 5; const ELEMENT_DIM: usize = 5;
DMatrix::from_column_slice(ELEMENT_DIM, ELEMENT_DIM, &[ DMatrix::from_column_slice(ELEMENT_DIM, ELEMENT_DIM, &[
@ -899,7 +933,7 @@ mod tests {
0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 0.0, 1.0,
]) ])
} }
// confirm that projection onto a configuration subspace is equivariant with // confirm that projection onto a configuration subspace is equivariant with
// respect to Euclidean motions // respect to Euclidean motions
#[test] #[test]
@ -913,13 +947,13 @@ mod tests {
problem_orig.gram.push_sym(0, 0, 1.0); problem_orig.gram.push_sym(0, 0, 1.0);
problem_orig.gram.push_sym(1, 1, 1.0); problem_orig.gram.push_sym(1, 1, 1.0);
problem_orig.gram.push_sym(0, 1, 0.5); problem_orig.gram.push_sym(0, 1, 0.5);
let Realization { result: result_orig, history: history_orig } = realize_gram( let Realization { result: result_orig, history: history_orig } =
&problem_orig, SCALED_TOL, 0.5, 0.9, 1.1, 200, 110 realize_gram(&problem_orig, SCALED_TOL, 0.5, 0.9, 1.1, 200, 110);
); let ConfigNeighborhood { config: config_orig, nbhd: tangent_orig } =
let ConfigNeighborhood { config: config_orig, nbhd: tangent_orig } = result_orig.unwrap(); result_orig.unwrap();
assert_eq!(config_orig, problem_orig.guess); assert_eq!(config_orig, problem_orig.guess);
assert_eq!(history_orig.scaled_loss.len(), 1); assert_eq!(history_orig.scaled_loss.len(), 1);
// find another pair of spheres that meet at 120°. we'll think of this // find another pair of spheres that meet at 120°. we'll think of this
// solution as a transformed version of the original one // solution as a transformed version of the original one
let guess_tfm = { let guess_tfm = {
@ -934,23 +968,24 @@ mod tests {
frozen: problem_orig.frozen, frozen: problem_orig.frozen,
guess: guess_tfm, guess: guess_tfm,
}; };
let Realization { result: result_tfm, history: history_tfm } = realize_gram( let Realization { result: result_tfm, history: history_tfm } =
&problem_tfm, SCALED_TOL, 0.5, 0.9, 1.1, 200, 110 realize_gram(&problem_tfm, SCALED_TOL, 0.5, 0.9, 1.1, 200, 110);
); let ConfigNeighborhood { config: config_tfm, nbhd: tangent_tfm } =
let ConfigNeighborhood { config: config_tfm, nbhd: tangent_tfm } = result_tfm.unwrap(); result_tfm.unwrap();
assert_eq!(config_tfm, problem_tfm.guess); assert_eq!(config_tfm, problem_tfm.guess);
assert_eq!(history_tfm.scaled_loss.len(), 1); assert_eq!(history_tfm.scaled_loss.len(), 1);
// project a nudge to the tangent space of the solution variety at the // project a nudge to the tangent space of the solution variety at the
// original solution // original solution
let motion_orig = DVector::from_column_slice(&[0.0, 0.0, 1.0, 0.0]); let motion_orig = DVector::from_column_slice(&[0.0, 0.0, 1.0, 0.0]);
let motion_orig_proj = tangent_orig.proj(&motion_orig.as_view(), 0); let motion_orig_proj = tangent_orig.proj(&motion_orig.as_view(), 0);
// project the equivalent nudge to the tangent space of the solution // project the equivalent nudge to the tangent space of the solution
// variety at the transformed solution // variety at the transformed solution
let motion_tfm = DVector::from_column_slice(&[FRAC_1_SQRT_2, 0.0, FRAC_1_SQRT_2, 0.0]); let motion_tfm = DVector::from_column_slice(
&[FRAC_1_SQRT_2, 0.0, FRAC_1_SQRT_2, 0.0]);
let motion_tfm_proj = tangent_tfm.proj(&motion_tfm.as_view(), 0); let motion_tfm_proj = tangent_tfm.proj(&motion_tfm.as_view(), 0);
// take the transformation that sends the original solution to the // take the transformation that sends the original solution to the
// transformed solution and apply it to the motion that the original // transformed solution and apply it to the motion that the original
// solution makes in response to the nudge // solution makes in response to the nudge
@ -964,12 +999,14 @@ mod tests {
]); ]);
let transl = translation(Vector3::new(0.0, 0.0, 7.0)); let transl = translation(Vector3::new(0.0, 0.0, 7.0));
let motion_proj_tfm = transl * rot * motion_orig_proj; let motion_proj_tfm = transl * rot * motion_orig_proj;
// confirm that the projection of the nudge is equivariant. we loosen // confirm that the projection of the nudge is equivariant. we loosen
// the comparison tolerance because the transformation seems to // the comparison tolerance because the transformation seems to
// introduce some numerical error // introduce some numerical error
const SCALED_TOL_TFM: f64 = 1.0e-9; const SCALED_TOL_TFM: f64 = 1.0e-9;
let tol_sq = ((problem_orig.guess.nrows() * problem_orig.guess.ncols()) as f64) * SCALED_TOL_TFM * SCALED_TOL_TFM; let tol_sq = ((problem_orig.guess.nrows()
* problem_orig.guess.ncols()) as f64)
* SCALED_TOL_TFM * SCALED_TOL_TFM;
assert!((motion_proj_tfm - motion_tfm_proj).norm_squared() < tol_sq); assert!((motion_proj_tfm - motion_tfm_proj).norm_squared() < tol_sq);
} }
} }

8
app-proto/src/main.rs
View file

@ -30,7 +30,7 @@ impl AppState {
selection: create_signal(BTreeSet::default()), selection: create_signal(BTreeSet::default()),
} }
} }
// in single-selection mode, select the given element. in multiple-selection // in single-selection mode, select the given element. in multiple-selection
// mode, toggle whether the given element is selected // mode, toggle whether the given element is selected
fn select(&self, element: &Rc<dyn Element>, multi: bool) { fn select(&self, element: &Rc<dyn Element>, multi: bool) {
@ -53,10 +53,10 @@ fn main() {
// set the console error panic hook // set the console error panic hook
#[cfg(feature = "console_error_panic_hook")] #[cfg(feature = "console_error_panic_hook")]
console_error_panic_hook::set_once(); console_error_panic_hook::set_once();
sycamore::render(|| { sycamore::render(|| {
provide_context(AppState::new()); provide_context(AppState::new());
view! { view! {
div(id = "sidebar") { div(id = "sidebar") {
AddRemove {} AddRemove {}
@ -66,4 +66,4 @@ fn main() {
Display {} Display {}
} }
}); });
} }

6
app-proto/src/specified.rs
View file

@ -20,7 +20,7 @@ impl SpecifiedValue {
pub fn from_empty_spec() -> Self { pub fn from_empty_spec() -> Self {
Self { spec: String::new(), value: None } Self { spec: String::new(), value: None }
} }
pub fn is_present(&self) -> bool { pub fn is_present(&self) -> bool {
matches!(self.value, Some(_)) matches!(self.value, Some(_))
} }
@ -42,7 +42,7 @@ impl From<Option<f64>> for SpecifiedValue {
// if the specification is properly formatted, and `Error` if not // if the specification is properly formatted, and `Error` if not
impl TryFrom<String> for SpecifiedValue { impl TryFrom<String> for SpecifiedValue {
type Error = ParseFloatError; type Error = ParseFloatError;
fn try_from(spec: String) -> Result<Self, Self::Error> { fn try_from(spec: String) -> Result<Self, Self::Error> {
if spec.is_empty() { if spec.is_empty() {
Ok(Self::from_empty_spec()) Ok(Self::from_empty_spec())
@ -52,4 +52,4 @@ impl TryFrom<String> for SpecifiedValue {
) )
} }
} }
} }