diff --git a/app-proto/examples/common/mod.rs b/app-proto/examples/common/mod.rs new file mode 100644 index 0000000..b1b91b5 --- /dev/null +++ b/app-proto/examples/common/mod.rs @@ -0,0 +1,36 @@ +#![allow(dead_code)] + +use nalgebra::DMatrix; + +use dyna3::engine::{Q, DescentHistory, RealizationResult}; + +pub fn print_title(title: &str) { + println!("─── {title} ───"); +} + +pub fn print_realization_diagnostics(realization_result: &RealizationResult) { + let RealizationResult { result, history } = realization_result; + println!(); + if let Err(ref msg) = result { + println!("❌️ {msg}"); + } else { + println!("✅️ Target accuracy achieved!"); + } + println!("Steps: {}", history.scaled_loss.len() - 1); + println!("Loss: {}", history.scaled_loss.last().unwrap()); +} + +pub fn print_gram_matrix(config: &DMatrix) { + println!("\nCompleted Gram matrix:{}", (config.tr_mul(&*Q) * config).to_string().trim_end()); +} + +pub fn print_config(config: &DMatrix) { + println!("\nConfiguration:{}", config.to_string().trim_end()); +} + +pub fn print_loss_history(history: &DescentHistory) { + println!("\nStep │ Loss\n─────┼────────────────────────────────"); + for (step, scaled_loss) in history.scaled_loss.iter().enumerate() { + println!("{:<4} │ {}", step, scaled_loss); + } +} \ No newline at end of file diff --git a/app-proto/examples/irisawa-hexlet.rs b/app-proto/examples/irisawa-hexlet.rs index 639a494..750a0d0 100644 --- a/app-proto/examples/irisawa-hexlet.rs +++ b/app-proto/examples/irisawa-hexlet.rs @@ -1,25 +1,28 @@ -use dyna3::engine::{Q, examples::realize_irisawa_hexlet}; +mod common; + +use common::{ + print_gram_matrix, + print_loss_history, + print_realization_diagnostics, + print_title +}; +use dyna3::engine::{Realization, examples::realize_irisawa_hexlet}; fn main() { const SCALED_TOL: f64 = 1.0e-12; - let (config, _, success, history) = realize_irisawa_hexlet(SCALED_TOL); - print!("\nCompleted Gram matrix:{}", config.tr_mul(&*Q) * &config); - if success { - println!("Target accuracy achieved!"); - } else { - println!("Failed to reach target accuracy"); - } - println!("Steps: {}", history.scaled_loss.len() - 1); - println!("Loss: {}", history.scaled_loss.last().unwrap()); - if success { + let realization_result = realize_irisawa_hexlet(SCALED_TOL); + print_title("Irisawa hexlet"); + print_realization_diagnostics(&realization_result); + if let Ok(Realization { config, .. }) = realization_result.result { + // print the diameters of the chain spheres println!("\nChain diameters:"); println!(" {} sun (given)", 1.0 / config[(3, 3)]); for k in 4..9 { println!(" {} sun", 1.0 / config[(3, k)]); } + + // print the completed Gram matrix + print_gram_matrix(&config); } - println!("\nStep │ Loss\n─────┼────────────────────────────────"); - for (step, scaled_loss) in history.scaled_loss.into_iter().enumerate() { - println!("{:<4} │ {}", step, scaled_loss); - } + print_loss_history(&realization_result.history); } \ No newline at end of file diff --git a/app-proto/examples/kaleidocycle.rs b/app-proto/examples/kaleidocycle.rs index 2779ab1..21f8187 100644 --- a/app-proto/examples/kaleidocycle.rs +++ b/app-proto/examples/kaleidocycle.rs @@ -1,30 +1,37 @@ +mod common; + use nalgebra::{DMatrix, DVector}; -use dyna3::engine::{Q, examples::realize_kaleidocycle}; +use common::{ + print_config, + print_gram_matrix, + print_realization_diagnostics, + print_title +}; +use dyna3::engine::{Realization, examples::realize_kaleidocycle}; fn main() { const SCALED_TOL: f64 = 1.0e-12; - let (config, tangent, success, history) = realize_kaleidocycle(SCALED_TOL); - print!("Completed Gram matrix:{}", config.tr_mul(&*Q) * &config); - print!("Configuration:{}", config); - if success { - println!("Target accuracy achieved!"); - } else { - println!("Failed to reach target accuracy"); + let realization_result = realize_kaleidocycle(SCALED_TOL); + print_title("Kaleidocycle"); + print_realization_diagnostics(&realization_result); + if let Ok(Realization { config, tangent }) = realization_result.result { + // print the completed Gram matrix and the realized configuration + print_gram_matrix(&config); + print_config(&config); + + // find the kaleidocycle's twist motion by projecting onto the tangent + // space + const N_POINTS: usize = 12; + let up = DVector::from_column_slice(&[0.0, 0.0, 1.0, 0.0]); + let down = -&up; + let twist_motion: DMatrix<_> = (0..N_POINTS).step_by(4).flat_map( + |n| [ + tangent.proj(&up.as_view(), n), + tangent.proj(&down.as_view(), n+1) + ] + ).sum(); + let normalization = 5.0 / twist_motion[(2, 0)]; + println!("\nTwist motion:{}", (normalization * twist_motion).to_string().trim_end()); } - println!("Steps: {}", history.scaled_loss.len() - 1); - println!("Loss: {}\n", history.scaled_loss.last().unwrap()); - - // find the kaleidocycle's twist motion by projecting onto the tangent space - const N_POINTS: usize = 12; - let up = DVector::from_column_slice(&[0.0, 0.0, 1.0, 0.0]); - let down = -&up; - let twist_motion: DMatrix<_> = (0..N_POINTS).step_by(4).flat_map( - |n| [ - tangent.proj(&up.as_view(), n), - tangent.proj(&down.as_view(), n+1) - ] - ).sum(); - let normalization = 5.0 / twist_motion[(2, 0)]; - print!("Twist motion:{}", normalization * twist_motion); } \ No newline at end of file diff --git a/app-proto/examples/point-on-sphere.rs b/app-proto/examples/point-on-sphere.rs index 880d7b0..774368e 100644 --- a/app-proto/examples/point-on-sphere.rs +++ b/app-proto/examples/point-on-sphere.rs @@ -1,4 +1,19 @@ -use dyna3::engine::{Q, point, realize_gram, sphere, ConstraintProblem}; +mod common; + +use common::{ + print_config, + print_gram_matrix, + print_loss_history, + print_realization_diagnostics, + print_title +}; +use dyna3::engine::{ + point, + realize_gram, + sphere, + ConstraintProblem, + Realization +}; fn main() { let mut problem = ConstraintProblem::from_guess(&[ @@ -11,21 +26,14 @@ fn main() { } } problem.frozen.push(3, 0, problem.guess[(3, 0)]); - println!(); - let (config, _, success, history) = realize_gram( + let realization_result = realize_gram( &problem, 1.0e-12, 0.5, 0.9, 1.1, 200, 110 ); - print!("\nCompleted Gram matrix:{}", config.tr_mul(&*Q) * &config); - print!("Configuration:{}", config); - if success { - println!("Target accuracy achieved!"); - } else { - println!("Failed to reach target accuracy"); - } - println!("Steps: {}", history.scaled_loss.len() - 1); - println!("Loss: {}", history.scaled_loss.last().unwrap()); - println!("\nStep │ Loss\n─────┼────────────────────────────────"); - for (step, scaled_loss) in history.scaled_loss.into_iter().enumerate() { - println!("{:<4} │ {}", step, scaled_loss); + print_title("Point on a sphere"); + print_realization_diagnostics(&realization_result); + if let Ok(Realization{ config, .. }) = realization_result.result { + print_gram_matrix(&config); + print_config(&config); } + print_loss_history(&realization_result.history); } \ No newline at end of file diff --git a/app-proto/examples/three-spheres.rs b/app-proto/examples/three-spheres.rs index 3f3cc44..dd2cdc0 100644 --- a/app-proto/examples/three-spheres.rs +++ b/app-proto/examples/three-spheres.rs @@ -1,4 +1,12 @@ -use dyna3::engine::{Q, realize_gram, sphere, ConstraintProblem}; +mod common; + +use common::{ + print_gram_matrix, + print_loss_history, + print_realization_diagnostics, + print_title +}; +use dyna3::engine::{realize_gram, sphere, ConstraintProblem, Realization}; fn main() { let mut problem = ConstraintProblem::from_guess({ @@ -14,20 +22,13 @@ fn main() { problem.gram.push_sym(j, k, if j == k { 1.0 } else { -1.0 }); } } - println!(); - let (config, _, success, history) = realize_gram( + let realization_result = realize_gram( &problem, 1.0e-12, 0.5, 0.9, 1.1, 200, 110 ); - print!("\nCompleted Gram matrix:{}", config.tr_mul(&*Q) * &config); - if success { - println!("Target accuracy achieved!"); - } else { - println!("Failed to reach target accuracy"); - } - println!("Steps: {}", history.scaled_loss.len() - 1); - println!("Loss: {}", history.scaled_loss.last().unwrap()); - println!("\nStep │ Loss\n─────┼────────────────────────────────"); - for (step, scaled_loss) in history.scaled_loss.into_iter().enumerate() { - println!("{:<4} │ {}", step, scaled_loss); + print_title("Three spheres"); + print_realization_diagnostics(&realization_result); + if let Ok(Realization{ config, .. }) = realization_result.result { + print_gram_matrix(&config); } + print_loss_history(&realization_result.history); } \ No newline at end of file diff --git a/app-proto/run-examples b/app-proto/run-examples index 52173b0..b3e3121 100755 --- a/app-proto/run-examples +++ b/app-proto/run-examples @@ -6,7 +6,7 @@ # http://xion.io/post/code/rust-examples.html # -cargo run --example irisawa-hexlet -cargo run --example three-spheres -cargo run --example point-on-sphere +cargo run --example irisawa-hexlet; echo +cargo run --example three-spheres; echo +cargo run --example point-on-sphere; echo cargo run --example kaleidocycle \ No newline at end of file diff --git a/app-proto/src/assembly.rs b/app-proto/src/assembly.rs index 6c91fc0..9f66055 100644 --- a/app-proto/src/assembly.rs +++ b/app-proto/src/assembly.rs @@ -24,7 +24,9 @@ use crate::{ realize_gram, sphere, ConfigSubspace, - ConstraintProblem + ConstraintProblem, + Realization, + RealizationResult }, outline::OutlineItem, specified::SpecifiedValue @@ -687,22 +689,25 @@ impl Assembly { console_log!("Old configuration:{:>8.3}", problem.guess); // look for a configuration with the given Gram matrix - let (config, tangent, success, history) = realize_gram( + let RealizationResult { result, history } = realize_gram( &problem, 1.0e-12, 0.5, 0.9, 1.1, 200, 110 ); /* DEBUG */ // report the outcome of the search - if success { - console_log!("Target accuracy achieved!") + if let Err(ref msg) = result { + console_log!("❌️ {msg}"); } else { - console_log!("Failed to reach target accuracy") + console_log!("✅️ Target accuracy achieved!"); } console_log!("Steps: {}", history.scaled_loss.len() - 1); - console_log!("Loss: {}", *history.scaled_loss.last().unwrap()); - console_log!("Tangent dimension: {}", tangent.dim()); + console_log!("Loss: {}", history.scaled_loss.last().unwrap()); - if success { + if let Ok(Realization { config, tangent }) = result { + /* DEBUG */ + // report the tangent dimension + console_log!("Tangent dimension: {}", tangent.dim()); + // read out the solution for elt in self.elements.get_clone_untracked() { elt.representation().update( diff --git a/app-proto/src/engine.rs b/app-proto/src/engine.rs index c5d7b00..2524fcd 100644 --- a/app-proto/src/engine.rs +++ b/app-proto/src/engine.rs @@ -393,6 +393,16 @@ fn seek_better_config( None } +pub struct Realization { + pub config: DMatrix, + pub tangent: ConfigSubspace +} + +pub struct RealizationResult { + pub result: Result, + pub history: DescentHistory +} + // seek a matrix `config` that matches the partial matrix `problem.frozen` and // has `config' * Q * config` matching the partial matrix `problem.gram`. start // at `problem.guess`, set the frozen entries to their desired values, and then @@ -405,7 +415,7 @@ pub fn realize_gram( reg_scale: f64, max_descent_steps: i32, max_backoff_steps: i32 -) -> (DMatrix, ConfigSubspace, bool, DescentHistory) { +) -> RealizationResult { // destructure the problem data let ConstraintProblem { gram, guess, frozen @@ -491,19 +501,20 @@ pub fn realize_gram( history.base_step.push(base_step.clone()); // use backtracking line search to find a better configuration - match seek_better_config( + if let Some((better_state, backoff_steps)) = seek_better_config( gram, &state, &base_step, neg_grad.dot(&base_step), min_efficiency, backoff, max_backoff_steps ) { - Some((better_state, backoff_steps)) => { - state = better_state; - history.backoff_steps.push(backoff_steps); - }, - None => return (state.config, ConfigSubspace::zero(assembly_dim), false, history) + state = better_state; + history.backoff_steps.push(backoff_steps); + } else { + return RealizationResult { + result: Err("Line search failed".to_string()), + history + } }; } - let success = state.loss < tol; - let tangent = if success { + let result = if state.loss < tol { // express the uniform basis in the standard basis const UNIFORM_DIM: usize = 4; let total_dim_unif = UNIFORM_DIM * assembly_dim; @@ -516,11 +527,13 @@ pub fn realize_gram( } // find the kernel of the Hessian. give it the uniform inner product - ConfigSubspace::symmetric_kernel(hess, unif_to_std, assembly_dim) + let tangent = ConfigSubspace::symmetric_kernel(hess, unif_to_std, assembly_dim); + + Ok(Realization { config: state.config, tangent }) } else { - ConfigSubspace::zero(assembly_dim) + Err("Failed to reach target accuracy".to_string()) }; - (state.config, tangent, success, history) + RealizationResult{ result, history } } // --- tests --- @@ -539,7 +552,7 @@ pub mod examples { // "Japan's 'Wasan' Mathematical Tradition", by Abe Haruki // https://www.nippon.com/en/japan-topics/c12801/ // - pub fn realize_irisawa_hexlet(scaled_tol: f64) -> (DMatrix, ConfigSubspace, bool, DescentHistory) { + pub fn realize_irisawa_hexlet(scaled_tol: f64) -> RealizationResult { let mut problem = ConstraintProblem::from_guess( [ sphere(0.0, 0.0, 0.0, 15.0), @@ -590,7 +603,7 @@ pub mod examples { // set up a kaleidocycle, made of points with fixed distances between them, // and find its tangent space - pub fn realize_kaleidocycle(scaled_tol: f64) -> (DMatrix, ConfigSubspace, bool, DescentHistory) { + pub fn realize_kaleidocycle(scaled_tol: f64) -> RealizationResult { const N_HINGES: usize = 6; let mut problem = ConstraintProblem::from_guess( (0..N_HINGES).step_by(2).flat_map( @@ -714,10 +727,10 @@ mod tests { } problem.frozen.push(3, 0, problem.guess[(3, 0)]); problem.frozen.push(3, 1, 0.5); - let (config, _, success, history) = realize_gram( + let RealizationResult { result, history } = realize_gram( &problem, 1.0e-12, 0.5, 0.9, 1.1, 200, 110 ); - assert_eq!(success, true); + let config = result.unwrap().config; for base_step in history.base_step.into_iter() { for &MatrixEntry { index, .. } in &problem.frozen { assert_eq!(base_step[index], 0.0); @@ -732,7 +745,7 @@ mod tests { fn irisawa_hexlet_test() { // solve Irisawa's problem const SCALED_TOL: f64 = 1.0e-12; - let (config, _, _, _) = realize_irisawa_hexlet(SCALED_TOL); + let config = realize_irisawa_hexlet(SCALED_TOL).result.unwrap().config; // check against Irisawa's solution let entry_tol = SCALED_TOL.sqrt(); @@ -759,11 +772,11 @@ mod tests { for n in 0..ELEMENT_DIM { problem.frozen.push(n, 0, problem.guess[(n, 0)]); } - let (config, tangent, success, history) = realize_gram( + let RealizationResult { result, history } = realize_gram( &problem, SCALED_TOL, 0.5, 0.9, 1.1, 200, 110 ); + let Realization { config, tangent } = result.unwrap(); assert_eq!(config, problem.guess); - assert_eq!(success, true); assert_eq!(history.scaled_loss.len(), 1); // list some motions that should form a basis for the tangent space of @@ -831,8 +844,8 @@ mod tests { fn tangent_test_kaleidocycle() { // set up a kaleidocycle and find its tangent space const SCALED_TOL: f64 = 1.0e-12; - let (config, tangent, success, history) = realize_kaleidocycle(SCALED_TOL); - assert_eq!(success, true); + let RealizationResult { result, history } = realize_kaleidocycle(SCALED_TOL); + let Realization { config, tangent } = result.unwrap(); assert_eq!(history.scaled_loss.len(), 1); // list some motions that should form a basis for the tangent space of @@ -920,11 +933,11 @@ mod tests { problem_orig.gram.push_sym(0, 0, 1.0); problem_orig.gram.push_sym(1, 1, 1.0); problem_orig.gram.push_sym(0, 1, 0.5); - let (config_orig, tangent_orig, success_orig, history_orig) = realize_gram( + let RealizationResult { result: result_orig, history: history_orig } = realize_gram( &problem_orig, SCALED_TOL, 0.5, 0.9, 1.1, 200, 110 ); + let Realization { config: config_orig, tangent: tangent_orig } = result_orig.unwrap(); assert_eq!(config_orig, problem_orig.guess); - assert_eq!(success_orig, true); assert_eq!(history_orig.scaled_loss.len(), 1); // find another pair of spheres that meet at 120°. we'll think of this @@ -941,11 +954,11 @@ mod tests { guess: guess_tfm, frozen: problem_orig.frozen }; - let (config_tfm, tangent_tfm, success_tfm, history_tfm) = realize_gram( + let RealizationResult { result: result_tfm, history: history_tfm } = realize_gram( &problem_tfm, SCALED_TOL, 0.5, 0.9, 1.1, 200, 110 ); + let Realization { config: config_tfm, tangent: tangent_tfm } = result_tfm.unwrap(); assert_eq!(config_tfm, problem_tfm.guess); - assert_eq!(success_tfm, true); assert_eq!(history_tfm.scaled_loss.len(), 1); // project a nudge to the tangent space of the solution variety at the