diff --git a/app-proto/Cargo.toml b/app-proto/Cargo.toml index e623b26..38205a7 100644 --- a/app-proto/Cargo.toml +++ b/app-proto/Cargo.toml @@ -6,6 +6,7 @@ edition = "2021" [features] default = ["console_error_panic_hook"] +dev = [] [dependencies] itertools = "0.13.0" @@ -36,7 +37,12 @@ features = [ 'WebGlVertexArrayObject' ] +# the self-dependency specifies features to use for tests and examples +# +# https://github.com/rust-lang/cargo/issues/2911#issuecomment-1483256987 +# [dev-dependencies] +dyna3 = { path = ".", default-features = false, features = ["dev"] } wasm-bindgen-test = "0.3.34" [profile.release] diff --git a/app-proto/examples/irisawa-hexlet.rs b/app-proto/examples/irisawa-hexlet.rs new file mode 100644 index 0000000..fc14f91 --- /dev/null +++ b/app-proto/examples/irisawa-hexlet.rs @@ -0,0 +1,25 @@ +use dyna3::engine::{Q, irisawa::realize_irisawa_hexlet}; + +fn main() { + const SCALED_TOL: f64 = 1.0e-12; + let (config, success, history) = realize_irisawa_hexlet(SCALED_TOL); + 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 { + println!("\nChain diameters:"); + println!(" {} sun (given)", 1.0 / config[(3, 3)]); + for k in 4..9 { + println!(" {} sun", 1.0 / config[(3, k)]); + } + } + println!("\nStep │ Loss\n─────┼────────────────────────────────"); + for (step, scaled_loss) in history.scaled_loss.into_iter().enumerate() { + println!("{:<4} │ {}", step, scaled_loss); + } +} \ No newline at end of file diff --git a/app-proto/examples/point-on-sphere.rs b/app-proto/examples/point-on-sphere.rs new file mode 100644 index 0000000..0e4765a --- /dev/null +++ b/app-proto/examples/point-on-sphere.rs @@ -0,0 +1,38 @@ +use nalgebra::DMatrix; + +use dyna3::engine::{Q, point, realize_gram, sphere, PartialMatrix}; + +fn main() { + let gram = { + let mut gram_to_be = PartialMatrix::new(); + for j in 0..2 { + for k in j..2 { + gram_to_be.push_sym(j, k, if (j, k) == (1, 1) { 1.0 } else { 0.0 }); + } + } + gram_to_be + }; + let guess = DMatrix::from_columns(&[ + point(0.0, 0.0, 2.0), + sphere(0.0, 0.0, 0.0, 1.0) + ]); + let frozen = [(3, 0)]; + println!(); + let (config, success, history) = realize_gram( + &gram, guess, &frozen, + 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); + } +} \ No newline at end of file diff --git a/app-proto/examples/three-spheres.rs b/app-proto/examples/three-spheres.rs new file mode 100644 index 0000000..d348b18 --- /dev/null +++ b/app-proto/examples/three-spheres.rs @@ -0,0 +1,40 @@ +use nalgebra::DMatrix; + +use dyna3::engine::{Q, realize_gram, sphere, PartialMatrix}; + +fn main() { + let gram = { + let mut gram_to_be = PartialMatrix::new(); + for j in 0..3 { + for k in j..3 { + gram_to_be.push_sym(j, k, if j == k { 1.0 } else { -1.0 }); + } + } + gram_to_be + }; + let guess = { + let a: f64 = 0.75_f64.sqrt(); + DMatrix::from_columns(&[ + sphere(1.0, 0.0, 0.0, 1.0), + sphere(-0.5, a, 0.0, 1.0), + sphere(-0.5, -a, 0.0, 1.0) + ]) + }; + println!(); + let (config, success, history) = realize_gram( + &gram, guess, &[], + 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); + } +} \ No newline at end of file diff --git a/app-proto/run-examples b/app-proto/run-examples index 6a5e3ae..9791a9d 100755 --- a/app-proto/run-examples +++ b/app-proto/run-examples @@ -1,8 +1,9 @@ -# based on "Enabling print statements in Cargo tests", by Jon Almeida +# run all Cargo examples, as described here: # -# https://jonalmeida.com/posts/2015/01/23/print-cargo/ +# Karol Kuczmarski. "Add examples to your Rust libraries" +# http://xion.io/post/code/rust-examples.html # -cargo test -- --nocapture engine::tests::irisawa_hexlet_test -cargo test -- --nocapture engine::tests::three_spheres_example -cargo test -- --nocapture engine::tests::point_on_sphere_example +cargo run --example irisawa-hexlet +cargo run --example three-spheres +cargo run --example point-on-sphere \ No newline at end of file diff --git a/app-proto/src/engine.rs b/app-proto/src/engine.rs index 343b96e..285a9c6 100644 --- a/app-proto/src/engine.rs +++ b/app-proto/src/engine.rs @@ -4,7 +4,7 @@ use web_sys::{console, wasm_bindgen::JsValue}; /* DEBUG */ // --- elements --- -#[cfg(test)] +#[cfg(feature = "dev")] pub fn point(x: f64, y: f64, z: f64) -> DVector { DVector::from_column_slice(&[x, y, z, 0.5, 0.5*(x*x + y*y + z*z)]) } @@ -113,7 +113,7 @@ impl DescentHistory { // the Lorentz form lazy_static! { - static ref Q: DMatrix = DMatrix::from_row_slice(5, 5, &[ + pub static ref Q: DMatrix = DMatrix::from_row_slice(5, 5, &[ 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, @@ -277,12 +277,79 @@ pub fn realize_gram( // --- tests --- -#[cfg(test)] -mod tests { +// this problem is from a sangaku by Irisawa Shintarō Hiroatsu. the article +// below includes a nice translation of the problem statement, which was +// recorded in Uchida Itsumi's book _Kokon sankan_ (_Mathematics, Past and +// Present_) +// +// "Japan's 'Wasan' Mathematical Tradition", by Abe Haruki +// https://www.nippon.com/en/japan-topics/c12801/ +// +#[cfg(feature = "dev")] +pub mod irisawa { use std::{array, f64::consts::PI}; use super::*; + pub fn realize_irisawa_hexlet(scaled_tol: f64) -> (DMatrix, bool, DescentHistory) { + let gram = { + let mut gram_to_be = PartialMatrix::new(); + for s in 0..9 { + // each sphere is represented by a spacelike vector + gram_to_be.push_sym(s, s, 1.0); + + // the circumscribing sphere is tangent to all of the other + // spheres, with matching orientation + if s > 0 { + gram_to_be.push_sym(0, s, 1.0); + } + + if s > 2 { + // each chain sphere is tangent to the "sun" and "moon" + // spheres, with opposing orientation + for n in 1..3 { + gram_to_be.push_sym(s, n, -1.0); + } + + // each chain sphere is tangent to the next chain sphere, + // with opposing orientation + let s_next = 3 + (s-2) % 6; + gram_to_be.push_sym(s, s_next, -1.0); + } + } + gram_to_be + }; + + let guess = DMatrix::from_columns( + [ + sphere(0.0, 0.0, 0.0, 15.0), + sphere(0.0, 0.0, -9.0, 5.0), + sphere(0.0, 0.0, 11.0, 3.0) + ].into_iter().chain( + (1..=6).map( + |k| { + let ang = (k as f64) * PI/3.0; + sphere(9.0 * ang.cos(), 9.0 * ang.sin(), 0.0, 2.5) + } + ) + ).collect::>().as_slice() + ); + + // the frozen entries fix the radii of the circumscribing sphere, the + // "sun" and "moon" spheres, and one of the chain spheres + let frozen: [(usize, usize); 4] = array::from_fn(|k| (3, k)); + + realize_gram( + &gram, guess, &frozen, + scaled_tol, 0.5, 0.9, 1.1, 200, 110 + ) + } +} + +#[cfg(test)] +mod tests { + use super::{*, irisawa::realize_irisawa_hexlet}; + #[test] fn sub_proj_test() { let target = PartialMatrix(vec![ @@ -328,182 +395,20 @@ mod tests { assert!(state.loss.abs() < f64::EPSILON); } - // this problem is from a sangaku by Irisawa Shintarō Hiroatsu. the article - // below includes a nice translation of the problem statement, which was - // recorded in Uchida Itsumi's book _Kokon sankan_ (_Mathematics, Past and - // Present_) - // - // "Japan's 'Wasan' Mathematical Tradition", by Abe Haruki - // https://www.nippon.com/en/japan-topics/c12801/ - // #[test] fn irisawa_hexlet_test() { - let gram = PartialMatrix({ - let mut entries = Vec::::new(); - for s in 0..9 { - // each sphere is represented by a spacelike vector - entries.push(MatrixEntry { index: (s, s), value: 1.0 }); - - // the circumscribing sphere is tangent to all of the other - // spheres, with matching orientation - if s > 0 { - entries.push(MatrixEntry { index: (0, s), value: 1.0 }); - entries.push(MatrixEntry { index: (s, 0), value: 1.0 }); - } - - if s > 2 { - // each chain sphere is tangent to the "sun" and "moon" - // spheres, with opposing orientation - for n in 1..3 { - entries.push(MatrixEntry { index: (s, n), value: -1.0 }); - entries.push(MatrixEntry { index: (n, s), value: -1.0 }); - } - - // each chain sphere is tangent to the next chain sphere, - // with opposing orientation - let s_next = 3 + (s-2) % 6; - entries.push(MatrixEntry { index: (s, s_next), value: -1.0 }); - entries.push(MatrixEntry { index: (s_next, s), value: -1.0 }); - } - } - entries - }); - let guess = DMatrix::from_columns( - [ - sphere(0.0, 0.0, 0.0, 15.0), - sphere(0.0, 0.0, -9.0, 5.0), - sphere(0.0, 0.0, 11.0, 3.0) - ].into_iter().chain( - (1..=6).map( - |k| { - let ang = (k as f64) * PI/3.0; - sphere(9.0 * ang.cos(), 9.0 * ang.sin(), 0.0, 2.5) - } - ) - ).collect::>().as_slice() - ); - let frozen: [(usize, usize); 4] = array::from_fn(|k| (3, k)); + // solve Irisawa's problem const SCALED_TOL: f64 = 1.0e-12; - let (config, success, history) = realize_gram( - &gram, guess, &frozen, - SCALED_TOL, 0.5, 0.9, 1.1, 200, 110 - ); + let (config, _, _) = realize_irisawa_hexlet(SCALED_TOL); + + // check against Irisawa's solution 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]; for (k, diam) in solution_diams.into_iter().enumerate() { assert!((config[(3, k)] - 1.0 / diam).abs() < entry_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 { - println!("\nChain diameters:"); - println!(" {} sun (given)", 1.0 / config[(3, 3)]); - for k in 4..9 { - println!(" {} sun", 1.0 / config[(3, k)]); - } - } - println!("\nStep │ Loss\n─────┼────────────────────────────────"); - for (step, scaled_loss) in history.scaled_loss.into_iter().enumerate() { - println!("{:<4} │ {}", step, scaled_loss); - } } - // --- process inspection examples --- - - // these tests are meant for human inspection, not automated use. run them - // one at a time in `--nocapture` mode and read through the results and - // optimization histories that they print out. the `run-examples` script - // will run all of them - - #[test] - fn three_spheres_example() { - let gram = PartialMatrix({ - let mut entries = Vec::::new(); - for j in 0..3 { - for k in 0..3 { - entries.push(MatrixEntry { - index: (j, k), - value: if j == k { 1.0 } else { -1.0 } - }); - } - } - entries - }); - let guess = { - let a: f64 = 0.75_f64.sqrt(); - DMatrix::from_columns(&[ - sphere(1.0, 0.0, 0.0, 1.0), - sphere(-0.5, a, 0.0, 1.0), - sphere(-0.5, -a, 0.0, 1.0) - ]) - }; - println!(); - let (config, success, history) = realize_gram( - &gram, guess, &[], - 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); - } - } - - #[test] - fn point_on_sphere_example() { - let gram = PartialMatrix({ - let mut entries = Vec::::new(); - for j in 0..2 { - for k in 0..2 { - entries.push(MatrixEntry { - index: (j, k), - value: if (j, k) == (1, 1) { 1.0 } else { 0.0 } - }); - } - } - entries - }); - let guess = DMatrix::from_columns(&[ - point(0.0, 0.0, 2.0), - sphere(0.0, 0.0, 0.0, 1.0) - ]); - let frozen = [(3, 0)]; - println!(); - let (config, success, history) = realize_gram( - &gram, guess, &frozen, - 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); - } - } - - /* TO DO */ - // --- new test placed here to avoid merge conflict --- - // at the frozen indices, the optimization steps should have exact zeros, // and the realized configuration should match the initial guess #[test] diff --git a/app-proto/src/lib.rs b/app-proto/src/lib.rs new file mode 100644 index 0000000..0d9bc4a --- /dev/null +++ b/app-proto/src/lib.rs @@ -0,0 +1 @@ +pub mod engine; \ No newline at end of file