diff --git a/.circleci/config.yml b/.circleci/config.yml index ee64d771..b21798d9 100644 --- a/.circleci/config.yml +++ b/.circleci/config.yml @@ -47,10 +47,10 @@ jobs: - checkout - run: name: test - command: cargo test --all-features + command: cargo test --features arbitrary --features serde-serialize --features abomonation-serialize --features sparse --features debug --features io --features compare --features libm - run: name: test nalgebra-glm - command: cargo test -p nalgebra-glm --all-features + command: cargo test -p nalgebra-glm --features arbitrary --features serde-serialize --features abomonation-serialize --features sparse --features debug --features io --features compare --features libm build-wasm: executor: rust-executor steps: diff --git a/.github/FUNDING.yml b/.github/FUNDING.yml index 356ddf23..5ee5fec9 100644 --- a/.github/FUNDING.yml +++ b/.github/FUNDING.yml @@ -1,7 +1,7 @@ # These are supported funding model platforms -github: [ "sebcrozet" ] # Replace with up to 4 GitHub Sponsors-enabled usernames e.g., [user1, user2] -patreon: sebcrozet # Replace with a single Patreon username +github: [ "dimforge" ] # Replace with up to 4 GitHub Sponsors-enabled usernames e.g., [user1, user2] +patreon: # Replace with a single Patreon username open_collective: # Replace with a single Open Collective username ko_fi: # Replace with a single Ko-fi username tidelift: # Replace with a single Tidelift platform-name/package-name e.g., npm/babel diff --git a/CHANGELOG.md b/CHANGELOG.md index 57f18f25..853529f0 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -4,15 +4,37 @@ documented here. This project adheres to [Semantic Versioning](https://semver.org/). -## [0.22.0] - WIP +## [0.23.0] - WIP ### Added + * The `.inverse_transform_unit_vector(v)` was added to `Rotation2/3`, `Isometry2/3`, `UnitQuaternion`, and `UnitComplex`. + It applies the corresponding rotation to a unit vector `Unit`. + * The `Point.map(f)` and `Point.apply(f)` to apply a function to each component of the point, similarly to `Vector.map(f)` + and `Vector.apply(f)`. + * The `Quaternion::from([N; 4])` conversion to build a quaternion from an array of four elements. + * The `Isometry::from(Translation)` conversion to build an isometry from a translation. + * The `Vector::ith_axis(i)` which build a unit vector, e.g., `Unit>` with its i-th component set to 1.0 and the + others set to zero. + * The `Isometry.lerp_slerp` and `Isometry.try_lerp_slerp` methods to interpolate between two isometries using linear + interpolation for the translational part, and spherical interpolation for the rotational part. + * The `Rotation2.slerp`, `Rotation3.slerp`, and `UnitQuaternion.slerp` method for + spherical interpolation. + +## [0.22.0] +In this release, we are using the new version 0.2 of simba. One major change of that version is that the +use of `libm` is now opt-in when building targetting `no-std` environment. If you are using floating-point +operations with nalgebra in a `no-std` environment, you will need to enable the new `libm` feature +of nalgebra for your code to compile again. + +### Added + * The `libm` feature that enables `libm` when building for `no-std` environment. + * The `libm-force` feature that enables `libm` even when building for a not `no-std` environment. * `Cholesky::new_unchecked` which build a Cholesky decomposition without checking that its input is positive-definite. It can be use with SIMD types. * The `Default` trait is now implemented for matrices, and quaternions. They are all filled with zeros, except for `UnitQuaternion` which is initialized with the identity. * Matrix exponential `matrix.exp()`. - + * The `Vector::ith(i, x)` that builds a vector filled with zeros except for the `i`-th component set to `x`. ## [0.21.0] In this release, we are no longer relying on traits from the __alga__ crate for our generic code. diff --git a/Cargo.toml b/Cargo.toml index 4adf2ab1..a36df272 100644 --- a/Cargo.toml +++ b/Cargo.toml @@ -1,6 +1,6 @@ [package] name = "nalgebra" -version = "0.21.1" +version = "0.23.0" authors = [ "Sébastien Crozet " ] description = "Linear algebra library with transformations and statically-sized or dynamically-sized matrices." @@ -10,7 +10,7 @@ repository = "https://github.com/rustsim/nalgebra" readme = "README.md" categories = [ "science" ] keywords = [ "linear", "algebra", "matrix", "vector", "math" ] -license = "BSD-3-Clause" +license = "Apache-2.0" edition = "2018" exclude = ["/ci/*", "/.travis.yml", "/Makefile"] @@ -24,32 +24,36 @@ default = [ "std" ] std = [ "matrixmultiply", "rand/std", "rand_distr", "simba/std" ] stdweb = [ "rand/stdweb" ] arbitrary = [ "quickcheck" ] -serde-serialize = [ "serde", "serde_derive", "num-complex/serde" ] +serde-serialize = [ "serde", "num-complex/serde" ] abomonation-serialize = [ "abomonation" ] sparse = [ ] debug = [ "approx/num-complex", "rand/std" ] alloc = [ ] io = [ "pest", "pest_derive" ] +compare = [ "matrixcompare-core" ] +libm = [ "simba/libm" ] +libm-force = [ "simba/libm_force" ] + [dependencies] -typenum = "1.11" -generic-array = "0.13" +typenum = "1.12" +generic-array = "0.14" rand = { version = "0.7", default-features = false } num-traits = { version = "0.2", default-features = false } -num-complex = { version = "0.2", default-features = false } -num-rational = { version = "0.2", default-features = false } -approx = { version = "0.3", default-features = false } -simba = { version = "0.1", default-features = false } +num-complex = { version = "0.3", default-features = false } +num-rational = { version = "0.3", default-features = false } +approx = { version = "0.4", default-features = false } +simba = { version = "0.3", default-features = false } alga = { version = "0.9", default-features = false, optional = true } -rand_distr = { version = "0.2", optional = true } +rand_distr = { version = "0.3", optional = true } matrixmultiply = { version = "0.2", optional = true } -serde = { version = "1.0", optional = true } -serde_derive = { version = "1.0", optional = true } +serde = { version = "1.0", features = [ "derive" ], optional = true } abomonation = { version = "0.7", optional = true } mint = { version = "0.5", optional = true } quickcheck = { version = "0.9", optional = true } -pest = { version = "2.0", optional = true } -pest_derive = { version = "2.0", optional = true } +pest = { version = "2", optional = true } +pest_derive = { version = "2", optional = true } +matrixcompare-core = { version = "0.1", optional = true } [dev-dependencies] serde_json = "1.0" @@ -61,6 +65,9 @@ rand_isaac = "0.2" ### https://github.com/rust-lang/cargo/issues/4866 #criterion = "0.2.10" +# For matrix comparison macro +matrixcompare = "0.1.3" + [workspace] members = [ "nalgebra-lapack", "nalgebra-glm" ] @@ -71,6 +78,3 @@ path = "benches/lib.rs" [profile.bench] lto = true - -#[patch.crates-io] -#simba = { path = "../simba" } \ No newline at end of file diff --git a/LICENSE b/LICENSE index f7f80044..97f4383a 100644 --- a/LICENSE +++ b/LICENSE @@ -1,27 +1,201 @@ -Copyright (c) 2013, Sébastien Crozet -All rights reserved. + Apache License + Version 2.0, January 2004 + http://www.apache.org/licenses/ -Redistribution and use in source and binary forms, with or without -modification, are permitted provided that the following conditions are met: + TERMS AND CONDITIONS FOR USE, REPRODUCTION, AND DISTRIBUTION -1. Redistributions of source code must retain the above copyright notice, this - list of conditions and the following disclaimer. + 1. Definitions. -2. Redistributions in binary form must reproduce the above copyright notice, - this list of conditions and the following disclaimer in the documentation - and/or other materials provided with the distribution. + "License" shall mean the terms and conditions for use, reproduction, + and distribution as defined by Sections 1 through 9 of this document. -3. Neither the name of the author nor the names of its contributors may be used - to endorse or promote products derived from this software without specific - prior written permission. + "Licensor" shall mean the copyright owner or entity authorized by + the copyright owner that is granting the License. -THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND -ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED -WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE -DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE -FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL -DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR -SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER -CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, -OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE -OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. + "Legal Entity" shall mean the union of the acting entity and all + other entities that control, are controlled by, or are under common + control with that entity. For the purposes of this definition, + "control" means (i) the power, direct or indirect, to cause the + direction or management of such entity, whether by contract or + otherwise, or (ii) ownership of fifty percent (50%) or more of the + outstanding shares, or (iii) beneficial ownership of such entity. + + "You" (or "Your") shall mean an individual or Legal Entity + exercising permissions granted by this License. + + "Source" form shall mean the preferred form for making modifications, + including but not limited to software source code, documentation + source, and configuration files. + + "Object" form shall mean any form resulting from mechanical + transformation or translation of a Source form, including but + not limited to compiled object code, generated documentation, + and conversions to other media types. + + "Work" shall mean the work of authorship, whether in Source or + Object form, made available under the License, as indicated by a + copyright notice that is included in or attached to the work + (an example is provided in the Appendix below). + + "Derivative Works" shall mean any work, whether in Source or Object + form, that is based on (or derived from) the Work and for which the + editorial revisions, annotations, elaborations, or other modifications + represent, as a whole, an original work of authorship. For the purposes + of this License, Derivative Works shall not include works that remain + separable from, or merely link (or bind by name) to the interfaces of, + the Work and Derivative Works thereof. + + "Contribution" shall mean any work of authorship, including + the original version of the Work and any modifications or additions + to that Work or Derivative Works thereof, that is intentionally + submitted to Licensor for inclusion in the Work by the copyright owner + or by an individual or Legal Entity authorized to submit on behalf of + the copyright owner. For the purposes of this definition, "submitted" + means any form of electronic, verbal, or written communication sent + to the Licensor or its representatives, including but not limited to + communication on electronic mailing lists, source code control systems, + and issue tracking systems that are managed by, or on behalf of, the + Licensor for the purpose of discussing and improving the Work, but + excluding communication that is conspicuously marked or otherwise + designated in writing by the copyright owner as "Not a Contribution." + + "Contributor" shall mean Licensor and any individual or Legal Entity + on behalf of whom a Contribution has been received by Licensor and + subsequently incorporated within the Work. + + 2. Grant of Copyright License. Subject to the terms and conditions of + this License, each Contributor hereby grants to You a perpetual, + worldwide, non-exclusive, no-charge, royalty-free, irrevocable + copyright license to reproduce, prepare Derivative Works of, + publicly display, publicly perform, sublicense, and distribute the + Work and such Derivative Works in Source or Object form. + + 3. Grant of Patent License. Subject to the terms and conditions of + this License, each Contributor hereby grants to You a perpetual, + worldwide, non-exclusive, no-charge, royalty-free, irrevocable + (except as stated in this section) patent license to make, have made, + use, offer to sell, sell, import, and otherwise transfer the Work, + where such license applies only to those patent claims licensable + by such Contributor that are necessarily infringed by their + Contribution(s) alone or by combination of their Contribution(s) + with the Work to which such Contribution(s) was submitted. If You + institute patent litigation against any entity (including a + cross-claim or counterclaim in a lawsuit) alleging that the Work + or a Contribution incorporated within the Work constitutes direct + or contributory patent infringement, then any patent licenses + granted to You under this License for that Work shall terminate + as of the date such litigation is filed. + + 4. Redistribution. You may reproduce and distribute copies of the + Work or Derivative Works thereof in any medium, with or without + modifications, and in Source or Object form, provided that You + meet the following conditions: + + (a) You must give any other recipients of the Work or + Derivative Works a copy of this License; and + + (b) You must cause any modified files to carry prominent notices + stating that You changed the files; and + + (c) You must retain, in the Source form of any Derivative Works + that You distribute, all copyright, patent, trademark, and + attribution notices from the Source form of the Work, + excluding those notices that do not pertain to any part of + the Derivative Works; and + + (d) If the Work includes a "NOTICE" text file as part of its + distribution, then any Derivative Works that You distribute must + include a readable copy of the attribution notices contained + within such NOTICE file, excluding those notices that do not + pertain to any part of the Derivative Works, in at least one + of the following places: within a NOTICE text file distributed + as part of the Derivative Works; within the Source form or + documentation, if provided along with the Derivative Works; or, + within a display generated by the Derivative Works, if and + wherever such third-party notices normally appear. The contents + of the NOTICE file are for informational purposes only and + do not modify the License. You may add Your own attribution + notices within Derivative Works that You distribute, alongside + or as an addendum to the NOTICE text from the Work, provided + that such additional attribution notices cannot be construed + as modifying the License. + + You may add Your own copyright statement to Your modifications and + may provide additional or different license terms and conditions + for use, reproduction, or distribution of Your modifications, or + for any such Derivative Works as a whole, provided Your use, + reproduction, and distribution of the Work otherwise complies with + the conditions stated in this License. + + 5. Submission of Contributions. Unless You explicitly state otherwise, + any Contribution intentionally submitted for inclusion in the Work + by You to the Licensor shall be under the terms and conditions of + this License, without any additional terms or conditions. + Notwithstanding the above, nothing herein shall supersede or modify + the terms of any separate license agreement you may have executed + with Licensor regarding such Contributions. + + 6. Trademarks. This License does not grant permission to use the trade + names, trademarks, service marks, or product names of the Licensor, + except as required for reasonable and customary use in describing the + origin of the Work and reproducing the content of the NOTICE file. + + 7. Disclaimer of Warranty. Unless required by applicable law or + agreed to in writing, Licensor provides the Work (and each + Contributor provides its Contributions) on an "AS IS" BASIS, + WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or + implied, including, without limitation, any warranties or conditions + of TITLE, NON-INFRINGEMENT, MERCHANTABILITY, or FITNESS FOR A + PARTICULAR PURPOSE. You are solely responsible for determining the + appropriateness of using or redistributing the Work and assume any + risks associated with Your exercise of permissions under this License. + + 8. Limitation of Liability. In no event and under no legal theory, + whether in tort (including negligence), contract, or otherwise, + unless required by applicable law (such as deliberate and grossly + negligent acts) or agreed to in writing, shall any Contributor be + liable to You for damages, including any direct, indirect, special, + incidental, or consequential damages of any character arising as a + result of this License or out of the use or inability to use the + Work (including but not limited to damages for loss of goodwill, + work stoppage, computer failure or malfunction, or any and all + other commercial damages or losses), even if such Contributor + has been advised of the possibility of such damages. + + 9. Accepting Warranty or Additional Liability. While redistributing + the Work or Derivative Works thereof, You may choose to offer, + and charge a fee for, acceptance of support, warranty, indemnity, + or other liability obligations and/or rights consistent with this + License. However, in accepting such obligations, You may act only + on Your own behalf and on Your sole responsibility, not on behalf + of any other Contributor, and only if You agree to indemnify, + defend, and hold each Contributor harmless for any liability + incurred by, or claims asserted against, such Contributor by reason + of your accepting any such warranty or additional liability. + + END OF TERMS AND CONDITIONS + + APPENDIX: How to apply the Apache License to your work. + + To apply the Apache License to your work, attach the following + boilerplate notice, with the fields enclosed by brackets "[]" + replaced with your own identifying information. (Don't include + the brackets!) The text should be enclosed in the appropriate + comment syntax for the file format. We also recommend that a + file or class name and description of purpose be included on the + same "printed page" as the copyright notice for easier + identification within third-party archives. + + Copyright 2020 Sébastien Crozet + + Licensed under the Apache License, Version 2.0 (the "License"); + you may not use this file except in compliance with the License. + You may obtain a copy of the License at + + http://www.apache.org/licenses/LICENSE-2.0 + + Unless required by applicable law or agreed to in writing, software + distributed under the License is distributed on an "AS IS" BASIS, + WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. + See the License for the specific language governing permissions and + limitations under the License. diff --git a/Makefile b/Makefile index 01707964..0af9f3f2 100644 --- a/Makefile +++ b/Makefile @@ -1,5 +1,5 @@ all: - cargo test --features "debug arbitrary serde-serialize abomonation-serialize" + cargo test --features "debug arbitrary serde-serialize abomonation-serialize compare" # cargo check --features "debug arbitrary serde-serialize" doc: @@ -9,4 +9,4 @@ bench: cargo bench test: - cargo test --features "debug arbitrary serde-serialize abomonation-serialize" + cargo test --features "debug arbitrary serde-serialize abomonation-serialize compare" diff --git a/README.md b/README.md index df3a82f3..72d058b5 100644 --- a/README.md +++ b/README.md @@ -5,14 +5,14 @@ - - Build status + + Build status crates.io - - + +

@@ -29,11 +29,3 @@

----- - -

-  Click this button if you wish to donate to support the development of nalgebra: -

- -

- Become a Patron! -

diff --git a/examples/matrixcompare.rs b/examples/matrixcompare.rs new file mode 100644 index 00000000..713be25f --- /dev/null +++ b/examples/matrixcompare.rs @@ -0,0 +1,82 @@ +extern crate nalgebra as na; + +use matrixcompare::comparators::{AbsoluteElementwiseComparator, ExactElementwiseComparator}; +use matrixcompare::compare_matrices; +use na::{MatrixMN, U3, U4}; + +fn compare_integers_fail() { + println!("Comparing two integer matrices."); + + #[rustfmt::skip] + let a = MatrixMN::<_, U3, U4>::from_row_slice(&[ + 0, 1, 2, 3, + 4, 5, 6, 7, + 8, 9, -2, 11 + ]); + + #[rustfmt::skip] + let b = MatrixMN::<_, U3, U4>::from_row_slice(&[ + 0, 1, 2, 3, + 4, 5, 6, 7, + 8, 9, 10, 11 + ]); + + if let Err(err) = compare_matrices(a, b, &ExactElementwiseComparator) { + println!("{}", err); + } +} + +fn compare_different_size() { + println!("Comparing matrices of different size."); + #[rustfmt::skip] + let a = MatrixMN::<_, U3, U3>::from_row_slice(&[ + 0, 1, 2, + 4, 5, 6, + 8, 9, 10, + ]); + + #[rustfmt::skip] + let b = MatrixMN::<_, U3, U4>::from_row_slice(&[ + 0, 1, 2, 3, + 4, 5, 6, 7, + 8, 9, 10, 11 + ]); + + if let Err(err) = compare_matrices(a, b, &ExactElementwiseComparator) { + println!("{}", err); + } +} + +fn compare_f64_abs_tol_fail() { + println!("Comparing two f64 matrices."); + + #[rustfmt::skip] + let a = MatrixMN::::from_row_slice(&[ + 0.0, 1.0, 2.0 + 1e-10, + 4.0, 5.0, 6.0, + 8.0, 9.0, 10.0, + ]); + + #[rustfmt::skip] + let b = MatrixMN::<_, U3, U3>::from_row_slice(&[ + 0.0, 1.0, 2.0, + 4.0, 5.0, 6.0, + 8.0, 9.0, 10.0 + ]); + + let cmp = AbsoluteElementwiseComparator { tol: 1e-12 }; + if let Err(err) = compare_matrices(a, b, &cmp) { + println!("{}", err); + } +} + +fn main() { + // This example mostly serves the purpose of demonstrating the kind of error messages + // that are given upon comparison failure. + // The more typical use case is using `assert_matrix_eq!` in tests. + compare_integers_fail(); + println!("======================================================"); + compare_f64_abs_tol_fail(); + println!("======================================================"); + compare_different_size(); +} diff --git a/examples/reshaping.rs b/examples/reshaping.rs new file mode 100644 index 00000000..60cb9bdd --- /dev/null +++ b/examples/reshaping.rs @@ -0,0 +1,54 @@ +#![cfg_attr(rustfmt, rustfmt_skip)] + +extern crate nalgebra as na; + +use na::{DMatrix, Dynamic, Matrix2x3, Matrix3x2, U2, U3}; + +fn main() { + // Matrices can be reshaped in-place without moving or copying values. + let m1 = Matrix2x3::new( + 1.1, 1.2, 1.3, + 2.1, 2.2, 2.3 + ); + let m2 = Matrix3x2::new( + 1.1, 2.2, + 2.1, 1.3, + 1.2, 2.3 + ); + + let m3 = m1.reshape_generic(U3, U2); + assert_eq!(m3, m2); + + // Note that, for statically sized matrices, invalid reshapes will not compile: + //let m4 = m3.reshape_generic(U3, U3); + + // If dynamically sized matrices are used, the reshaping is checked at run-time. + let dm1 = DMatrix::from_row_slice( + 4, + 3, + &[ + 1.0, 0.0, 0.0, + 0.0, 0.0, 1.0, + 0.0, 0.0, 0.0, + 0.0, 1.0, 0.0 + ], + ); + let dm2 = DMatrix::from_row_slice( + 6, + 2, + &[ + 1.0, 0.0, + 0.0, 1.0, + 0.0, 0.0, + 0.0, 1.0, + 0.0, 0.0, + 0.0, 0.0, + ], + ); + + let dm3 = dm1.reshape_generic(Dynamic::new(6), Dynamic::new(2)); + assert_eq!(dm3, dm2); + + // Invalid reshapings of dynamic matrices will panic at run-time. + //let dm4 = dm3.reshape_generic(Dynamic::new(6), Dynamic::new(6)); +} diff --git a/nalgebra-glm/Cargo.toml b/nalgebra-glm/Cargo.toml index b72e2988..f4e12f8e 100644 --- a/nalgebra-glm/Cargo.toml +++ b/nalgebra-glm/Cargo.toml @@ -1,6 +1,6 @@ [package] name = "nalgebra-glm" -version = "0.7.0" +version = "0.9.0" authors = ["sebcrozet "] description = "A computer-graphics oriented API for nalgebra, inspired by the C++ GLM library." @@ -23,6 +23,6 @@ abomonation-serialize = [ "nalgebra/abomonation-serialize" ] [dependencies] num-traits = { version = "0.2", default-features = false } -approx = { version = "0.3", default-features = false } -simba = { version = "0.1", default-features = false } -nalgebra = { path = "..", version = "0.21", default-features = false } +approx = { version = "0.4", default-features = false } +simba = { version = "0.3", default-features = false } +nalgebra = { path = "..", version = "0.23", default-features = false } diff --git a/nalgebra-glm/src/ext/matrix_clip_space.rs b/nalgebra-glm/src/ext/matrix_clip_space.rs index 260c5196..5ee725ae 100644 --- a/nalgebra-glm/src/ext/matrix_clip_space.rs +++ b/nalgebra-glm/src/ext/matrix_clip_space.rs @@ -530,7 +530,7 @@ pub fn perspective_lh_no(aspect: N, fovy: N, near: N, far: N) -> T ); assert!( !relative_eq!(aspect, N::zero()), - "The apsect ratio must not be zero." + "The aspect ratio must not be zero." ); let one = N::one(); @@ -566,7 +566,7 @@ pub fn perspective_lh_zo(aspect: N, fovy: N, near: N, far: N) -> T ); assert!( !relative_eq!(aspect, N::zero()), - "The apsect ratio must not be zero." + "The aspect ratio must not be zero." ); let one = N::one(); @@ -632,7 +632,7 @@ pub fn perspective_rh_no(aspect: N, fovy: N, near: N, far: N) -> T ); assert!( !relative_eq!(aspect, N::zero()), - "The apsect ratio must not be zero." + "The aspect ratio must not be zero." ); let negone = -N::one(); @@ -669,7 +669,7 @@ pub fn perspective_rh_zo(aspect: N, fovy: N, near: N, far: N) -> T ); assert!( !relative_eq!(aspect, N::zero()), - "The apsect ratio must not be zero." + "The aspect ratio must not be zero." ); let negone = -N::one(); diff --git a/nalgebra-lapack/Cargo.toml b/nalgebra-lapack/Cargo.toml index b262e8b9..e70359cc 100644 --- a/nalgebra-lapack/Cargo.toml +++ b/nalgebra-lapack/Cargo.toml @@ -1,6 +1,6 @@ [package] name = "nalgebra-lapack" -version = "0.13.0" +version = "0.14.0" authors = [ "Sébastien Crozet ", "Andrew Straw " ] description = "Linear algebra library with transformations and satically-sized or dynamically-sized matrices." @@ -23,10 +23,10 @@ accelerate = ["lapack-src/accelerate"] intel-mkl = ["lapack-src/intel-mkl"] [dependencies] -nalgebra = { version = "0.21", path = ".." } +nalgebra = { version = "0.22" } # , path = ".." } num-traits = "0.2" num-complex = { version = "0.2", default-features = false } -simba = "0.1" +simba = "0.2" serde = { version = "1.0", optional = true } serde_derive = { version = "1.0", optional = true } lapack = { version = "0.16", default-features = false } @@ -34,7 +34,7 @@ lapack-src = { version = "0.5", default-features = false } # clippy = "*" [dev-dependencies] -nalgebra = { version = "0.21", path = "..", features = [ "arbitrary" ] } +nalgebra = { version = "0.22", features = [ "arbitrary" ] } # path = ".." } quickcheck = "0.9" approx = "0.3" rand = "0.7" diff --git a/rustfmt.toml b/rustfmt.toml index 8dc8e61a..e69de29b 100644 --- a/rustfmt.toml +++ b/rustfmt.toml @@ -1,3 +0,0 @@ -unstable_features = true -indent_style = "Block" -where_single_line = true \ No newline at end of file diff --git a/src/base/array_storage.rs b/src/base/array_storage.rs index 1743f06b..b86e3a80 100644 --- a/src/base/array_storage.rs +++ b/src/base/array_storage.rs @@ -24,7 +24,9 @@ use typenum::Prod; use crate::base::allocator::Allocator; use crate::base::default_allocator::DefaultAllocator; use crate::base::dimension::{DimName, U1}; -use crate::base::storage::{ContiguousStorage, ContiguousStorageMut, Owned, Storage, StorageMut}; +use crate::base::storage::{ + ContiguousStorage, ContiguousStorageMut, Owned, ReshapableStorage, Storage, StorageMut, +}; use crate::base::Scalar; /* @@ -267,6 +269,25 @@ where { } +impl ReshapableStorage for ArrayStorage +where + N: Scalar, + R1: DimName, + C1: DimName, + R1::Value: Mul, + Prod: ArrayLength, + R2: DimName, + C2: DimName, + R2::Value: Mul>, + Prod: ArrayLength, +{ + type Output = ArrayStorage; + + fn reshape_generic(self, _: R2, _: C2) -> Self::Output { + ArrayStorage { data: self.data } + } +} + /* * * Allocation-less serde impls. diff --git a/src/base/blas.rs b/src/base/blas.rs index add17af2..f245f24a 100644 --- a/src/base/blas.rs +++ b/src/base/blas.rs @@ -284,7 +284,16 @@ where { assert!( self.nrows() == rhs.nrows(), - "Dot product dimensions mismatch." + "Dot product dimensions mismatch for shapes {:?} and {:?}: left rows != right rows.", + self.shape(), + rhs.shape(), + ); + + assert!( + self.ncols() == rhs.ncols(), + "Dot product dimensions mismatch for shapes {:?} and {:?}: left cols != right cols.", + self.shape(), + rhs.shape(), ); // So we do some special cases for common fixed-size vectors of dimension lower than 8 @@ -361,8 +370,8 @@ where while self.nrows() - i >= 8 { acc0 += unsafe { - conjugate(self.get_unchecked((i + 0, j)).inlined_clone()) - * rhs.get_unchecked((i + 0, j)).inlined_clone() + conjugate(self.get_unchecked((i, j)).inlined_clone()) + * rhs.get_unchecked((i, j)).inlined_clone() }; acc1 += unsafe { conjugate(self.get_unchecked((i + 1, j)).inlined_clone()) @@ -496,8 +505,9 @@ where ShapeConstraint: DimEq + DimEq, { let (nrows, ncols) = self.shape(); - assert!( - (ncols, nrows) == rhs.shape(), + assert_eq!( + (ncols, nrows), + rhs.shape(), "Transposed dot product dimension mismatch." ); diff --git a/src/base/cg.rs b/src/base/cg.rs index 6a7db4bd..935388af 100644 --- a/src/base/cg.rs +++ b/src/base/cg.rs @@ -11,11 +11,12 @@ use crate::base::allocator::Allocator; use crate::base::dimension::{DimName, DimNameDiff, DimNameSub, U1}; use crate::base::storage::{Storage, StorageMut}; use crate::base::{ - DefaultAllocator, Matrix3, Matrix4, MatrixN, Scalar, SquareMatrix, Unit, Vector, Vector3, - VectorN, + DefaultAllocator, Matrix3, Matrix4, MatrixN, Scalar, SquareMatrix, Unit, Vector, Vector2, + Vector3, VectorN, }; use crate::geometry::{ - Isometry, IsometryMatrix3, Orthographic3, Perspective3, Point, Point3, Rotation2, Rotation3, + Isometry, IsometryMatrix3, Orthographic3, Perspective3, Point, Point2, Point3, Rotation2, + Rotation3, }; use simba::scalar::{ClosedAdd, ClosedMul, RealField}; @@ -70,6 +71,26 @@ impl Matrix3 { pub fn new_rotation(angle: N) -> Self { Rotation2::new(angle).to_homogeneous() } + + /// Creates a new homogeneous matrix that applies a scaling factor for each dimension with respect to point. + /// + /// Can be used to implement "zoom_to" functionality. + #[inline] + pub fn new_nonuniform_scaling_wrt_point(scaling: &Vector2, pt: &Point2) -> Self { + let _0 = N::zero(); + let _1 = N::one(); + Matrix3::new( + scaling.x, + _0, + pt.x - pt.x * scaling.x, + _0, + scaling.y, + pt.y - pt.y * scaling.y, + _0, + _0, + _1, + ) + } } impl Matrix4 { @@ -90,6 +111,33 @@ impl Matrix4 { Isometry::rotation_wrt_point(rot, pt).to_homogeneous() } + /// Creates a new homogeneous matrix that applies a scaling factor for each dimension with respect to point. + /// + /// Can be used to implement "zoom_to" functionality. + #[inline] + pub fn new_nonuniform_scaling_wrt_point(scaling: &Vector3, pt: &Point3) -> Self { + let _0 = N::zero(); + let _1 = N::one(); + Matrix4::new( + scaling.x, + _0, + _0, + pt.x - pt.x * scaling.x, + _0, + scaling.y, + _0, + pt.y - pt.y * scaling.y, + _0, + _0, + scaling.z, + pt.z - pt.z * scaling.z, + _0, + _0, + _0, + _1, + ) + } + /// Builds a 3D homogeneous rotation matrix from an axis and an angle (multiplied together). /// /// Returns the identity matrix if the given argument is zero. diff --git a/src/base/construction.rs b/src/base/construction.rs index ae8c10d4..aa4cf956 100644 --- a/src/base/construction.rs +++ b/src/base/construction.rs @@ -1011,6 +1011,20 @@ where N: Scalar + Zero + One, DefaultAllocator: Allocator, { + /// The column vector with `val` as its i-th component. + #[inline] + pub fn ith(i: usize, val: N) -> Self { + let mut res = Self::zeros(); + res[i] = val; + res + } + + /// The column unit vector with `N::one()` as its i-th component. + #[inline] + pub fn ith_axis(i: usize) -> Unit { + Unit::new_unchecked(Self::ith(i, N::one())) + } + /// The column vector with a 1 as its first component, and zero elsewhere. #[inline] pub fn x() -> Self diff --git a/src/base/dimension.rs b/src/base/dimension.rs index c171ff23..3f77998e 100644 --- a/src/base/dimension.rs +++ b/src/base/dimension.rs @@ -23,7 +23,7 @@ impl Dynamic { /// A dynamic size equal to `value`. #[inline] pub fn new(value: usize) -> Self { - Self { value: value } + Self { value } } } diff --git a/src/base/edition.rs b/src/base/edition.rs index b3133648..d69582d9 100644 --- a/src/base/edition.rs +++ b/src/base/edition.rs @@ -13,7 +13,7 @@ use crate::base::dimension::Dynamic; use crate::base::dimension::{ Dim, DimAdd, DimDiff, DimMin, DimMinimum, DimName, DimSub, DimSum, U1, }; -use crate::base::storage::{Storage, StorageMut}; +use crate::base::storage::{ReshapableStorage, Storage, StorageMut}; #[cfg(any(feature = "std", feature = "alloc"))] use crate::base::DMatrix; use crate::base::{DefaultAllocator, Matrix, MatrixMN, RowVector, Scalar, Vector}; @@ -745,7 +745,7 @@ impl> Matrix { self.resize_generic(R2::name(), C2::name(), val) } - /// Resizes `self` such that it has dimensions `new_nrows × now_ncols`. + /// Resizes `self` such that it has dimensions `new_nrows × new_ncols`. /// /// The values are copied such that `self[(i, j)] == result[(i, j)]`. If the result has more /// rows and/or columns than `self`, then the extra rows or columns are filled with `val`. @@ -813,6 +813,80 @@ impl> Matrix { } } +impl Matrix +where + N: Scalar, + R: Dim, + C: Dim, +{ + /// Reshapes `self` such that it has dimensions `new_nrows × new_ncols`. + /// + /// This will reinterpret `self` as if it is a matrix with `new_nrows` rows and `new_ncols` + /// columns. The arrangements of the component in the output matrix are the same as what + /// would be obtained by `Matrix::from_slice_generic(self.as_slice(), new_nrows, new_ncols)`. + /// + /// If `self` is a dynamically-sized matrix, then its components are neither copied nor moved. + /// If `self` is staticyll-sized, then a copy may happen in some situations. + /// This function will panic if the given dimensions are such that the number of elements of + /// the input matrix are not equal to the number of elements of the output matrix. + /// + /// # Examples + /// + /// ``` + /// # use nalgebra::{Matrix3x2, Matrix2x3, DMatrix, U2, U3, Dynamic}; + /// + /// let m1 = Matrix2x3::new( + /// 1.1, 1.2, 1.3, + /// 2.1, 2.2, 2.3 + /// ); + /// let m2 = Matrix3x2::new( + /// 1.1, 2.2, + /// 2.1, 1.3, + /// 1.2, 2.3 + /// ); + /// let reshaped = m1.reshape_generic(U3, U2); + /// assert_eq!(reshaped, m2); + /// + /// let dm1 = DMatrix::from_row_slice( + /// 4, + /// 3, + /// &[ + /// 1.0, 0.0, 0.0, + /// 0.0, 0.0, 1.0, + /// 0.0, 0.0, 0.0, + /// 0.0, 1.0, 0.0 + /// ], + /// ); + /// let dm2 = DMatrix::from_row_slice( + /// 6, + /// 2, + /// &[ + /// 1.0, 0.0, + /// 0.0, 1.0, + /// 0.0, 0.0, + /// 0.0, 1.0, + /// 0.0, 0.0, + /// 0.0, 0.0, + /// ], + /// ); + /// let reshaped = dm1.reshape_generic(Dynamic::new(6), Dynamic::new(2)); + /// assert_eq!(reshaped, dm2); + /// ``` + pub fn reshape_generic( + self, + new_nrows: R2, + new_ncols: C2, + ) -> Matrix + where + R2: Dim, + C2: Dim, + S: ReshapableStorage, + { + let data = self.data.reshape_generic(new_nrows, new_ncols); + Matrix::from_data(data) + } +} + #[cfg(any(feature = "std", feature = "alloc"))] impl DMatrix { /// Resizes this matrix in-place. diff --git a/src/base/matrix.rs b/src/base/matrix.rs index e4821bf8..ca7fa77b 100644 --- a/src/base/matrix.rs +++ b/src/base/matrix.rs @@ -159,13 +159,39 @@ impl Abomonation for Matrix> matrixcompare_core::Matrix + for Matrix +{ + fn rows(&self) -> usize { + self.nrows() + } + + fn cols(&self) -> usize { + self.ncols() + } + + fn access(&self) -> matrixcompare_core::Access { + matrixcompare_core::Access::Dense(self) + } +} + +#[cfg(feature = "compare")] +impl> matrixcompare_core::DenseAccess + for Matrix +{ + fn fetch_single(&self, row: usize, col: usize) -> N { + self.index((row, col)).clone() + } +} + impl Matrix { /// Creates a new matrix with the given data without statically checking that the matrix /// dimension matches the storage dimension. #[inline] pub unsafe fn from_data_statically_unchecked(data: S) -> Matrix { Matrix { - data: data, + data, _phantoms: PhantomData, } } @@ -538,8 +564,9 @@ impl> Matrix { let mut res = unsafe { MatrixMN::new_uninitialized_generic(nrows, ncols) }; - assert!( - (nrows.value(), ncols.value()) == rhs.shape(), + assert_eq!( + (nrows.value(), ncols.value()), + rhs.shape(), "Matrix simultaneous traversal error: dimension mismatch." ); @@ -578,9 +605,14 @@ impl> Matrix { let mut res = unsafe { MatrixMN::new_uninitialized_generic(nrows, ncols) }; - assert!( - (nrows.value(), ncols.value()) == b.shape() - && (nrows.value(), ncols.value()) == c.shape(), + assert_eq!( + (nrows.value(), ncols.value()), + b.shape(), + "Matrix simultaneous traversal error: dimension mismatch." + ); + assert_eq!( + (nrows.value(), ncols.value()), + c.shape(), "Matrix simultaneous traversal error: dimension mismatch." ); @@ -636,8 +668,9 @@ impl> Matrix { let mut res = init; - assert!( - (nrows.value(), ncols.value()) == rhs.shape(), + assert_eq!( + (nrows.value(), ncols.value()), + rhs.shape(), "Matrix simultaneous traversal error: dimension mismatch." ); @@ -884,8 +917,9 @@ impl> Matrix { { let (nrows, ncols) = self.shape(); - assert!( - (nrows, ncols) == rhs.shape(), + assert_eq!( + (nrows, ncols), + rhs.shape(), "Matrix simultaneous traversal error: dimension mismatch." ); @@ -922,12 +956,14 @@ impl> Matrix { { let (nrows, ncols) = self.shape(); - assert!( - (nrows, ncols) == b.shape(), + assert_eq!( + (nrows, ncols), + b.shape(), "Matrix simultaneous traversal error: dimension mismatch." ); - assert!( - (nrows, ncols) == c.shape(), + assert_eq!( + (nrows, ncols), + c.shape(), "Matrix simultaneous traversal error: dimension mismatch." ); @@ -1427,8 +1463,9 @@ where #[inline] fn lt(&self, right: &Self) -> bool { - assert!( - self.shape() == right.shape(), + assert_eq!( + self.shape(), + right.shape(), "Matrix comparison error: dimensions mismatch." ); self.iter().zip(right.iter()).all(|(a, b)| a.lt(b)) @@ -1436,8 +1473,9 @@ where #[inline] fn le(&self, right: &Self) -> bool { - assert!( - self.shape() == right.shape(), + assert_eq!( + self.shape(), + right.shape(), "Matrix comparison error: dimensions mismatch." ); self.iter().zip(right.iter()).all(|(a, b)| a.le(b)) @@ -1445,8 +1483,9 @@ where #[inline] fn gt(&self, right: &Self) -> bool { - assert!( - self.shape() == right.shape(), + assert_eq!( + self.shape(), + right.shape(), "Matrix comparison error: dimensions mismatch." ); self.iter().zip(right.iter()).all(|(a, b)| a.gt(b)) @@ -1454,8 +1493,9 @@ where #[inline] fn ge(&self, right: &Self) -> bool { - assert!( - self.shape() == right.shape(), + assert_eq!( + self.shape(), + right.shape(), "Matrix comparison error: dimensions mismatch." ); self.iter().zip(right.iter()).all(|(a, b)| a.ge(b)) @@ -1602,7 +1642,11 @@ impl + SameNumberOfColumns, { - assert!(self.shape() == (2, 1), "2D perpendicular product "); + assert!( + self.shape() == (2, 1), + "2D perpendicular product requires (2, 1) vector but found {:?}", + self.shape() + ); unsafe { self.get_unchecked((0, 0)).inlined_clone() * b.get_unchecked((1, 0)).inlined_clone() @@ -1626,13 +1670,11 @@ impl + SameNumberOfColumns, { let shape = self.shape(); - assert!( - shape == b.shape(), - "Vector cross product dimension mismatch." - ); + assert_eq!(shape, b.shape(), "Vector cross product dimension mismatch."); assert!( (shape.0 == 3 && shape.1 == 1) || (shape.0 == 1 && shape.1 == 3), - "Vector cross product dimension mismatch." + "Vector cross product dimension mismatch: must be (3, 1) or (1, 3) but found {:?}.", + shape ); if shape.0 == 3 { diff --git a/src/base/ops.rs b/src/base/ops.rs index 12b26d1a..7f0c6e4e 100644 --- a/src/base/ops.rs +++ b/src/base/ops.rs @@ -154,8 +154,8 @@ macro_rules! componentwise_binop_impl( out: &mut Matrix) where SB: Storage, SC: StorageMut { - assert!(self.shape() == rhs.shape(), "Matrix addition/subtraction dimensions mismatch."); - assert!(self.shape() == out.shape(), "Matrix addition/subtraction output dimensions mismatch."); + assert_eq!(self.shape(), rhs.shape(), "Matrix addition/subtraction dimensions mismatch."); + assert_eq!(self.shape(), out.shape(), "Matrix addition/subtraction output dimensions mismatch."); // This is the most common case and should be deduced at compile-time. // FIXME: use specialization instead? @@ -188,7 +188,7 @@ macro_rules! componentwise_binop_impl( C2: Dim, SA: StorageMut, SB: Storage { - assert!(self.shape() == rhs.shape(), "Matrix addition/subtraction dimensions mismatch."); + assert_eq!(self.shape(), rhs.shape(), "Matrix addition/subtraction dimensions mismatch."); // This is the most common case and should be deduced at compile-time. // FIXME: use specialization instead? @@ -218,7 +218,7 @@ macro_rules! componentwise_binop_impl( where R2: Dim, C2: Dim, SB: StorageMut { - assert!(self.shape() == rhs.shape(), "Matrix addition/subtraction dimensions mismatch."); + assert_eq!(self.shape(), rhs.shape(), "Matrix addition/subtraction dimensions mismatch."); // This is the most common case and should be deduced at compile-time. // FIXME: use specialization instead? @@ -277,7 +277,7 @@ macro_rules! componentwise_binop_impl( #[inline] fn $method(self, rhs: &'b Matrix) -> Self::Output { - assert!(self.shape() == rhs.shape(), "Matrix addition/subtraction dimensions mismatch."); + assert_eq!(self.shape(), rhs.shape(), "Matrix addition/subtraction dimensions mismatch."); let mut res = self.into_owned_sum::(); res.$method_assign_statically_unchecked(rhs); res @@ -296,7 +296,7 @@ macro_rules! componentwise_binop_impl( #[inline] fn $method(self, rhs: Matrix) -> Self::Output { let mut rhs = rhs.into_owned_sum::(); - assert!(self.shape() == rhs.shape(), "Matrix addition/subtraction dimensions mismatch."); + assert_eq!(self.shape(), rhs.shape(), "Matrix addition/subtraction dimensions mismatch."); self.$method_assign_statically_unchecked_rhs(&mut rhs); rhs } @@ -728,11 +728,21 @@ where assert!( nrows1 == nrows2, - "Matrix multiplication dimensions mismatch." + "Matrix multiplication dimensions mismatch {:?} and {:?}: left rows != right rows.", + self.shape(), + rhs.shape() ); assert!( - nrows3 == ncols1 && ncols3 == ncols2, - "Matrix multiplication output dimensions mismatch." + ncols1 == nrows3, + "Matrix multiplication output dimensions mismatch {:?} and {:?}: left cols != right rows.", + self.shape(), + out.shape() + ); + assert!( + ncols2 == ncols3, + "Matrix multiplication output dimensions mismatch {:?} and {:?}: left cols != right cols", + rhs.shape(), + out.shape() ); for i in 0..ncols1 { diff --git a/src/base/storage.rs b/src/base/storage.rs index 9f039ba2..598cb061 100644 --- a/src/base/storage.rs +++ b/src/base/storage.rs @@ -171,7 +171,7 @@ pub unsafe trait StorageMut: Storage { /// A matrix storage that is stored contiguously in memory. /// -/// The storage requirement means that for any value of `i` in `[0, nrows * ncols[`, the value +/// The storage requirement means that for any value of `i` in `[0, nrows * ncols - 1]`, the value /// `.get_unchecked_linear` returns one of the matrix component. This trait is unsafe because /// failing to comply to this may cause Undefined Behaviors. pub unsafe trait ContiguousStorage: @@ -181,10 +181,26 @@ pub unsafe trait ContiguousStorage: /// A mutable matrix storage that is stored contiguously in memory. /// -/// The storage requirement means that for any value of `i` in `[0, nrows * ncols[`, the value +/// The storage requirement means that for any value of `i` in `[0, nrows * ncols - 1]`, the value /// `.get_unchecked_linear` returns one of the matrix component. This trait is unsafe because /// failing to comply to this may cause Undefined Behaviors. pub unsafe trait ContiguousStorageMut: ContiguousStorage + StorageMut { } + +/// A matrix storage that can be reshaped in-place. +pub trait ReshapableStorage: Storage +where + N: Scalar, + R1: Dim, + C1: Dim, + R2: Dim, + C2: Dim, +{ + /// The reshaped storage type. + type Output: Storage; + + /// Reshapes the storage into the output storage type. + fn reshape_generic(self, nrows: R2, ncols: C2) -> Self::Output; +} diff --git a/src/base/vec_storage.rs b/src/base/vec_storage.rs index 26948693..37f31213 100644 --- a/src/base/vec_storage.rs +++ b/src/base/vec_storage.rs @@ -8,7 +8,9 @@ use crate::base::allocator::Allocator; use crate::base::constraint::{SameNumberOfRows, ShapeConstraint}; use crate::base::default_allocator::DefaultAllocator; use crate::base::dimension::{Dim, DimName, Dynamic, U1}; -use crate::base::storage::{ContiguousStorage, ContiguousStorageMut, Owned, Storage, StorageMut}; +use crate::base::storage::{ + ContiguousStorage, ContiguousStorageMut, Owned, ReshapableStorage, Storage, StorageMut, +}; use crate::base::{Scalar, Vector}; #[cfg(feature = "abomonation-serialize")] @@ -41,11 +43,7 @@ impl VecStorage { nrows.value() * ncols.value() == data.len(), "Data storage buffer dimension mismatch." ); - Self { - data: data, - nrows: nrows, - ncols: ncols, - } + Self { data, nrows, ncols } } /// The underlying data storage. @@ -229,6 +227,42 @@ unsafe impl ContiguousStorageMut for VecStorag { } +impl ReshapableStorage for VecStorage +where + N: Scalar, + C1: Dim, + C2: Dim, +{ + type Output = VecStorage; + + fn reshape_generic(self, nrows: Dynamic, ncols: C2) -> Self::Output { + assert_eq!(nrows.value() * ncols.value(), self.data.len()); + VecStorage { + data: self.data, + nrows, + ncols, + } + } +} + +impl ReshapableStorage for VecStorage +where + N: Scalar, + C1: Dim, + R2: DimName, +{ + type Output = VecStorage; + + fn reshape_generic(self, nrows: R2, ncols: Dynamic) -> Self::Output { + assert_eq!(nrows.value() * ncols.value(), self.data.len()); + VecStorage { + data: self.data, + nrows, + ncols, + } + } +} + unsafe impl StorageMut for VecStorage where DefaultAllocator: Allocator, @@ -244,6 +278,42 @@ where } } +impl ReshapableStorage for VecStorage +where + N: Scalar, + R1: DimName, + C2: Dim, +{ + type Output = VecStorage; + + fn reshape_generic(self, nrows: Dynamic, ncols: C2) -> Self::Output { + assert_eq!(nrows.value() * ncols.value(), self.data.len()); + VecStorage { + data: self.data, + nrows, + ncols, + } + } +} + +impl ReshapableStorage for VecStorage +where + N: Scalar, + R1: DimName, + R2: DimName, +{ + type Output = VecStorage; + + fn reshape_generic(self, nrows: R2, ncols: Dynamic) -> Self::Output { + assert_eq!(nrows.value() * ncols.value(), self.data.len()); + VecStorage { + data: self.data, + nrows, + ncols, + } + } +} + #[cfg(feature = "abomonation-serialize")] impl Abomonation for VecStorage { unsafe fn entomb(&self, writer: &mut W) -> IOResult<()> { diff --git a/src/geometry/abstract_rotation.rs b/src/geometry/abstract_rotation.rs index ff057a09..52471851 100644 --- a/src/geometry/abstract_rotation.rs +++ b/src/geometry/abstract_rotation.rs @@ -1,6 +1,6 @@ use crate::allocator::Allocator; use crate::geometry::{Rotation, UnitComplex, UnitQuaternion}; -use crate::{DefaultAllocator, DimName, Point, Scalar, SimdRealField, VectorN, U2, U3}; +use crate::{DefaultAllocator, DimName, Point, Scalar, SimdRealField, Unit, VectorN, U2, U3}; use simba::scalar::ClosedMul; @@ -24,6 +24,13 @@ pub trait AbstractRotation: PartialEq + ClosedMul + Clone fn inverse_transform_vector(&self, v: &VectorN) -> VectorN where DefaultAllocator: Allocator; + /// Apply the inverse rotation to the given unit vector. + fn inverse_transform_unit_vector(&self, v: &Unit>) -> Unit> + where + DefaultAllocator: Allocator, + { + Unit::new_unchecked(self.inverse_transform_vector(&**v)) + } /// Apply the inverse rotation to the given point. fn inverse_transform_point(&self, p: &Point) -> Point where @@ -74,6 +81,14 @@ where self.inverse_transform_vector(v) } + #[inline] + fn inverse_transform_unit_vector(&self, v: &Unit>) -> Unit> + where + DefaultAllocator: Allocator, + { + self.inverse_transform_unit_vector(v) + } + #[inline] fn inverse_transform_point(&self, p: &Point) -> Point where diff --git a/src/geometry/isometry.rs b/src/geometry/isometry.rs index 2e4b7cd5..edb27b5a 100755 --- a/src/geometry/isometry.rs +++ b/src/geometry/isometry.rs @@ -14,10 +14,12 @@ use simba::scalar::{RealField, SubsetOf}; use simba::simd::SimdRealField; use crate::base::allocator::Allocator; -use crate::base::dimension::{DimName, DimNameAdd, DimNameSum, U1}; +use crate::base::dimension::{DimName, DimNameAdd, DimNameSum, U1, U2, U3}; use crate::base::storage::Owned; -use crate::base::{DefaultAllocator, MatrixN, Scalar, VectorN}; -use crate::geometry::{AbstractRotation, Point, Translation}; +use crate::base::{DefaultAllocator, MatrixN, Scalar, Unit, VectorN}; +use crate::geometry::{ + AbstractRotation, Point, Rotation2, Rotation3, Translation, UnitComplex, UnitQuaternion, +}; /// A direct isometry, i.e., a rotation followed by a translation, aka. a rigid-body motion, aka. an element of a Special Euclidean (SE) group. #[repr(C)] @@ -82,21 +84,23 @@ where } } -impl + Copy> Copy - for Isometry +impl Copy for Isometry where DefaultAllocator: Allocator, Owned: Copy, { } -impl + Clone> Clone for Isometry +impl Clone for Isometry where DefaultAllocator: Allocator, { #[inline] fn clone(&self) -> Self { - Self::from_parts(self.translation.clone(), self.rotation.clone()) + Self { + rotation: self.rotation.clone(), + translation: self.translation.clone(), + } } } @@ -348,6 +352,237 @@ where pub fn inverse_transform_vector(&self, v: &VectorN) -> VectorN { self.rotation.inverse_transform_vector(v) } + + /// Transform the given unit vector by the inverse of this isometry, ignoring the + /// translation component of the isometry. This may be + /// less expensive than computing the entire isometry inverse and then + /// transforming the point. + /// + /// # Example + /// + /// ``` + /// # #[macro_use] extern crate approx; + /// # use std::f32; + /// # use nalgebra::{Isometry3, Translation3, UnitQuaternion, Vector3}; + /// let tra = Translation3::new(0.0, 0.0, 3.0); + /// let rot = UnitQuaternion::from_scaled_axis(Vector3::z() * f32::consts::FRAC_PI_2); + /// let iso = Isometry3::from_parts(tra, rot); + /// + /// let transformed_point = iso.inverse_transform_unit_vector(&Vector3::x_axis()); + /// assert_relative_eq!(transformed_point, -Vector3::y_axis(), epsilon = 1.0e-6); + /// ``` + #[inline] + pub fn inverse_transform_unit_vector(&self, v: &Unit>) -> Unit> { + self.rotation.inverse_transform_unit_vector(v) + } +} + +impl Isometry> { + /// Interpolates between two isometries using a linear interpolation for the translation part, + /// and a spherical interpolation for the rotation part. + /// + /// Panics if the angle between both rotations is 180 degrees (in which case the interpolation + /// is not well-defined). Use `.try_lerp_slerp` instead to avoid the panic. + /// + /// # Examples: + /// + /// ``` + /// # use nalgebra::{Vector3, Translation3, Isometry3, UnitQuaternion}; + /// + /// let t1 = Translation3::new(1.0, 2.0, 3.0); + /// let t2 = Translation3::new(4.0, 8.0, 12.0); + /// let q1 = UnitQuaternion::from_euler_angles(std::f32::consts::FRAC_PI_4, 0.0, 0.0); + /// let q2 = UnitQuaternion::from_euler_angles(-std::f32::consts::PI, 0.0, 0.0); + /// let iso1 = Isometry3::from_parts(t1, q1); + /// let iso2 = Isometry3::from_parts(t2, q2); + /// + /// let iso3 = iso1.lerp_slerp(&iso2, 1.0 / 3.0); + /// + /// assert_eq!(iso3.translation.vector, Vector3::new(2.0, 4.0, 6.0)); + /// assert_eq!(iso3.rotation.euler_angles(), (std::f32::consts::FRAC_PI_2, 0.0, 0.0)); + /// ``` + #[inline] + pub fn lerp_slerp(&self, other: &Self, t: N) -> Self + where + N: RealField, + { + let tr = self.translation.vector.lerp(&other.translation.vector, t); + let rot = self.rotation.slerp(&other.rotation, t); + Self::from_parts(tr.into(), rot) + } + + /// Attempts to interpolate between two isometries using a linear interpolation for the translation part, + /// and a spherical interpolation for the rotation part. + /// + /// Retuns `None` if the angle between both rotations is 180 degrees (in which case the interpolation + /// is not well-defined). + /// + /// # Examples: + /// + /// ``` + /// # use nalgebra::{Vector3, Translation3, Isometry3, UnitQuaternion}; + /// + /// let t1 = Translation3::new(1.0, 2.0, 3.0); + /// let t2 = Translation3::new(4.0, 8.0, 12.0); + /// let q1 = UnitQuaternion::from_euler_angles(std::f32::consts::FRAC_PI_4, 0.0, 0.0); + /// let q2 = UnitQuaternion::from_euler_angles(-std::f32::consts::PI, 0.0, 0.0); + /// let iso1 = Isometry3::from_parts(t1, q1); + /// let iso2 = Isometry3::from_parts(t2, q2); + /// + /// let iso3 = iso1.lerp_slerp(&iso2, 1.0 / 3.0); + /// + /// assert_eq!(iso3.translation.vector, Vector3::new(2.0, 4.0, 6.0)); + /// assert_eq!(iso3.rotation.euler_angles(), (std::f32::consts::FRAC_PI_2, 0.0, 0.0)); + /// ``` + #[inline] + pub fn try_lerp_slerp(&self, other: &Self, t: N, epsilon: N) -> Option + where + N: RealField, + { + let tr = self.translation.vector.lerp(&other.translation.vector, t); + let rot = self.rotation.try_slerp(&other.rotation, t, epsilon)?; + Some(Self::from_parts(tr.into(), rot)) + } +} + +impl Isometry> { + /// Interpolates between two isometries using a linear interpolation for the translation part, + /// and a spherical interpolation for the rotation part. + /// + /// Panics if the angle between both rotations is 180 degrees (in which case the interpolation + /// is not well-defined). Use `.try_lerp_slerp` instead to avoid the panic. + /// + /// # Examples: + /// + /// ``` + /// # use nalgebra::{Vector3, Translation3, Rotation3, IsometryMatrix3}; + /// + /// let t1 = Translation3::new(1.0, 2.0, 3.0); + /// let t2 = Translation3::new(4.0, 8.0, 12.0); + /// let q1 = Rotation3::from_euler_angles(std::f32::consts::FRAC_PI_4, 0.0, 0.0); + /// let q2 = Rotation3::from_euler_angles(-std::f32::consts::PI, 0.0, 0.0); + /// let iso1 = IsometryMatrix3::from_parts(t1, q1); + /// let iso2 = IsometryMatrix3::from_parts(t2, q2); + /// + /// let iso3 = iso1.lerp_slerp(&iso2, 1.0 / 3.0); + /// + /// assert_eq!(iso3.translation.vector, Vector3::new(2.0, 4.0, 6.0)); + /// assert_eq!(iso3.rotation.euler_angles(), (std::f32::consts::FRAC_PI_2, 0.0, 0.0)); + /// ``` + #[inline] + pub fn lerp_slerp(&self, other: &Self, t: N) -> Self + where + N: RealField, + { + let tr = self.translation.vector.lerp(&other.translation.vector, t); + let rot = self.rotation.slerp(&other.rotation, t); + Self::from_parts(tr.into(), rot) + } + + /// Attempts to interpolate between two isometries using a linear interpolation for the translation part, + /// and a spherical interpolation for the rotation part. + /// + /// Retuns `None` if the angle between both rotations is 180 degrees (in which case the interpolation + /// is not well-defined). + /// + /// # Examples: + /// + /// ``` + /// # use nalgebra::{Vector3, Translation3, Rotation3, IsometryMatrix3}; + /// + /// let t1 = Translation3::new(1.0, 2.0, 3.0); + /// let t2 = Translation3::new(4.0, 8.0, 12.0); + /// let q1 = Rotation3::from_euler_angles(std::f32::consts::FRAC_PI_4, 0.0, 0.0); + /// let q2 = Rotation3::from_euler_angles(-std::f32::consts::PI, 0.0, 0.0); + /// let iso1 = IsometryMatrix3::from_parts(t1, q1); + /// let iso2 = IsometryMatrix3::from_parts(t2, q2); + /// + /// let iso3 = iso1.lerp_slerp(&iso2, 1.0 / 3.0); + /// + /// assert_eq!(iso3.translation.vector, Vector3::new(2.0, 4.0, 6.0)); + /// assert_eq!(iso3.rotation.euler_angles(), (std::f32::consts::FRAC_PI_2, 0.0, 0.0)); + /// ``` + #[inline] + pub fn try_lerp_slerp(&self, other: &Self, t: N, epsilon: N) -> Option + where + N: RealField, + { + let tr = self.translation.vector.lerp(&other.translation.vector, t); + let rot = self.rotation.try_slerp(&other.rotation, t, epsilon)?; + Some(Self::from_parts(tr.into(), rot)) + } +} + +impl Isometry> { + /// Interpolates between two isometries using a linear interpolation for the translation part, + /// and a spherical interpolation for the rotation part. + /// + /// Panics if the angle between both rotations is 180 degrees (in which case the interpolation + /// is not well-defined). Use `.try_lerp_slerp` instead to avoid the panic. + /// + /// # Examples: + /// + /// ``` + /// # #[macro_use] extern crate approx; + /// # use nalgebra::{Vector2, Translation2, UnitComplex, Isometry2}; + /// + /// let t1 = Translation2::new(1.0, 2.0); + /// let t2 = Translation2::new(4.0, 8.0); + /// let q1 = UnitComplex::new(std::f32::consts::FRAC_PI_4); + /// let q2 = UnitComplex::new(-std::f32::consts::PI); + /// let iso1 = Isometry2::from_parts(t1, q1); + /// let iso2 = Isometry2::from_parts(t2, q2); + /// + /// let iso3 = iso1.lerp_slerp(&iso2, 1.0 / 3.0); + /// + /// assert_eq!(iso3.translation.vector, Vector2::new(2.0, 4.0)); + /// assert_relative_eq!(iso3.rotation.angle(), std::f32::consts::FRAC_PI_2); + /// ``` + #[inline] + pub fn lerp_slerp(&self, other: &Self, t: N) -> Self + where + N: RealField, + { + let tr = self.translation.vector.lerp(&other.translation.vector, t); + let rot = self.rotation.slerp(&other.rotation, t); + Self::from_parts(tr.into(), rot) + } +} + +impl Isometry> { + /// Interpolates between two isometries using a linear interpolation for the translation part, + /// and a spherical interpolation for the rotation part. + /// + /// Panics if the angle between both rotations is 180 degrees (in which case the interpolation + /// is not well-defined). Use `.try_lerp_slerp` instead to avoid the panic. + /// + /// # Examples: + /// + /// ``` + /// # #[macro_use] extern crate approx; + /// # use nalgebra::{Vector2, Translation2, Rotation2, IsometryMatrix2}; + /// + /// let t1 = Translation2::new(1.0, 2.0); + /// let t2 = Translation2::new(4.0, 8.0); + /// let q1 = Rotation2::new(std::f32::consts::FRAC_PI_4); + /// let q2 = Rotation2::new(-std::f32::consts::PI); + /// let iso1 = IsometryMatrix2::from_parts(t1, q1); + /// let iso2 = IsometryMatrix2::from_parts(t2, q2); + /// + /// let iso3 = iso1.lerp_slerp(&iso2, 1.0 / 3.0); + /// + /// assert_eq!(iso3.translation.vector, Vector2::new(2.0, 4.0)); + /// assert_relative_eq!(iso3.rotation.angle(), std::f32::consts::FRAC_PI_2); + /// ``` + #[inline] + pub fn lerp_slerp(&self, other: &Self, t: N) -> Self + where + N: RealField, + { + let tr = self.translation.vector.lerp(&other.translation.vector, t); + let rot = self.rotation.slerp(&other.rotation, t); + Self::from_parts(tr.into(), rot) + } } // NOTE: we don't require `R: Rotation<...>` here because this is not useful for the implementation diff --git a/src/geometry/isometry_alias.rs b/src/geometry/isometry_alias.rs index ba9a69e7..69d33cb4 100644 --- a/src/geometry/isometry_alias.rs +++ b/src/geometry/isometry_alias.rs @@ -13,3 +13,15 @@ pub type IsometryMatrix2 = Isometry>; /// A 3-dimensional direct isometry using a rotation matrix for its rotational part. Also known as a rigid-body motion, or as an element of SE(3). pub type IsometryMatrix3 = Isometry>; + +// This tests that the types correctly implement `Copy`, without having to run tests +// (when targeting no-std for example). +#[allow(dead_code)] +fn ensure_copy() { + fn is_copy() {} + + is_copy::>(); + is_copy::>(); + is_copy::>(); + is_copy::>(); +} diff --git a/src/geometry/isometry_conversion.rs b/src/geometry/isometry_conversion.rs index 7cb6fe85..bdb128f9 100644 --- a/src/geometry/isometry_conversion.rs +++ b/src/geometry/isometry_conversion.rs @@ -151,6 +151,17 @@ where } } +impl> From> + for Isometry +where + DefaultAllocator: Allocator, +{ + #[inline] + fn from(tra: Translation) -> Self { + Self::from_parts(tra, R::identity()) + } +} + impl From> for MatrixN> where D: DimNameAdd, diff --git a/src/geometry/isometry_ops.rs b/src/geometry/isometry_ops.rs index 8b00fcda..e676d492 100644 --- a/src/geometry/isometry_ops.rs +++ b/src/geometry/isometry_ops.rs @@ -149,6 +149,7 @@ isometry_binop_impl_all!( [ref ref] => { let shift = self.rotation.transform_vector(&rhs.translation.vector); + #[allow(clippy::suspicious_arithmetic_impl)] Isometry::from_parts(Translation::from(&self.translation.vector + shift), self.rotation.clone() * rhs.rotation.clone()) // FIXME: too bad we have to clone. }; @@ -157,10 +158,10 @@ isometry_binop_impl_all!( isometry_binop_impl_all!( Div, div; self: Isometry, rhs: Isometry, Output = Isometry; - [val val] => self * rhs.inverse(); - [ref val] => self * rhs.inverse(); - [val ref] => self * rhs.inverse(); - [ref ref] => self * rhs.inverse(); + [val val] => #[allow(clippy::suspicious_arithmetic_impl)] { self * rhs.inverse() }; + [ref val] => #[allow(clippy::suspicious_arithmetic_impl)] { self * rhs.inverse() }; + [val ref] => #[allow(clippy::suspicious_arithmetic_impl)] { self * rhs.inverse() }; + [ref ref] => #[allow(clippy::suspicious_arithmetic_impl)] { self * rhs.inverse() }; ); // Isometry ×= Translation @@ -289,6 +290,7 @@ isometry_binop_impl_all!( [ref val] => self * &right; [val ref] => &self * right; [ref ref] => { + #[allow(clippy::suspicious_arithmetic_impl)] let new_tr = &self.translation.vector + self.rotation.transform_vector(&right.vector); Isometry::from_parts(Translation::from(new_tr), self.rotation.clone()) }; @@ -427,10 +429,10 @@ isometry_from_composition_impl_all!( self: Rotation, right: Isometry>, Output = Isometry>; // FIXME: don't call inverse explicitly? - [val val] => self * right.inverse(); - [ref val] => self * right.inverse(); - [val ref] => self * right.inverse(); - [ref ref] => self * right.inverse(); + [val val] => #[allow(clippy::suspicious_arithmetic_impl)] { self * right.inverse() }; + [ref val] => #[allow(clippy::suspicious_arithmetic_impl)] { self * right.inverse() }; + [val ref] => #[allow(clippy::suspicious_arithmetic_impl)] { self * right.inverse() }; + [ref ref] => #[allow(clippy::suspicious_arithmetic_impl)] { self * right.inverse() }; ); // Isometry × UnitQuaternion @@ -479,10 +481,10 @@ isometry_from_composition_impl_all!( self: UnitQuaternion, right: Isometry>, Output = Isometry>; // FIXME: don't call inverse explicitly? - [val val] => self * right.inverse(); - [ref val] => self * right.inverse(); - [val ref] => self * right.inverse(); - [ref ref] => self * right.inverse(); + [val val] => #[allow(clippy::suspicious_arithmetic_impl)] { self * right.inverse() }; + [ref val] => #[allow(clippy::suspicious_arithmetic_impl)] { self * right.inverse() }; + [val ref] => #[allow(clippy::suspicious_arithmetic_impl)] { self * right.inverse() }; + [ref ref] => #[allow(clippy::suspicious_arithmetic_impl)] { self * right.inverse() }; ); // Translation × Rotation diff --git a/src/geometry/orthographic.rs b/src/geometry/orthographic.rs index a7d5372f..9f0eadaf 100644 --- a/src/geometry/orthographic.rs +++ b/src/geometry/orthographic.rs @@ -139,7 +139,7 @@ impl Orthographic3 { /// ``` #[inline] pub fn from_matrix_unchecked(matrix: Matrix4) -> Self { - Self { matrix: matrix } + Self { matrix } } /// Creates a new orthographic projection matrix from an aspect ratio and the vertical field of view. @@ -151,7 +151,7 @@ impl Orthographic3 { ); assert!( !relative_eq!(aspect, N::zero()), - "The apsect ratio must not be zero." + "The aspect ratio must not be zero." ); let half: N = crate::convert(0.5); diff --git a/src/geometry/perspective.rs b/src/geometry/perspective.rs index 325a450c..f550bcc6 100644 --- a/src/geometry/perspective.rs +++ b/src/geometry/perspective.rs @@ -75,7 +75,7 @@ impl Perspective3 { ); assert!( !relative_eq!(aspect, N::zero()), - "The apsect ratio must not be zero." + "The aspect ratio must not be zero." ); let matrix = Matrix4::identity(); @@ -97,7 +97,7 @@ impl Perspective3 { /// projection. #[inline] pub fn from_matrix_unchecked(matrix: Matrix4) -> Self { - Self { matrix: matrix } + Self { matrix } } /// Retrieves the inverse of the underlying homogeneous matrix. @@ -294,7 +294,7 @@ impl Arbitrary for Perspective3 { impl From> for Matrix4 { #[inline] - fn from(orth: Perspective3) -> Self { - orth.into_inner() + fn from(pers: Perspective3) -> Self { + pers.into_inner() } } diff --git a/src/geometry/point.rs b/src/geometry/point.rs index ab26589b..e0c26054 100644 --- a/src/geometry/point.rs +++ b/src/geometry/point.rs @@ -102,6 +102,45 @@ impl Point where DefaultAllocator: Allocator, { + /// Returns a point containing the result of `f` applied to each of its entries. + /// + /// # Example + /// ``` + /// # use nalgebra::{Point2, Point3}; + /// let p = Point2::new(1.0, 2.0); + /// assert_eq!(p.map(|e| e * 10.0), Point2::new(10.0, 20.0)); + /// + /// // This works in any dimension. + /// let p = Point3::new(1.1, 2.1, 3.1); + /// assert_eq!(p.map(|e| e as u32), Point3::new(1, 2, 3)); + /// ``` + #[inline] + pub fn map N2>(&self, f: F) -> Point + where + DefaultAllocator: Allocator, + { + self.coords.map(f).into() + } + + /// Replaces each component of `self` by the result of a closure `f` applied on it. + /// + /// # Example + /// ``` + /// # use nalgebra::{Point2, Point3}; + /// let mut p = Point2::new(1.0, 2.0); + /// p.apply(|e| e * 10.0); + /// assert_eq!(p, Point2::new(10.0, 20.0)); + /// + /// // This works in any dimension. + /// let mut p = Point3::new(1.0, 2.0, 3.0); + /// p.apply(|e| e * 10.0); + /// assert_eq!(p, Point3::new(10.0, 20.0, 30.0)); + /// ``` + #[inline] + pub fn apply N>(&mut self, f: F) { + self.coords.apply(f) + } + /// Converts this point into a vector in homogeneous coordinates, i.e., appends a `1` at the /// end of it. /// @@ -135,7 +174,7 @@ where #[deprecated(note = "Use Point::from(vector) instead.")] #[inline] pub fn from_coordinates(coords: VectorN) -> Self { - Self { coords: coords } + Self { coords } } /// The dimension of this point. diff --git a/src/geometry/point_coordinates.rs b/src/geometry/point_coordinates.rs index b56e120e..268c9b0d 100644 --- a/src/geometry/point_coordinates.rs +++ b/src/geometry/point_coordinates.rs @@ -1,4 +1,3 @@ -use std::mem; use std::ops::{Deref, DerefMut}; use crate::base::allocator::Allocator; @@ -22,7 +21,7 @@ macro_rules! deref_impl( #[inline] fn deref(&self) -> &Self::Target { - unsafe { mem::transmute(self) } + &*self.coords } } @@ -30,7 +29,7 @@ macro_rules! deref_impl( where DefaultAllocator: Allocator { #[inline] fn deref_mut(&mut self) -> &mut Self::Target { - unsafe { mem::transmute(self) } + &mut *self.coords } } } diff --git a/src/geometry/quaternion.rs b/src/geometry/quaternion.rs index 0c996ede..69d54967 100755 --- a/src/geometry/quaternion.rs +++ b/src/geometry/quaternion.rs @@ -389,6 +389,7 @@ where /// ``` #[inline] pub fn outer(&self, other: &Self) -> Self { + #[allow(clippy::eq_op)] (self * other - other * self).half() } @@ -1542,6 +1543,26 @@ where pub fn inverse_transform_vector(&self, v: &Vector3) -> Vector3 { self.inverse() * v } + + /// Rotate a vector by the inverse of this unit quaternion. This may be + /// cheaper than inverting the unit quaternion and transforming the + /// vector. + /// + /// # Example + /// + /// ``` + /// # #[macro_use] extern crate approx; + /// # use std::f32; + /// # use nalgebra::{UnitQuaternion, Vector3}; + /// let rot = UnitQuaternion::from_axis_angle(&Vector3::z_axis(), f32::consts::FRAC_PI_2); + /// let transformed_vector = rot.inverse_transform_unit_vector(&Vector3::x_axis()); + /// + /// assert_relative_eq!(transformed_vector, -Vector3::y_axis(), epsilon = 1.0e-6); + /// ``` + #[inline] + pub fn inverse_transform_unit_vector(&self, v: &Unit>) -> Unit> { + self.inverse() * v + } } impl Default for UnitQuaternion { diff --git a/src/geometry/quaternion_construction.rs b/src/geometry/quaternion_construction.rs index 175275d4..33e05867 100644 --- a/src/geometry/quaternion_construction.rs +++ b/src/geometry/quaternion_construction.rs @@ -374,7 +374,8 @@ where } /// The unit quaternion needed to make `a` and `b` be collinear and point toward the same - /// direction. + /// direction. Returns `None` if both `a` and `b` are collinear and point to opposite directions, as then the + /// rotation desired is not unique. /// /// # Example /// ``` @@ -491,18 +492,18 @@ where // The cosinus may be out of [-1, 1] because of inaccuracies. if cos <= -N::one() { - return None; + None } else if cos >= N::one() { - return Some(Self::identity()); + Some(Self::identity()) } else { - return Some(Self::from_axis_angle(&axis, cos.acos() * s)); + Some(Self::from_axis_angle(&axis, cos.acos() * s)) } } else if na.dot(&nb) < N::zero() { // PI // // The rotation axis is undefined but the angle not zero. This is not a // simple rotation. - return None; + None } else { // Zero Some(Self::identity()) diff --git a/src/geometry/quaternion_conversion.rs b/src/geometry/quaternion_conversion.rs index 37a8bebd..330c97e4 100644 --- a/src/geometry/quaternion_conversion.rs +++ b/src/geometry/quaternion_conversion.rs @@ -265,6 +265,15 @@ impl From> for Quaternion { } } +impl From<[N; 4]> for Quaternion { + #[inline] + fn from(coords: [N; 4]) -> Self { + Self { + coords: coords.into(), + } + } +} + impl From<[Quaternion; 2]> for Quaternion where N: From<[::Element; 2]>, diff --git a/src/geometry/quaternion_ops.rs b/src/geometry/quaternion_ops.rs index e18c1b16..d91d7fd8 100644 --- a/src/geometry/quaternion_ops.rs +++ b/src/geometry/quaternion_ops.rs @@ -121,11 +121,10 @@ quaternion_op_impl!( 'b); quaternion_op_impl!( -Add, add; -(U4, U1), (U4, U1); -self: Quaternion, rhs: Quaternion, Output = Quaternion; -Quaternion::from(self.coords + rhs.coords); -); + Add, add; + (U4, U1), (U4, U1); + self: Quaternion, rhs: Quaternion, Output = Quaternion; + Quaternion::from(self.coords + rhs.coords); ); // Quaternion - Quaternion quaternion_op_impl!( @@ -150,11 +149,10 @@ quaternion_op_impl!( 'b); quaternion_op_impl!( -Sub, sub; -(U4, U1), (U4, U1); -self: Quaternion, rhs: Quaternion, Output = Quaternion; -Quaternion::from(self.coords - rhs.coords); -); + Sub, sub; + (U4, U1), (U4, U1); + self: Quaternion, rhs: Quaternion, Output = Quaternion; + Quaternion::from(self.coords - rhs.coords); ); // Quaternion × Quaternion quaternion_op_impl!( @@ -183,11 +181,10 @@ quaternion_op_impl!( 'b); quaternion_op_impl!( -Mul, mul; -(U4, U1), (U4, U1); -self: Quaternion, rhs: Quaternion, Output = Quaternion; -&self * &rhs; -); + Mul, mul; + (U4, U1), (U4, U1); + self: Quaternion, rhs: Quaternion, Output = Quaternion; + &self * &rhs; ); // UnitQuaternion × UnitQuaternion quaternion_op_impl!( @@ -212,18 +209,17 @@ quaternion_op_impl!( 'b); quaternion_op_impl!( -Mul, mul; -(U4, U1), (U4, U1); -self: UnitQuaternion, rhs: UnitQuaternion, Output = UnitQuaternion; -&self * &rhs; -); + Mul, mul; + (U4, U1), (U4, U1); + self: UnitQuaternion, rhs: UnitQuaternion, Output = UnitQuaternion; + &self * &rhs; ); // UnitQuaternion ÷ UnitQuaternion quaternion_op_impl!( Div, div; (U4, U1), (U4, U1); self: &'a UnitQuaternion, rhs: &'b UnitQuaternion, Output = UnitQuaternion; - self * rhs.inverse(); + #[allow(clippy::suspicious_arithmetic_impl)] { self * rhs.inverse() }; 'a, 'b); quaternion_op_impl!( @@ -241,11 +237,10 @@ quaternion_op_impl!( 'b); quaternion_op_impl!( -Div, div; -(U4, U1), (U4, U1); -self: UnitQuaternion, rhs: UnitQuaternion, Output = UnitQuaternion; -&self / &rhs; -); + Div, div; + (U4, U1), (U4, U1); + self: UnitQuaternion, rhs: UnitQuaternion, Output = UnitQuaternion; + &self / &rhs; ); // UnitQuaternion × Rotation quaternion_op_impl!( @@ -274,12 +269,11 @@ quaternion_op_impl!( 'b); quaternion_op_impl!( -Mul, mul; -(U4, U1), (U3, U3); -self: UnitQuaternion, rhs: Rotation, -Output = UnitQuaternion => U3, U3; -self * UnitQuaternion::::from_rotation_matrix(&rhs); -); + Mul, mul; + (U4, U1), (U3, U3); + self: UnitQuaternion, rhs: Rotation, + Output = UnitQuaternion => U3, U3; + self * UnitQuaternion::::from_rotation_matrix(&rhs); ); // UnitQuaternion ÷ Rotation quaternion_op_impl!( @@ -308,12 +302,11 @@ quaternion_op_impl!( 'b); quaternion_op_impl!( -Div, div; -(U4, U1), (U3, U3); -self: UnitQuaternion, rhs: Rotation, -Output = UnitQuaternion => U3, U3; -self / UnitQuaternion::::from_rotation_matrix(&rhs); -); + Div, div; + (U4, U1), (U3, U3); + self: UnitQuaternion, rhs: Rotation, + Output = UnitQuaternion => U3, U3; + self / UnitQuaternion::::from_rotation_matrix(&rhs); ); // Rotation × UnitQuaternion quaternion_op_impl!( @@ -342,12 +335,11 @@ quaternion_op_impl!( 'b); quaternion_op_impl!( -Mul, mul; -(U3, U3), (U4, U1); -self: Rotation, rhs: UnitQuaternion, -Output = UnitQuaternion => U3, U3; -UnitQuaternion::::from_rotation_matrix(&self) * rhs; -); + Mul, mul; + (U3, U3), (U4, U1); + self: Rotation, rhs: UnitQuaternion, + Output = UnitQuaternion => U3, U3; + UnitQuaternion::::from_rotation_matrix(&self) * rhs; ); // Rotation ÷ UnitQuaternion quaternion_op_impl!( @@ -376,12 +368,11 @@ quaternion_op_impl!( 'b); quaternion_op_impl!( -Div, div; -(U3, U3), (U4, U1); -self: Rotation, rhs: UnitQuaternion, -Output = UnitQuaternion => U3, U3; -UnitQuaternion::::from_rotation_matrix(&self) / rhs; -); + Div, div; + (U3, U3), (U4, U1); + self: Rotation, rhs: UnitQuaternion, + Output = UnitQuaternion => U3, U3; + UnitQuaternion::::from_rotation_matrix(&self) / rhs; ); // UnitQuaternion × Vector quaternion_op_impl!( @@ -415,12 +406,11 @@ quaternion_op_impl!( 'b); quaternion_op_impl!( -Mul, mul; -(U4, U1), (U3, U1) for SB: Storage ; -self: UnitQuaternion, rhs: Vector, -Output = Vector3 => U3, U4; -&self * &rhs; -); + Mul, mul; + (U4, U1), (U3, U1) for SB: Storage ; + self: UnitQuaternion, rhs: Vector, + Output = Vector3 => U3, U4; + &self * &rhs; ); // UnitQuaternion × Point quaternion_op_impl!( @@ -448,12 +438,11 @@ quaternion_op_impl!( 'b); quaternion_op_impl!( -Mul, mul; -(U4, U1), (U3, U1); -self: UnitQuaternion, rhs: Point3, -Output = Point3 => U3, U4; -Point3::from(self * rhs.coords); -); + Mul, mul; + (U4, U1), (U3, U1); + self: UnitQuaternion, rhs: Point3, + Output = Point3 => U3, U4; + Point3::from(self * rhs.coords); ); // UnitQuaternion × Unit quaternion_op_impl!( @@ -481,12 +470,11 @@ quaternion_op_impl!( 'b); quaternion_op_impl!( -Mul, mul; -(U4, U1), (U3, U1) for SB: Storage ; -self: UnitQuaternion, rhs: Unit>, -Output = Unit> => U3, U4; -Unit::new_unchecked(self * rhs.into_inner()); -); + Mul, mul; + (U4, U1), (U3, U1) for SB: Storage ; + self: UnitQuaternion, rhs: Unit>, + Output = Unit> => U3, U4; + Unit::new_unchecked(self * rhs.into_inner()); ); macro_rules! scalar_op_impl( ($($Op: ident, $op: ident, $OpAssign: ident, $op_assign: ident);* $(;)*) => {$( diff --git a/src/geometry/rotation.rs b/src/geometry/rotation.rs index c6cd3d86..c30a3a97 100755 --- a/src/geometry/rotation.rs +++ b/src/geometry/rotation.rs @@ -19,7 +19,7 @@ use simba::simd::SimdRealField; use crate::base::allocator::Allocator; use crate::base::dimension::{DimName, DimNameAdd, DimNameSum, U1}; -use crate::base::{DefaultAllocator, MatrixN, Scalar, VectorN}; +use crate::base::{DefaultAllocator, MatrixN, Scalar, Unit, VectorN}; use crate::geometry::Point; /// A rotation matrix. @@ -255,7 +255,7 @@ where "Unable to create a rotation from a non-square matrix." ); - Self { matrix: matrix } + Self { matrix } } /// Transposes `self`. @@ -441,6 +441,25 @@ where pub fn inverse_transform_vector(&self, v: &VectorN) -> VectorN { self.matrix().tr_mul(v) } + + /// Rotate the given vector by the inverse of this rotation. This may be + /// cheaper than inverting the rotation and then transforming the given + /// vector. + /// + /// # Example + /// ``` + /// # #[macro_use] extern crate approx; + /// # use std::f32; + /// # use nalgebra::{Rotation2, Rotation3, UnitQuaternion, Vector3}; + /// let rot = Rotation3::new(Vector3::z() * f32::consts::FRAC_PI_2); + /// let transformed_vector = rot.inverse_transform_unit_vector(&Vector3::x_axis()); + /// + /// assert_relative_eq!(transformed_vector, -Vector3::y_axis(), epsilon = 1.0e-6); + /// ``` + #[inline] + pub fn inverse_transform_unit_vector(&self, v: &Unit>) -> Unit> { + Unit::new_unchecked(self.inverse_transform_vector(&**v)) + } } impl Eq for Rotation where DefaultAllocator: Allocator {} diff --git a/src/geometry/rotation_ops.rs b/src/geometry/rotation_ops.rs index b3330877..9b90c856 100644 --- a/src/geometry/rotation_ops.rs +++ b/src/geometry/rotation_ops.rs @@ -59,10 +59,10 @@ md_impl_all!( Div, div; (D, D), (D, D) for D: DimName; self: Rotation, right: Rotation, Output = Rotation; - [val val] => self * right.inverse(); - [ref val] => self * right.inverse(); - [val ref] => self * right.inverse(); - [ref ref] => self * right.inverse(); + [val val] => #[allow(clippy::suspicious_arithmetic_impl)] { self * right.inverse() }; + [ref val] => #[allow(clippy::suspicious_arithmetic_impl)] { self * right.inverse() }; + [val ref] => #[allow(clippy::suspicious_arithmetic_impl)] { self * right.inverse() }; + [ref ref] => #[allow(clippy::suspicious_arithmetic_impl)] { self * right.inverse() }; ); // Rotation × Matrix @@ -98,10 +98,10 @@ md_impl_all!( where DefaultAllocator: Allocator where ShapeConstraint: AreMultipliable; self: Matrix, right: Rotation, Output = MatrixMN; - [val val] => self * right.inverse(); - [ref val] => self * right.inverse(); - [val ref] => self * right.inverse(); - [ref ref] => self * right.inverse(); + [val val] => #[allow(clippy::suspicious_arithmetic_impl)] { self * right.inverse() }; + [ref val] => #[allow(clippy::suspicious_arithmetic_impl)] { self * right.inverse() }; + [val ref] => #[allow(clippy::suspicious_arithmetic_impl)] { self * right.inverse() }; + [ref ref] => #[allow(clippy::suspicious_arithmetic_impl)] { self * right.inverse() }; ); // Rotation × Point diff --git a/src/geometry/rotation_specialization.rs b/src/geometry/rotation_specialization.rs index 1ee4a9df..085cb2d8 100644 --- a/src/geometry/rotation_specialization.rs +++ b/src/geometry/rotation_specialization.rs @@ -236,6 +236,31 @@ impl Rotation2 { pub fn scaled_axis(&self) -> VectorN { Vector1::new(self.angle()) } + + /// Spherical linear interpolation between two rotation matrices. + /// + /// # Examples: + /// + /// ``` + /// # #[macro_use] extern crate approx; + /// # use nalgebra::geometry::Rotation2; + /// + /// let rot1 = Rotation2::new(std::f32::consts::FRAC_PI_4); + /// let rot2 = Rotation2::new(-std::f32::consts::PI); + /// + /// let rot = rot1.slerp(&rot2, 1.0 / 3.0); + /// + /// assert_relative_eq!(rot.angle(), std::f32::consts::FRAC_PI_2); + /// ``` + #[inline] + pub fn slerp(&self, other: &Self, t: N) -> Self + where + N::Element: SimdRealField, + { + let c1 = UnitComplex::from(*self); + let c2 = UnitComplex::from(*other); + c1.slerp(&c2, t).into() + } } impl Distribution> for Standard @@ -862,6 +887,53 @@ where Self::identity() } } + + /// Spherical linear interpolation between two rotation matrices. + /// + /// Panics if the angle between both rotations is 180 degrees (in which case the interpolation + /// is not well-defined). Use `.try_slerp` instead to avoid the panic. + /// + /// # Examples: + /// + /// ``` + /// # use nalgebra::geometry::Rotation3; + /// + /// let q1 = Rotation3::from_euler_angles(std::f32::consts::FRAC_PI_4, 0.0, 0.0); + /// let q2 = Rotation3::from_euler_angles(-std::f32::consts::PI, 0.0, 0.0); + /// + /// let q = q1.slerp(&q2, 1.0 / 3.0); + /// + /// assert_eq!(q.euler_angles(), (std::f32::consts::FRAC_PI_2, 0.0, 0.0)); + /// ``` + #[inline] + pub fn slerp(&self, other: &Self, t: N) -> Self + where + N: RealField, + { + let q1 = UnitQuaternion::from(*self); + let q2 = UnitQuaternion::from(*other); + q1.slerp(&q2, t).into() + } + + /// Computes the spherical linear interpolation between two rotation matrices or returns `None` + /// if both rotations are approximately 180 degrees apart (in which case the interpolation is + /// not well-defined). + /// + /// # Arguments + /// * `self`: the first rotation to interpolate from. + /// * `other`: the second rotation to interpolate toward. + /// * `t`: the interpolation parameter. Should be between 0 and 1. + /// * `epsilon`: the value below which the sinus of the angle separating both rotations + /// must be to return `None`. + #[inline] + pub fn try_slerp(&self, other: &Self, t: N, epsilon: N) -> Option + where + N: RealField, + { + let q1 = Rotation3::from(*self); + let q2 = Rotation3::from(*other); + q1.try_slerp(&q2, t, epsilon).map(|q| q.into()) + } } impl Distribution> for Standard diff --git a/src/geometry/similarity_ops.rs b/src/geometry/similarity_ops.rs index 9043a667..d011871d 100644 --- a/src/geometry/similarity_ops.rs +++ b/src/geometry/similarity_ops.rs @@ -158,15 +158,14 @@ similarity_binop_impl_all!( similarity_binop_impl_all!( Div, div; self: Similarity, rhs: Similarity, Output = Similarity; - [val val] => self * rhs.inverse(); - [ref val] => self * rhs.inverse(); - [val ref] => self * rhs.inverse(); - [ref ref] => self * rhs.inverse(); + [val val] => #[allow(clippy::suspicious_arithmetic_impl)] { self * rhs.inverse() }; + [ref val] => #[allow(clippy::suspicious_arithmetic_impl)] { self * rhs.inverse() }; + [val ref] => #[allow(clippy::suspicious_arithmetic_impl)] { self * rhs.inverse() }; + [ref ref] => #[allow(clippy::suspicious_arithmetic_impl)] { self * rhs.inverse() }; ); // Similarity ×= Translation similarity_binop_assign_impl_all!( - MulAssign, mul_assign; self: Similarity, rhs: Translation; [val] => *self *= &rhs; @@ -281,6 +280,7 @@ similarity_binop_impl_all!( [ref ref] => { let shift = self.isometry.rotation.transform_vector(&rhs.translation.vector) * self.scaling(); Similarity::from_parts( + #[allow(clippy::suspicious_arithmetic_impl)] Translation::from(&self.isometry.translation.vector + shift), self.isometry.rotation.clone() * rhs.rotation.clone(), self.scaling()) @@ -290,10 +290,10 @@ similarity_binop_impl_all!( similarity_binop_impl_all!( Div, div; self: Similarity, rhs: Isometry, Output = Similarity; - [val val] => self * rhs.inverse(); - [ref val] => self * rhs.inverse(); - [val ref] => self * rhs.inverse(); - [ref ref] => self * rhs.inverse(); + [val val] => #[allow(clippy::suspicious_arithmetic_impl)] { self * rhs.inverse() }; + [ref val] => #[allow(clippy::suspicious_arithmetic_impl)] { self * rhs.inverse() }; + [val ref] => #[allow(clippy::suspicious_arithmetic_impl)] { self * rhs.inverse() }; + [ref ref] => #[allow(clippy::suspicious_arithmetic_impl)] { self * rhs.inverse() }; ); // Isometry × Similarity @@ -322,10 +322,10 @@ similarity_binop_impl_all!( similarity_binop_impl_all!( Div, div; self: Isometry, rhs: Similarity, Output = Similarity; - [val val] => self * rhs.inverse(); - [ref val] => self * rhs.inverse(); - [val ref] => self * rhs.inverse(); - [ref ref] => self * rhs.inverse(); + [val val] => #[allow(clippy::suspicious_arithmetic_impl)] { self * rhs.inverse() }; + [ref val] => #[allow(clippy::suspicious_arithmetic_impl)] { self * rhs.inverse() }; + [val ref] => #[allow(clippy::suspicious_arithmetic_impl)] { self * rhs.inverse() }; + [ref ref] => #[allow(clippy::suspicious_arithmetic_impl)] { self * rhs.inverse() }; ); // Similarity × Point @@ -364,6 +364,7 @@ similarity_binop_impl_all!( [ref ref] => { let shift = self.isometry.rotation.transform_vector(&right.vector) * self.scaling(); Similarity::from_parts( + #[allow(clippy::suspicious_arithmetic_impl)] Translation::from(&self.isometry.translation.vector + shift), self.isometry.rotation.clone(), self.scaling()) @@ -495,10 +496,10 @@ similarity_from_composition_impl_all!( self: Rotation, right: Similarity>, Output = Similarity>; // FIXME: don't call inverse explicitly? - [val val] => self * right.inverse(); - [ref val] => self * right.inverse(); - [val ref] => self * right.inverse(); - [ref ref] => self * right.inverse(); + [val val] => #[allow(clippy::suspicious_arithmetic_impl)] { self * right.inverse() }; + [ref val] => #[allow(clippy::suspicious_arithmetic_impl)] { self * right.inverse() }; + [val ref] => #[allow(clippy::suspicious_arithmetic_impl)] { self * right.inverse() }; + [ref ref] => #[allow(clippy::suspicious_arithmetic_impl)] { self * right.inverse() }; ); // Similarity × UnitQuaternion @@ -556,10 +557,10 @@ similarity_from_composition_impl_all!( self: UnitQuaternion, right: Similarity>, Output = Similarity>; // FIXME: don't call inverse explicitly? - [val val] => self * right.inverse(); - [ref val] => self * right.inverse(); - [val ref] => self * right.inverse(); - [ref ref] => self * right.inverse(); + [val val] => #[allow(clippy::suspicious_arithmetic_impl)] { self * right.inverse() }; + [ref val] => #[allow(clippy::suspicious_arithmetic_impl)] { self * right.inverse() }; + [val ref] => #[allow(clippy::suspicious_arithmetic_impl)] { self * right.inverse() }; + [ref ref] => #[allow(clippy::suspicious_arithmetic_impl)] { self * right.inverse() }; ); // Similarity × UnitComplex diff --git a/src/geometry/transform.rs b/src/geometry/transform.rs index 42586509..80b72771 100755 --- a/src/geometry/transform.rs +++ b/src/geometry/transform.rs @@ -245,7 +245,7 @@ where #[inline] pub fn from_matrix_unchecked(matrix: MatrixN>) -> Self { Transform { - matrix: matrix, + matrix, _phantom: PhantomData, } } diff --git a/src/geometry/transform_ops.rs b/src/geometry/transform_ops.rs index 5773782c..6986a90a 100644 --- a/src/geometry/transform_ops.rs +++ b/src/geometry/transform_ops.rs @@ -143,6 +143,7 @@ md_impl_all!( if C::has_normalizer() { let normalizer = self.matrix().fixed_slice::(D::dim(), 0); + #[allow(clippy::suspicious_arithmetic_impl)] let n = normalizer.tr_dot(&rhs.coords) + unsafe { *self.matrix().get_unchecked((D::dim(), D::dim())) }; if !n.is_zero() { @@ -293,10 +294,10 @@ md_impl_all!( Div, div where N: RealField; (DimNameSum, DimNameSum), (DimNameSum, DimNameSum) for D: DimNameAdd, CA: TCategoryMul, CB: SubTCategoryOf; self: Transform, rhs: Transform, Output = Transform; - [val val] => self * rhs.inverse(); - [ref val] => self * rhs.inverse(); - [val ref] => self * rhs.clone().inverse(); - [ref ref] => self * rhs.clone().inverse(); + [val val] => #[allow(clippy::suspicious_arithmetic_impl)] { self * rhs.inverse() }; + [ref val] => #[allow(clippy::suspicious_arithmetic_impl)] { self * rhs.inverse() }; + [val ref] => #[allow(clippy::suspicious_arithmetic_impl)] { self * rhs.clone().inverse() }; + [ref ref] => #[allow(clippy::suspicious_arithmetic_impl)] { self * rhs.clone().inverse() }; ); // Transform ÷ Rotation @@ -304,10 +305,10 @@ md_impl_all!( Div, div where N: RealField; (DimNameSum, DimNameSum), (D, D) for D: DimNameAdd, C: TCategoryMul; self: Transform, rhs: Rotation, Output = Transform; - [val val] => self * rhs.inverse(); - [ref val] => self * rhs.inverse(); - [val ref] => self * rhs.inverse(); - [ref ref] => self * rhs.inverse(); + [val val] => #[allow(clippy::suspicious_arithmetic_impl)] { self * rhs.inverse() }; + [ref val] => #[allow(clippy::suspicious_arithmetic_impl)] { self * rhs.inverse() }; + [val ref] => #[allow(clippy::suspicious_arithmetic_impl)] { self * rhs.inverse() }; + [ref ref] => #[allow(clippy::suspicious_arithmetic_impl)] { self * rhs.inverse() }; ); // Rotation ÷ Transform @@ -315,10 +316,10 @@ md_impl_all!( Div, div where N: RealField; (D, D), (DimNameSum, DimNameSum) for D: DimNameAdd, C: TCategoryMul; self: Rotation, rhs: Transform, Output = Transform; - [val val] => self.inverse() * rhs; - [ref val] => self.inverse() * rhs; - [val ref] => self.inverse() * rhs; - [ref ref] => self.inverse() * rhs; + [val val] => #[allow(clippy::suspicious_arithmetic_impl)] { self.inverse() * rhs }; + [ref val] => #[allow(clippy::suspicious_arithmetic_impl)] { self.inverse() * rhs }; + [val ref] => #[allow(clippy::suspicious_arithmetic_impl)] { self.inverse() * rhs }; + [ref ref] => #[allow(clippy::suspicious_arithmetic_impl)] { self.inverse() * rhs }; ); // Transform ÷ UnitQuaternion @@ -326,10 +327,10 @@ md_impl_all!( Div, div where N: RealField; (U4, U4), (U4, U1) for C: TCategoryMul; self: Transform, rhs: UnitQuaternion, Output = Transform; - [val val] => self * rhs.inverse(); - [ref val] => self * rhs.inverse(); - [val ref] => self * rhs.inverse(); - [ref ref] => self * rhs.inverse(); + [val val] => #[allow(clippy::suspicious_arithmetic_impl)] { self * rhs.inverse() }; + [ref val] => #[allow(clippy::suspicious_arithmetic_impl)] { self * rhs.inverse() }; + [val ref] => #[allow(clippy::suspicious_arithmetic_impl)] { self * rhs.inverse() }; + [ref ref] => #[allow(clippy::suspicious_arithmetic_impl)] { self * rhs.inverse() }; ); // UnitQuaternion ÷ Transform @@ -337,10 +338,10 @@ md_impl_all!( Div, div where N: RealField; (U4, U1), (U4, U4) for C: TCategoryMul; self: UnitQuaternion, rhs: Transform, Output = Transform; - [val val] => self.inverse() * rhs; - [ref val] => self.inverse() * rhs; - [val ref] => self.inverse() * rhs; - [ref ref] => self.inverse() * rhs; + [val val] => #[allow(clippy::suspicious_arithmetic_impl)] { self.inverse() * rhs }; + [ref val] => #[allow(clippy::suspicious_arithmetic_impl)] { self.inverse() * rhs }; + [val ref] => #[allow(clippy::suspicious_arithmetic_impl)] { self.inverse() * rhs }; + [ref ref] => #[allow(clippy::suspicious_arithmetic_impl)] { self.inverse() * rhs }; ); // // Transform ÷ Isometry @@ -402,10 +403,10 @@ md_impl_all!( Div, div where N: RealField; (DimNameSum, DimNameSum), (D, U1) for D: DimNameAdd, C: TCategoryMul; self: Transform, rhs: Translation, Output = Transform; - [val val] => self * rhs.inverse(); - [ref val] => self * rhs.inverse(); - [val ref] => self * rhs.inverse(); - [ref ref] => self * rhs.inverse(); + [val val] => #[allow(clippy::suspicious_arithmetic_impl)] { self * rhs.inverse() }; + [ref val] => #[allow(clippy::suspicious_arithmetic_impl)] { self * rhs.inverse() }; + [val ref] => #[allow(clippy::suspicious_arithmetic_impl)] { self * rhs.inverse() }; + [ref ref] => #[allow(clippy::suspicious_arithmetic_impl)] { self * rhs.inverse() }; ); // Translation ÷ Transform @@ -414,10 +415,10 @@ md_impl_all!( (D, U1), (DimNameSum, DimNameSum) for D: DimNameAdd, C: TCategoryMul; self: Translation, rhs: Transform, Output = Transform; - [val val] => self.inverse() * rhs; - [ref val] => self.inverse() * rhs; - [val ref] => self.inverse() * rhs; - [ref ref] => self.inverse() * rhs; + [val val] => #[allow(clippy::suspicious_arithmetic_impl)] { self.inverse() * rhs }; + [ref val] => #[allow(clippy::suspicious_arithmetic_impl)] { self.inverse() * rhs }; + [val ref] => #[allow(clippy::suspicious_arithmetic_impl)] { self.inverse() * rhs }; + [ref ref] => #[allow(clippy::suspicious_arithmetic_impl)] { self.inverse() * rhs }; ); // Transform ×= Transform diff --git a/src/geometry/translation.rs b/src/geometry/translation.rs index 88d52ce4..26bbfd12 100755 --- a/src/geometry/translation.rs +++ b/src/geometry/translation.rs @@ -119,7 +119,7 @@ where #[inline] #[deprecated(note = "Use `::from` instead.")] pub fn from_vector(vector: VectorN) -> Translation { - Translation { vector: vector } + Translation { vector } } /// Inverts `self`. diff --git a/src/geometry/translation_ops.rs b/src/geometry/translation_ops.rs index c00b25bc..8d543c20 100644 --- a/src/geometry/translation_ops.rs +++ b/src/geometry/translation_ops.rs @@ -13,44 +13,50 @@ use crate::geometry::{Point, Translation}; add_sub_impl!(Mul, mul, ClosedAdd; (D, U1), (D, U1) -> (D) for D: DimName; self: &'a Translation, right: &'b Translation, Output = Translation; - Translation::from(&self.vector + &right.vector); 'a, 'b); + #[allow(clippy::suspicious_arithmetic_impl)] { Translation::from(&self.vector + &right.vector) }; + 'a, 'b); add_sub_impl!(Mul, mul, ClosedAdd; (D, U1), (D, U1) -> (D) for D: DimName; self: &'a Translation, right: Translation, Output = Translation; - Translation::from(&self.vector + right.vector); 'a); + #[allow(clippy::suspicious_arithmetic_impl)] { Translation::from(&self.vector + right.vector) }; + 'a); add_sub_impl!(Mul, mul, ClosedAdd; (D, U1), (D, U1) -> (D) for D: DimName; self: Translation, right: &'b Translation, Output = Translation; - Translation::from(self.vector + &right.vector); 'b); + #[allow(clippy::suspicious_arithmetic_impl)] { Translation::from(self.vector + &right.vector) }; + 'b); add_sub_impl!(Mul, mul, ClosedAdd; (D, U1), (D, U1) -> (D) for D: DimName; self: Translation, right: Translation, Output = Translation; - Translation::from(self.vector + right.vector); ); + #[allow(clippy::suspicious_arithmetic_impl)] { Translation::from(self.vector + right.vector) }; ); // Translation ÷ Translation // FIXME: instead of calling inverse explicitly, could we just add a `mul_tr` or `mul_inv` method? add_sub_impl!(Div, div, ClosedSub; (D, U1), (D, U1) -> (D) for D: DimName; self: &'a Translation, right: &'b Translation, Output = Translation; - Translation::from(&self.vector - &right.vector); 'a, 'b); + #[allow(clippy::suspicious_arithmetic_impl)] { Translation::from(&self.vector - &right.vector) }; + 'a, 'b); add_sub_impl!(Div, div, ClosedSub; (D, U1), (D, U1) -> (D) for D: DimName; self: &'a Translation, right: Translation, Output = Translation; - Translation::from(&self.vector - right.vector); 'a); + #[allow(clippy::suspicious_arithmetic_impl)] { Translation::from(&self.vector - right.vector) }; + 'a); add_sub_impl!(Div, div, ClosedSub; (D, U1), (D, U1) -> (D) for D: DimName; self: Translation, right: &'b Translation, Output = Translation; - Translation::from(self.vector - &right.vector); 'b); + #[allow(clippy::suspicious_arithmetic_impl)] { Translation::from(self.vector - &right.vector) }; + 'b); add_sub_impl!(Div, div, ClosedSub; (D, U1), (D, U1) -> (D) for D: DimName; self: Translation, right: Translation, Output = Translation; - Translation::from(self.vector - right.vector); ); + #[allow(clippy::suspicious_arithmetic_impl)] { Translation::from(self.vector - right.vector) }; ); // Translation × Point // FIXME: we don't handle properly non-zero origins here. Do we want this to be the intended @@ -58,107 +64,45 @@ add_sub_impl!(Div, div, ClosedSub; add_sub_impl!(Mul, mul, ClosedAdd; (D, U1), (D, U1) -> (D) for D: DimName; self: &'a Translation, right: &'b Point, Output = Point; - right + &self.vector; 'a, 'b); + #[allow(clippy::suspicious_arithmetic_impl)] { right + &self.vector }; + 'a, 'b); add_sub_impl!(Mul, mul, ClosedAdd; (D, U1), (D, U1) -> (D) for D: DimName; self: &'a Translation, right: Point, Output = Point; - right + &self.vector; 'a); + #[allow(clippy::suspicious_arithmetic_impl)] { right + &self.vector }; + 'a); add_sub_impl!(Mul, mul, ClosedAdd; (D, U1), (D, U1) -> (D) for D: DimName; self: Translation, right: &'b Point, Output = Point; - right + self.vector; 'b); + #[allow(clippy::suspicious_arithmetic_impl)] { right + self.vector }; + 'b); add_sub_impl!(Mul, mul, ClosedAdd; (D, U1), (D, U1) -> (D) for D: DimName; self: Translation, right: Point, Output = Point; - right + self.vector; ); + #[allow(clippy::suspicious_arithmetic_impl)] { right + self.vector }; ); // Translation *= Translation add_sub_assign_impl!(MulAssign, mul_assign, ClosedAdd; (D, U1), (D, U1) for D: DimName; self: Translation, right: &'b Translation; - self.vector += &right.vector; 'b); + #[allow(clippy::suspicious_arithmetic_impl)] { self.vector += &right.vector }; + 'b); add_sub_assign_impl!(MulAssign, mul_assign, ClosedAdd; (D, U1), (D, U1) for D: DimName; self: Translation, right: Translation; - self.vector += right.vector; ); + #[allow(clippy::suspicious_arithmetic_impl)] { self.vector += right.vector }; ); add_sub_assign_impl!(DivAssign, div_assign, ClosedSub; (D, U1), (D, U1) for D: DimName; self: Translation, right: &'b Translation; - self.vector -= &right.vector; 'b); + #[allow(clippy::suspicious_arithmetic_impl)] { self.vector -= &right.vector }; + 'b); add_sub_assign_impl!(DivAssign, div_assign, ClosedSub; (D, U1), (D, U1) for D: DimName; self: Translation, right: Translation; - self.vector -= right.vector; ); - -/* -// Translation × Matrix -add_sub_impl!(Mul, mul; - (D1, D1), (R2, C2) for D1, R2, C2; - self: &'a Translation, right: &'b Matrix, Output = MatrixMN; - self.vector() * right; 'a, 'b); - -add_sub_impl!(Mul, mul; - (D1, D1), (R2, C2) for D1, R2, C2; - self: &'a Translation, right: Matrix, Output = MatrixMN; - self.vector() * right; 'a); - -add_sub_impl!(Mul, mul; - (D1, D1), (R2, C2) for D1, R2, C2; - self: Translation, right: &'b Matrix, Output = MatrixMN; - self.unwrap() * right; 'b); - -add_sub_impl!(Mul, mul; - (D1, D1), (R2, C2) for D1, R2, C2; - self: Translation, right: Matrix, Output = MatrixMN; - self.unwrap() * right; ); - -// Matrix × Translation -add_sub_impl!(Mul, mul; - (R1, C1), (D2, D2) for R1, C1, D2; - self: &'a Matrix, right: &'b Translation, Output = MatrixMN; - self * right.vector(); 'a, 'b); - -add_sub_impl!(Mul, mul; - (R1, C1), (D2, D2) for R1, C1, D2; - self: &'a Matrix, right: Translation, Output = MatrixMN; - self * right.unwrap(); 'a); - -add_sub_impl!(Mul, mul; - (R1, C1), (D2, D2) for R1, C1, D2; - self: Matrix, right: &'b Translation, Output = MatrixMN; - self * right.vector(); 'b); - - -add_sub_impl!(Mul, mul; - (R1, C1), (D2, D2) for R1, C1, D2; - self: Matrix, right: Translation, Output = MatrixMN; - self * right.unwrap(); ); - -// Matrix *= Translation -md_assign_impl!(MulAssign, mul_assign; - (R1, C1), (C1, C1) for R1, C1; - self: Matrix, right: &'b Translation; - self.mul_assign(right.vector()); 'b); - -md_assign_impl!(MulAssign, mul_assign; - (R1, C1), (C1, C1) for R1, C1; - self: Matrix, right: Translation; - self.mul_assign(right.unwrap()); ); - - -md_assign_impl!(DivAssign, div_assign; - (R1, C1), (C1, C1) for R1, C1; - self: Matrix, right: &'b Translation; - self.mul_assign(right.inverse().vector()); 'b); - -md_assign_impl!(DivAssign, div_assign; - (R1, C1), (C1, C1) for R1, C1; - self: Matrix, right: Translation; - self.mul_assign(right.inverse().unwrap()); ); -*/ + #[allow(clippy::suspicious_arithmetic_impl)] { self.vector -= right.vector }; ); diff --git a/src/geometry/unit_complex.rs b/src/geometry/unit_complex.rs index a8ac5bcd..a31f3ddf 100755 --- a/src/geometry/unit_complex.rs +++ b/src/geometry/unit_complex.rs @@ -360,6 +360,43 @@ where pub fn inverse_transform_vector(&self, v: &Vector2) -> Vector2 { self.inverse() * v } + + /// Rotate the given vector by the inverse of this unit complex number. + /// + /// # Example + /// ``` + /// # #[macro_use] extern crate approx; + /// # use nalgebra::{UnitComplex, Vector2}; + /// # use std::f32; + /// let rot = UnitComplex::new(f32::consts::FRAC_PI_2); + /// let transformed_vector = rot.inverse_transform_unit_vector(&Vector2::x_axis()); + /// assert_relative_eq!(transformed_vector, -Vector2::y_axis(), epsilon = 1.0e-6); + /// ``` + #[inline] + pub fn inverse_transform_unit_vector(&self, v: &Unit>) -> Unit> { + self.inverse() * v + } + + /// Spherical linear interpolation between two rotations represented as unit complex numbers. + /// + /// # Examples: + /// + /// ``` + /// # #[macro_use] extern crate approx; + /// # use nalgebra::geometry::UnitComplex; + /// + /// let rot1 = UnitComplex::new(std::f32::consts::FRAC_PI_4); + /// let rot2 = UnitComplex::new(-std::f32::consts::PI); + /// + /// let rot = rot1.slerp(&rot2, 1.0 / 3.0); + /// + /// assert_relative_eq!(rot.angle(), std::f32::consts::FRAC_PI_2); + /// ``` + + #[inline] + pub fn slerp(&self, other: &Self, t: N) -> Self { + Self::new(self.angle() * (N::one() - t) + other.angle() * t) + } } impl fmt::Display for UnitComplex { diff --git a/src/geometry/unit_complex_ops.rs b/src/geometry/unit_complex_ops.rs index f7af7088..af1a32f1 100644 --- a/src/geometry/unit_complex_ops.rs +++ b/src/geometry/unit_complex_ops.rs @@ -96,6 +96,7 @@ where #[inline] fn div(self, rhs: Self) -> Self::Output { + #[allow(clippy::suspicious_arithmetic_impl)] Unit::new_unchecked(self.into_inner() * rhs.conjugate().into_inner()) } } @@ -108,6 +109,7 @@ where #[inline] fn div(self, rhs: UnitComplex) -> Self::Output { + #[allow(clippy::suspicious_arithmetic_impl)] Unit::new_unchecked(self.complex() * rhs.conjugate().into_inner()) } } @@ -120,6 +122,7 @@ where #[inline] fn div(self, rhs: &'b UnitComplex) -> Self::Output { + #[allow(clippy::suspicious_arithmetic_impl)] Unit::new_unchecked(self.into_inner() * rhs.conjugate().into_inner()) } } @@ -132,6 +135,7 @@ where #[inline] fn div(self, rhs: &'b UnitComplex) -> Self::Output { + #[allow(clippy::suspicious_arithmetic_impl)] Unit::new_unchecked(self.complex() * rhs.conjugate().into_inner()) } } @@ -206,7 +210,7 @@ complex_op_impl_all!( [val val] => &self / &rhs; [ref val] => self / &rhs; [val ref] => &self / rhs; - [ref ref] => self * UnitComplex::from_rotation_matrix(rhs).inverse(); + [ref ref] => #[allow(clippy::suspicious_arithmetic_impl)] { self * UnitComplex::from_rotation_matrix(rhs).inverse() }; ); // Rotation × UnitComplex @@ -228,7 +232,7 @@ complex_op_impl_all!( [val val] => &self / &rhs; [ref val] => self / &rhs; [val ref] => &self / rhs; - [ref ref] => UnitComplex::from_rotation_matrix(self) * rhs.inverse(); + [ref ref] => #[allow(clippy::suspicious_arithmetic_impl)] { UnitComplex::from_rotation_matrix(self) * rhs.inverse() }; ); // UnitComplex × Point diff --git a/src/lib.rs b/src/lib.rs index 9db176bf..ac954f06 100644 --- a/src/lib.rs +++ b/src/lib.rs @@ -15,7 +15,7 @@ Simply add the following to your `Cargo.toml` file: ```.ignore [dependencies] -nalgebra = "0.18" +nalgebra = "0.21" ``` @@ -90,11 +90,9 @@ an optimized set of tools for computer graphics and physics. Those features incl #[cfg(feature = "arbitrary")] extern crate quickcheck; -#[cfg(feature = "serde")] -extern crate serde; -#[cfg(feature = "serde")] +#[cfg(feature = "serde-serialize")] #[macro_use] -extern crate serde_derive; +extern crate serde; #[cfg(feature = "abomonation-serialize")] extern crate abomonation; diff --git a/src/linalg/cholesky.rs b/src/linalg/cholesky.rs index d9f33f27..17187fdf 100644 --- a/src/linalg/cholesky.rs +++ b/src/linalg/cholesky.rs @@ -8,7 +8,7 @@ use simba::simd::SimdComplexField; use crate::allocator::Allocator; use crate::base::{DefaultAllocator, Matrix, MatrixMN, MatrixN, SquareMatrix, Vector}; use crate::constraint::{SameNumberOfRows, ShapeConstraint}; -use crate::dimension::{Dim, DimAdd, DimDiff, DimSub, DimSum, Dynamic, U1}; +use crate::dimension::{Dim, DimAdd, DimDiff, DimSub, DimSum, U1}; use crate::storage::{Storage, StorageMut}; /// The Cholesky decomposition of a symmetric-definite-positive matrix. @@ -364,7 +364,7 @@ where } } -impl, S: Storage> SquareMatrix +impl> SquareMatrix where DefaultAllocator: Allocator, { diff --git a/src/linalg/exp.rs b/src/linalg/exp.rs index c11fd757..6c43568d 100644 --- a/src/linalg/exp.rs +++ b/src/linalg/exp.rs @@ -10,10 +10,12 @@ use crate::{ convert, try_convert, ComplexField, MatrixN, RealField, }; +use crate::num::Zero; + // https://github.com/scipy/scipy/blob/c1372d8aa90a73d8a52f135529293ff4edb98fc8/scipy/sparse/linalg/matfuncs.py struct ExpmPadeHelper where - N: RealField, + N: ComplexField, D: DimMin, DefaultAllocator: Allocator + Allocator<(usize, usize), DimMinimum>, { @@ -27,20 +29,20 @@ where a8: Option>, a10: Option>, - d4_exact: Option, - d6_exact: Option, - d8_exact: Option, - d10_exact: Option, + d4_exact: Option, + d6_exact: Option, + d8_exact: Option, + d10_exact: Option, - d4_approx: Option, - d6_approx: Option, - d8_approx: Option, - d10_approx: Option, + d4_approx: Option, + d6_approx: Option, + d8_approx: Option, + d10_approx: Option, } impl ExpmPadeHelper where - N: RealField, + N: ComplexField, D: DimMin, DefaultAllocator: Allocator + Allocator<(usize, usize), DimMinimum>, { @@ -110,7 +112,7 @@ where } } - fn d4_tight(&mut self) -> N { + fn d4_tight(&mut self) -> N::RealField { if self.d4_exact.is_none() { self.calc_a4(); self.d4_exact = Some(one_norm(self.a4.as_ref().unwrap()).powf(convert(0.25))); @@ -118,7 +120,7 @@ where self.d4_exact.unwrap() } - fn d6_tight(&mut self) -> N { + fn d6_tight(&mut self) -> N::RealField { if self.d6_exact.is_none() { self.calc_a6(); self.d6_exact = Some(one_norm(self.a6.as_ref().unwrap()).powf(convert(1.0 / 6.0))); @@ -126,7 +128,7 @@ where self.d6_exact.unwrap() } - fn d8_tight(&mut self) -> N { + fn d8_tight(&mut self) -> N::RealField { if self.d8_exact.is_none() { self.calc_a8(); self.d8_exact = Some(one_norm(self.a8.as_ref().unwrap()).powf(convert(1.0 / 8.0))); @@ -134,7 +136,7 @@ where self.d8_exact.unwrap() } - fn d10_tight(&mut self) -> N { + fn d10_tight(&mut self) -> N::RealField { if self.d10_exact.is_none() { self.calc_a10(); self.d10_exact = Some(one_norm(self.a10.as_ref().unwrap()).powf(convert(1.0 / 10.0))); @@ -142,7 +144,7 @@ where self.d10_exact.unwrap() } - fn d4_loose(&mut self) -> N { + fn d4_loose(&mut self) -> N::RealField { if self.use_exact_norm { return self.d4_tight(); } @@ -159,7 +161,7 @@ where self.d4_approx.unwrap() } - fn d6_loose(&mut self) -> N { + fn d6_loose(&mut self) -> N::RealField { if self.use_exact_norm { return self.d6_tight(); } @@ -176,7 +178,7 @@ where self.d6_approx.unwrap() } - fn d8_loose(&mut self) -> N { + fn d8_loose(&mut self) -> N::RealField { if self.use_exact_norm { return self.d8_tight(); } @@ -193,7 +195,7 @@ where self.d8_approx.unwrap() } - fn d10_loose(&mut self) -> N { + fn d10_loose(&mut self) -> N::RealField { if self.use_exact_norm { return self.d10_tight(); } @@ -359,15 +361,20 @@ where fn ell(a: &MatrixN, m: u64) -> u64 where - N: RealField, + N: ComplexField, D: Dim, - DefaultAllocator: Allocator + Allocator, + DefaultAllocator: Allocator + + Allocator + + Allocator + + Allocator, { // 2m choose m = (2m)!/(m! * (2m-m)!) - let a_abs_onenorm = onenorm_matrix_power_nonm(&a.abs(), 2 * m + 1); + let a_abs = a.map(|x| x.abs()); - if a_abs_onenorm == N::zero() { + let a_abs_onenorm = onenorm_matrix_power_nonm(&a_abs, 2 * m + 1); + + if a_abs_onenorm == ::RealField::zero() { return 0; } @@ -399,27 +406,33 @@ where q.lu().solve(&p).unwrap() } -fn one_norm(m: &MatrixN) -> N +fn one_norm(m: &MatrixN) -> N::RealField where - N: RealField, + N: ComplexField, D: Dim, DefaultAllocator: Allocator, { - let mut max = N::zero(); + let mut max = ::RealField::zero(); for i in 0..m.ncols() { let col = m.column(i); - max = max.max(col.iter().fold(N::zero(), |a, b| a + b.abs())); + max = max.max( + col.iter() + .fold(::RealField::zero(), |a, b| a + b.abs()), + ); } max } -impl MatrixN +impl MatrixN where D: DimMin, - DefaultAllocator: - Allocator + Allocator<(usize, usize), DimMinimum> + Allocator, + DefaultAllocator: Allocator + + Allocator<(usize, usize), DimMinimum> + + Allocator + + Allocator + + Allocator, { /// Computes exponential of this matrix pub fn exp(&self) -> Self { @@ -430,19 +443,19 @@ where let mut h = ExpmPadeHelper::new(self.clone(), true); - let eta_1 = N::max(h.d4_loose(), h.d6_loose()); + let eta_1 = N::RealField::max(h.d4_loose(), h.d6_loose()); if eta_1 < convert(1.495585217958292e-002) && ell(&h.a, 3) == 0 { let (u, v) = h.pade3(); return solve_p_q(u, v); } - let eta_2 = N::max(h.d4_tight(), h.d6_loose()); + let eta_2 = N::RealField::max(h.d4_tight(), h.d6_loose()); if eta_2 < convert(2.539398330063230e-001) && ell(&h.a, 5) == 0 { let (u, v) = h.pade5(); return solve_p_q(u, v); } - let eta_3 = N::max(h.d6_tight(), h.d8_loose()); + let eta_3 = N::RealField::max(h.d6_tight(), h.d8_loose()); if eta_3 < convert(9.504178996162932e-001) && ell(&h.a, 7) == 0 { let (u, v) = h.pade7(); return solve_p_q(u, v); @@ -452,11 +465,11 @@ where return solve_p_q(u, v); } - let eta_4 = N::max(h.d8_loose(), h.d10_loose()); - let eta_5 = N::min(eta_3, eta_4); + let eta_4 = N::RealField::max(h.d8_loose(), h.d10_loose()); + let eta_5 = N::RealField::min(eta_3, eta_4); let theta_13 = convert(4.25); - let mut s = if eta_5 == N::zero() { + let mut s = if eta_5 == N::RealField::zero() { 0 } else { let l2 = try_convert((eta_5 / theta_13).log2().ceil()).unwrap(); diff --git a/src/linalg/full_piv_lu.rs b/src/linalg/full_piv_lu.rs index bfa6854f..146b2cb7 100644 --- a/src/linalg/full_piv_lu.rs +++ b/src/linalg/full_piv_lu.rs @@ -60,11 +60,7 @@ where let mut q = PermutationSequence::identity_generic(min_nrows_ncols); if min_nrows_ncols.value() == 0 { - return Self { - lu: matrix, - p: p, - q: q, - }; + return Self { lu: matrix, p, q }; } for i in 0..min_nrows_ncols.value() { @@ -90,11 +86,7 @@ where } } - Self { - lu: matrix, - p: p, - q: q, - } + Self { lu: matrix, p, q } } #[doc(hidden)] diff --git a/src/linalg/lu.rs b/src/linalg/lu.rs index b3acaaa0..6aa84d62 100644 --- a/src/linalg/lu.rs +++ b/src/linalg/lu.rs @@ -96,7 +96,7 @@ where let mut p = PermutationSequence::identity_generic(min_nrows_ncols); if min_nrows_ncols.value() == 0 { - return LU { lu: matrix, p: p }; + return LU { lu: matrix, p }; } for i in 0..min_nrows_ncols.value() { @@ -117,7 +117,7 @@ where } } - LU { lu: matrix, p: p } + LU { lu: matrix, p } } #[doc(hidden)] diff --git a/src/linalg/mod.rs b/src/linalg/mod.rs index f96cef0c..4a1edbcf 100644 --- a/src/linalg/mod.rs +++ b/src/linalg/mod.rs @@ -5,6 +5,10 @@ mod bidiagonal; mod cholesky; mod convolution; mod determinant; +// FIXME: this should not be needed. However, the exp uses +// explicit float operations on `f32` and `f64`. We need to +// get rid of these to allow exp to be used on a no-std context. +#[cfg(feature = "std")] mod exp; mod full_piv_lu; pub mod givens; @@ -27,6 +31,7 @@ mod symmetric_tridiagonal; pub use self::bidiagonal::*; pub use self::cholesky::*; pub use self::convolution::*; +#[cfg(feature = "std")] pub use self::exp::*; pub use self::full_piv_lu::*; pub use self::hessenberg::*; diff --git a/src/linalg/qr.rs b/src/linalg/qr.rs index f218d47b..48d7a896 100644 --- a/src/linalg/qr.rs +++ b/src/linalg/qr.rs @@ -57,20 +57,14 @@ where let mut diag = unsafe { MatrixMN::new_uninitialized_generic(min_nrows_ncols, U1) }; if min_nrows_ncols.value() == 0 { - return QR { - qr: matrix, - diag: diag, - }; + return QR { qr: matrix, diag }; } for ite in 0..min_nrows_ncols.value() { householder::clear_column_unchecked(&mut matrix, &mut diag[ite], ite, 0, None); } - QR { - qr: matrix, - diag: diag, - } + QR { qr: matrix, diag } } /// Retrieves the upper trapezoidal submatrix `R` of this decomposition. diff --git a/src/linalg/schur.rs b/src/linalg/schur.rs index 99e73ea1..8d8be72d 100644 --- a/src/linalg/schur.rs +++ b/src/linalg/schur.rs @@ -73,10 +73,8 @@ where pub fn try_new(m: MatrixN, eps: N::RealField, max_niter: usize) -> Option { let mut work = unsafe { VectorN::new_uninitialized_generic(m.data.shape().0, U1) }; - Self::do_decompose(m, &mut work, eps, max_niter, true).map(|(q, t)| Schur { - q: q.unwrap(), - t: t, - }) + Self::do_decompose(m, &mut work, eps, max_niter, true) + .map(|(q, t)| Schur { q: q.unwrap(), t }) } fn do_decompose( @@ -138,9 +136,9 @@ where let m = end - 1; let n = end; - let h11 = t[(start + 0, start + 0)]; - let h12 = t[(start + 0, start + 1)]; - let h21 = t[(start + 1, start + 0)]; + let h11 = t[(start, start)]; + let h12 = t[(start, start + 1)]; + let h21 = t[(start + 1, start)]; let h22 = t[(start + 1, start + 1)]; let h32 = t[(start + 2, start + 1)]; @@ -163,7 +161,7 @@ where if not_zero { if k > start { - t[(k + 0, k - 1)] = norm; + t[(k, k - 1)] = norm; t[(k + 1, k - 1)] = N::zero(); t[(k + 2, k - 1)] = N::zero(); } diff --git a/src/linalg/svd.rs b/src/linalg/svd.rs index 4be0ecdb..94ec4d96 100644 --- a/src/linalg/svd.rs +++ b/src/linalg/svd.rs @@ -218,9 +218,9 @@ where } } - diagonal[k + 0] = subm[(0, 0)]; + diagonal[k] = subm[(0, 0)]; diagonal[k + 1] = subm[(1, 1)]; - off_diagonal[k + 0] = subm[(0, 1)]; + off_diagonal[k] = subm[(0, 1)]; if k != n - 1 { off_diagonal[k + 1] = subm[(1, 2)]; @@ -244,7 +244,7 @@ where let u2 = u2.map(|u2| GivensRotation::new_unchecked(u2.c(), N::from_real(u2.s()))); let v2 = v2.map(|v2| GivensRotation::new_unchecked(v2.c(), N::from_real(v2.s()))); - diagonal[start + 0] = s[0]; + diagonal[start] = s[0]; diagonal[start + 1] = s[1]; off_diagonal[start] = N::RealField::zero(); diff --git a/src/linalg/symmetric_eigen.rs b/src/linalg/symmetric_eigen.rs index 70f57cdf..2d62816e 100644 --- a/src/linalg/symmetric_eigen.rs +++ b/src/linalg/symmetric_eigen.rs @@ -192,7 +192,7 @@ where let eigvals = m.eigenvalues().unwrap(); let basis = Vector2::new(eigvals.x - diag[start + 1], off_diag[start]); - diag[start + 0] = eigvals[0]; + diag[start] = eigvals[0]; diag[start + 1] = eigvals[1]; if let Some(ref mut q) = q { diff --git a/tests/core/cg.rs b/tests/core/cg.rs new file mode 100644 index 00000000..b061efbb --- /dev/null +++ b/tests/core/cg.rs @@ -0,0 +1,59 @@ +use na::{Matrix3, Matrix4, Point2, Point3, Vector2, Vector3}; + +/// See Example 3.4 of "Graphics and Visualization: Principles & Algorithms" +/// by Theoharis, Papaioannou, Platis, Patrikalakis. +#[test] +fn test_scaling_wrt_point_1() { + let a = Point2::new(0.0, 0.0); + let b = Point2::new(1.0, 1.0); + let c = Point2::new(5.0, 2.0); + + let scaling = Vector2::new(2.0, 2.0); + let scale_about = Matrix3::new_nonuniform_scaling_wrt_point(&scaling, &c); + + let expected_a = Point2::new(-5.0, -2.0); + let expected_b = Point2::new(-3.0, 0.0); + let result_a = scale_about.transform_point(&a); + let result_b = scale_about.transform_point(&b); + let result_c = scale_about.transform_point(&c); + + assert!(expected_a == result_a); + assert!(expected_b == result_b); + assert!(c == result_c); +} + +/// Based on the same example as the test above. +#[test] +fn test_scaling_wrt_point_2() { + let a = Point3::new(0.0, 0.0, 1.0); + let b = Point3::new(1.0, 1.0, 1.0); + let c = Point3::new(5.0, 2.0, 1.0); + + let scaling = Vector3::new(2.0, 2.0, 1.0); + let scale_about = Matrix4::new_nonuniform_scaling_wrt_point(&scaling, &c); + + let expected_a = Point3::new(-5.0, -2.0, 1.0); + let expected_b = Point3::new(-3.0, 0.0, 1.0); + + let result_a = scale_about.transform_point(&a); + let result_b = scale_about.transform_point(&b); + let result_c = scale_about.transform_point(&c); + + assert!(expected_a == result_a); + assert!(expected_b == result_b); + assert!(c == result_c); +} + +/// Based on https://github.com/emlowry/AiE/blob/50bae4068edb686cf8ffacdf6fab8e7cb22e7eb1/Year%201%20Classwork/MathTest/Matrix4x4TestGroup.cpp#L145 +#[test] +fn test_scaling_wrt_point_3() { + let about = Point3::new(2.0, 1.0, -2.0); + let scale = Vector3::new(2.0, 0.5, -1.0); + let pt = Point3::new(1.0, 2.0, 3.0); + let scale_about = Matrix4::new_nonuniform_scaling_wrt_point(&scale, &about); + + let expected = Point3::new(0.0, 1.5, -7.0); + let result = scale_about.transform_point(&pt); + + assert!(result == expected); +} diff --git a/tests/core/matrixcompare.rs b/tests/core/matrixcompare.rs new file mode 100644 index 00000000..df112173 --- /dev/null +++ b/tests/core/matrixcompare.rs @@ -0,0 +1,85 @@ +//! Tests for `matrixcompare` integration. +//! +//! The `matrixcompare` crate itself is responsible for testing the actual comparison. +//! The tests here only check that the necessary trait implementations are correctly implemented, +//! in addition to some sanity checks with example input. + +use nalgebra::{DMatrix, MatrixMN, U4, U5}; + +use matrixcompare::{assert_matrix_eq, DenseAccess}; + +#[cfg(feature = "arbitrary")] +quickcheck! { + fn fetch_single_is_equivalent_to_index_f64(matrix: DMatrix) -> bool { + for i in 0 .. matrix.nrows() { + for j in 0 .. matrix.ncols() { + if matrix.fetch_single(i, j) != *matrix.index((i, j)) { + return false; + } + } + } + + true + } + + fn matrixcompare_shape_agrees_with_matrix(matrix: DMatrix) -> bool { + matrix.nrows() == as matrixcompare::Matrix>::rows(&matrix) + && + matrix.ncols() == as matrixcompare::Matrix>::cols(&matrix) + } +} + +#[test] +fn assert_matrix_eq_dense_positive_comparison() { + #[rustfmt::skip] + let a = MatrixMN::<_, U4, U5>::from_row_slice(&[ + 1210, 1320, 1430, 1540, 1650, + 2310, 2420, 2530, 2640, 2750, + 3410, 3520, 3630, 3740, 3850, + 4510, 4620, 4730, 4840, 4950, + ]); + + #[rustfmt::skip] + let b = MatrixMN::<_, U4, U5>::from_row_slice(&[ + 1210, 1320, 1430, 1540, 1650, + 2310, 2420, 2530, 2640, 2750, + 3410, 3520, 3630, 3740, 3850, + 4510, 4620, 4730, 4840, 4950, + ]); + + // Test matrices of static size + assert_matrix_eq!(a, b); + assert_matrix_eq!(&a, b); + assert_matrix_eq!(a, &b); + assert_matrix_eq!(&a, &b); + + // Test matrices of dynamic size + let a_dyn = a.index((0..4, 0..5)); + let b_dyn = b.index((0..4, 0..5)); + assert_matrix_eq!(a_dyn, b_dyn); + assert_matrix_eq!(a_dyn, &b_dyn); + assert_matrix_eq!(&a_dyn, b_dyn); + assert_matrix_eq!(&a_dyn, &b_dyn); +} + +#[test] +#[should_panic] +fn assert_matrix_eq_dense_negative_comparison() { + #[rustfmt::skip] + let a = MatrixMN::<_, U4, U5>::from_row_slice(&[ + 1210, 1320, 1430, 1540, 1650, + 2310, 2420, 2530, 2640, 2750, + 3410, 3520, 3630, 3740, 3850, + 4510, 4620, -4730, 4840, 4950, + ]); + + #[rustfmt::skip] + let b = MatrixMN::<_, U4, U5>::from_row_slice(&[ + 1210, 1320, 1430, 1540, 1650, + 2310, 2420, 2530, 2640, 2750, + 3410, 3520, 3630, 3740, 3850, + 4510, 4620, 4730, 4840, 4950, + ]); + + assert_matrix_eq!(a, b); +} diff --git a/tests/core/mod.rs b/tests/core/mod.rs index ec1c4e3e..4f9bc745 100644 --- a/tests/core/mod.rs +++ b/tests/core/mod.rs @@ -1,6 +1,7 @@ #[cfg(feature = "abomonation-serialize")] mod abomonation; mod blas; +mod cg; mod conversion; mod edition; mod empty; @@ -10,5 +11,8 @@ mod matrix_slice; mod mint; mod serde; +#[cfg(feature = "compare")] +mod matrixcompare; + #[cfg(feature = "arbitrary")] pub mod helper; diff --git a/tests/linalg/exp.rs b/tests/linalg/exp.rs index 75122107..f5b5243a 100644 --- a/tests/linalg/exp.rs +++ b/tests/linalg/exp.rs @@ -126,4 +126,51 @@ mod tests { assert!(relative_eq!(f, m.exp(), epsilon = 1.0e-7)); } + + #[test] + fn exp_complex() { + use nalgebra::{Complex, DMatrix, DVector, Matrix2, RealField}; + + { + let z = Matrix2::>::zeros(); + + let identity = Matrix2::>::identity(); + + assert!((z.exp() - identity).norm() < 1e-7); + } + + { + let a = Matrix2::>::new( + Complex::::new(0.0, 1.0), + Complex::::new(0.0, 2.0), + Complex::::new(0.0, -1.0), + Complex::::new(0.0, 3.0), + ); + + let b = Matrix2::>::new( + Complex::::new(0.42645929666726, 1.89217550966333), + Complex::::new(-2.13721484276556, -0.97811251808259), + Complex::::new(1.06860742138278, 0.48905625904129), + Complex::::new(-1.7107555460983, 0.91406299158075), + ); + + assert!((a.exp() - b).norm() < 1.0e-07); + } + + { + let d1 = Complex::::new(0.0, ::pi()); + let d2 = Complex::::new(0.0, ::frac_pi_2()); + let d3 = Complex::::new(0.0, ::frac_pi_4()); + + let m = DMatrix::>::from_diagonal(&DVector::from_row_slice(&[d1, d2, d3])); + + let res = DMatrix::>::from_diagonal(&DVector::from_row_slice(&[ + d1.exp(), + d2.exp(), + d3.exp(), + ])); + + assert!((m.exp() - res).norm() < 1e-07); + } + } }