diff --git a/include/xo/unit/Quantity.hpp b/include/xo/unit/Quantity.hpp index 15265046..98d086d8 100644 --- a/include/xo/unit/Quantity.hpp +++ b/include/xo/unit/Quantity.hpp @@ -1,4 +1,4 @@ -/** @file Quantity2.hpp +/** @file Quantity.hpp * * Author: Roland Conybeare **/ @@ -34,7 +34,7 @@ namespace xo { public: constexpr Quantity(Repr scale, - const natural_unit & unit) + const natural_unit & unit) : scale_{scale}, unit_{unit} {} constexpr const repr_type & scale() const { return scale_; } @@ -63,6 +63,29 @@ namespace xo { rr.natural_unit_); } + template + static constexpr + auto divide(const Quantity & x, const Quantity2 & y) { + using r_repr_type = std::common_type_t; + using r_int_type = std::common_type_t; + + auto rr = detail::nu_ratio(x.unit(), y.unit()); + + /* note: nu_ratio() reports multiplicative outer scaling factors, + * so multiply is correct here + */ + r_repr_type r_scale = (::sqrt(rr.outer_scale_sq_) + * rr.outer_scale_exact_.template to() + * static_cast(x.scale()) + / static_cast(y.scale())); + + return Quantity(r_scale, + rr.natural_unit_); + } + + private: /** @brief quantity represents this multiple of a unit amount **/ Repr scale_ = Repr{}; @@ -89,8 +112,17 @@ namespace xo { { return Quantity::multiply(x, y); } + + /** note: won't have constexpr result until c++26 (when ::sqrt(), ::pow() are constexpr) + **/ + template + constexpr auto + operator/ (const Quantity & x, const OtherQuantity & y) + { + return Quantity::divide(x, y); + } } /*namespace qty*/ } /*namespace xo*/ -/** end Quantity2.hpp **/ +/** end Quantity.hpp **/ diff --git a/include/xo/unit/natural_unit.hpp b/include/xo/unit/natural_unit.hpp index 414266ac..15a7f65c 100644 --- a/include/xo/unit/natural_unit.hpp +++ b/include/xo/unit/natural_unit.hpp @@ -16,10 +16,11 @@ namespace xo { /** @class natural_unit * @brief an array representing the cartesian product of distinct basis-power-units * - * 1. Each bpu in the array represents a power of a basis dimension, e.g. "meter" or "second^2". - * 2. Each bpu in an array has a different dimension id. + * 1. Quantities are represented as a multiple of a natural unit + * 2. Each bpu in the array represents a power of a basis dimension, e.g. "meter" or "second^2". + * 3. Each bpu in an array has a different dimension id. * For example dim::time, if present, appears once. - * 3. Basis dimensions can appear in any order. + * 4. Basis dimensions can appear in any order. * Order used for constructing abbreviations: will get @c "kg.m" or @c "m.kg" * depending on the orderin of @c dim::distance and @c dim::mass in @c bpu_v_ **/ @@ -46,6 +47,14 @@ namespace xo { return retval; } + /** @brief remove bpu at position @p p **/ + constexpr void remove_bpu(size_t p) { + for (std::size_t i = p; i+1 < n_bpu_; ++i) + bpu_v_[i] = bpu_v_[i+1]; + + --n_bpu_; + } + constexpr void push_back(const bpu & bpu) { if (n_bpu_ < n_dim) bpu_v_[n_bpu_++] = bpu; @@ -187,8 +196,8 @@ namespace xo { p_target_bpu->scalefactor()); *p_target_bpu = bpu(p_target_bpu->native_dim(), - p_target_bpu->scalefactor(), - p_target_bpu->power() + rhs_bpu_orig.power()); + p_target_bpu->scalefactor(), + p_target_bpu->power() + rhs_bpu_orig.power()); return outer_scalefactor_result(rhs_bpu_rr.outer_scale_exact_, rhs_bpu_rr.outer_scale_sq_); @@ -197,21 +206,49 @@ namespace xo { template constexpr outer_scalefactor_result - bpu_array_product_inplace(natural_unit * p_target, - const bpu & bpu) + bpu_ratio_inplace(bpu * p_target_bpu, + const bpu & rhs_bpu_orig) + { + assert(rhs_bpu_orig.native_dim() == p_target_bpu->native_dim()); + + bpu2_rescale_result rhs_bpu_rr = bpu2_rescale(rhs_bpu_orig, + p_target_bpu->scalefactor()); + + *p_target_bpu = bpu(p_target_bpu->native_dim(), + p_target_bpu->scalefactor(), + p_target_bpu->power() - rhs_bpu_orig.power()); + + return outer_scalefactor_result(power_ratio_type(1,1) / rhs_bpu_rr.outer_scale_exact_, + 1.0 / rhs_bpu_rr.outer_scale_sq_); + } + + template + constexpr + outer_scalefactor_result + nu_product_inplace(natural_unit * p_target, + const bpu & bpu) { std::size_t i = 0; for (; i < p_target->n_bpu(); ++i) { - if ((*p_target)[i].native_dim() == bpu.native_dim()) { - outer_scalefactor_result retval = bpu_product_inplace(&((*p_target)[i]), bpu); + auto * p_target_bpu = &((*p_target)[i]); - /* TODO: strip 0 power */ + if (p_target_bpu->native_dim() == bpu.native_dim()) { + outer_scalefactor_result retval = bpu_product_inplace(p_target_bpu, bpu); + + if (p_target_bpu->power().is_zero()) { + /* dimension assoc'd with *p_target_bpu has been cancelled */ + p_target->remove_bpu(i); + } return retval; } } - /* control here: i=p_target->n_bpu() */ + /* control here: i=p_target->n_bpu() + * Dimension represented by bpu does not already appear in *p_target. + * Adopt bpu's scalefactor + */ + p_target->push_back(bpu); return outer_scalefactor_result @@ -219,6 +256,40 @@ namespace xo { 1.0 /*outer_scale_sq*/); } + template + constexpr + outer_scalefactor_result + nu_ratio_inplace(natural_unit * p_target, + const bpu & bpu) + { + std::size_t i = 0; + for (; i < p_target->n_bpu(); ++i) { + auto * p_target_bpu = &((*p_target)[i]); + + if (p_target_bpu->native_dim() == bpu.native_dim()) { + outer_scalefactor_result retval = bpu_ratio_inplace(p_target_bpu, bpu); + + if (p_target_bpu->power().is_zero()) { + /* dimension assoc'd with *p_target_bpu has been cancelled */ + p_target->remove_bpu(i); + } + + return retval; + } + } + + /* here: i=p_target->n_bpu() + * Dimension represented by bpu does not already appear in *p_target. + * Adopt bpu's scalefactor + */ + + p_target->push_back(bpu.reciprocal()); + + return outer_scalefactor_result + (ratio::ratio(1, 1) /*outer_scale_exact*/, + 1.0 /*outer_scale_sq*/); + } + template constexpr natural_unit nu_reciprocal(const natural_unit & nu) diff --git a/include/xo/unit/scaled_unit.hpp b/include/xo/unit/scaled_unit.hpp index 16766130..6f0d9268 100644 --- a/include/xo/unit/scaled_unit.hpp +++ b/include/xo/unit/scaled_unit.hpp @@ -67,15 +67,15 @@ namespace xo { template constexpr scaled_unit - nu_product(const natural_unit & lhs_bpu_array, - const bpu & rhs_bpu) + nu_bpu_product(const natural_unit & lhs_bpu_array, + const bpu & rhs_bpu) { natural_unit prod = lhs_bpu_array; - auto rr = bpu_array_product_inplace(&prod, rhs_bpu); + auto rr = nu_product_inplace(&prod, rhs_bpu); return scaled_unit(prod, - rr.outer_scale_exact_, - rr.outer_scale_sq_); + rr.outer_scale_exact_, + rr.outer_scale_sq_); }; template @@ -95,7 +95,7 @@ namespace xo { 1.0 /*outer_scale_sq*/)); for (std::size_t i = 0; i < rhs_bpu_array.n_bpu(); ++i) { - auto sfr2 = bpu_array_product_inplace(&prod, rhs_bpu_array[i]); + auto sfr2 = nu_product_inplace(&prod, rhs_bpu_array[i]); sfr.outer_scale_exact_ = sfr.outer_scale_exact_ * sfr2.outer_scale_exact_; sfr.outer_scale_sq_ *= sfr2.outer_scale_sq_; @@ -106,6 +106,36 @@ namespace xo { sfr.outer_scale_sq_); } + template + constexpr + scaled_unit + nu_ratio(const natural_unit & nu_lhs, + const natural_unit & nu_rhs) + { + natural_unit ratio = nu_lhs; + + /* accumulate product of scalefactors spun off by rescaling + * any basis-units in rhs_bpu_array that conflict with the same dimension + * in lh_bpu_array + */ + auto sfr = (detail::outer_scalefactor_result + (scalefactor_ratio_type(1, 1) /*outer_scale_exact*/, + 1.0 /*outer_scale_sq*/)); + + for (std::size_t i = 0; i < nu_rhs.n_bpu(); ++i) { + auto sfr2 = nu_ratio_inplace(&ratio, nu_rhs[i]); + + /* note: nu_ratio_inplace() reports multiplicative outer scaling factors, + * so multiply is correct here + */ + sfr.outer_scale_exact_ = sfr.outer_scale_exact_ * sfr2.outer_scale_exact_; + sfr.outer_scale_sq_ *= sfr2.outer_scale_sq_; + } + + return scaled_unit(ratio, + sfr.outer_scale_exact_, + sfr.outer_scale_sq_); + } } template @@ -114,7 +144,21 @@ namespace xo { const scaled_unit & y_unit) { auto rr = detail::nu_product(x_unit.natural_unit_, - y_unit.natural_unit_); + y_unit.natural_unit_); + + return (scaled_unit + (rr.natural_unit_, + rr.outer_scale_exact_ * x_unit.outer_scale_exact_ * y_unit.outer_scale_exact_, + rr.outer_scale_sq_ * x_unit.outer_scale_sq_ * y_unit.outer_scale_sq_)); + } + + template + inline constexpr scaled_unit + operator/ (const scaled_unit & x_unit, + const scaled_unit & y_unit) + { + auto rr = detail::nu_ratio(x_unit.natural_unit_, + y_unit.natural_unit_); return (scaled_unit (rr.natural_unit_, diff --git a/utest/unit.test.cpp b/utest/unit.test.cpp index 7d0d6856..a0404279 100644 --- a/utest/unit.test.cpp +++ b/utest/unit.test.cpp @@ -46,6 +46,9 @@ namespace xo { using xo::qty::natural_unit; using xo::qty::bpu_array_maker; using xo::qty::detail::nu_product; + using xo::qty::detail::nu_bpu_product; + using xo::qty::detail::nu_ratio_inplace; + using xo::qty::detail::nu_ratio; using xo::qty::unit_qty; using xo::print::ccs; @@ -595,7 +598,7 @@ namespace xo { static_assert(bpu.power() == power_ratio_type(-1, 2)); - constexpr auto prod_rr = nu_product(v, bpu); + constexpr auto prod_rr = nu_bpu_product(v, bpu); log && log(xtag("prod_rr.bpu_array", prod_rr.natural_unit_)); log && log(xtag("prod_rr.outer_scale_exact", prod_rr.outer_scale_exact_)); @@ -642,7 +645,7 @@ namespace xo { static_assert(bpu.power() == power_ratio_type(-1, 2)); - constexpr auto prod_rr = nu_product(v, bpu); + constexpr auto prod_rr = nu_bpu_product(v, bpu); log && log(xtag("prod_rr.bpu_array", prod_rr.natural_unit_)); log && log(xtag("prod_rr.outer_scale_exact", prod_rr.outer_scale_exact_)); @@ -686,7 +689,7 @@ namespace xo { static_assert(bpu.power() == power_ratio_type(-1, 1)); - constexpr auto prod_rr = nu_product(v, bpu); + constexpr auto prod_rr = nu_bpu_product(v, bpu); log && log(xtag("prod_rr.bpu_array", prod_rr.natural_unit_)); log && log(xtag("prod_rr.outer_scale_exact", prod_rr.outer_scale_exact_)); @@ -832,6 +835,104 @@ namespace xo { } } /*TEST_CASE(natural_unit2)*/ + TEST_CASE("natural_unit3", "[natural_unit]") { + constexpr bool c_debug_flag = true; + + // can get bits from /dev/random by uncommenting the 2nd line below + //uint64_t seed = xxx; + //rng::Seed seed; + + //auto rng = xo::rng::xoshiro256ss(seed); + + scope log(XO_DEBUG2(c_debug_flag, "TEST_CASE.natural_unit3")); + //log && log("(A)", xtag("foo", foo)); + + { + constexpr natural_unit v + = (bpu_array_maker::make_bpu_array + (bpu(dim::mass, scalefactor_ratio_type(1000, 1), power_ratio_type(1, 1)), + bpu(dim::distance, scalefactor_ratio_type(1, 1), power_ratio_type(1, 1)))); + + static_assert(v.n_bpu() == 2); + + log && log(xtag("v.abbrev", v.abbrev())); + + static_assert(v.abbrev().size() > 0); + static_assert(v.abbrev() == flatstring("kg.m")); + + { + natural_unit w = v; + + nu_ratio_inplace(&w, + bpu(dim::mass, scalefactor_ratio_type(1000, 1), power_ratio_type(1, 1))); + + log && log(xtag("w.abbrev", w.abbrev())); + + REQUIRE(w.n_bpu() == 1); + REQUIRE(w[0].native_dim() == dim::distance); + REQUIRE(w.abbrev() == flatstring("m")); + } + + { + constexpr natural_unit w + = (bpu_array_maker::make_bpu_array + (bpu(dim::mass, scalefactor_ratio_type(1000, 1), power_ratio_type(1, 1)))); + + static_assert(w.n_bpu() == 1); + + log && log(xtag("w.abbrev", w.abbrev())); + + constexpr auto rr = nu_ratio(v, w); + + log && log(xtag("rr", rr)); + + REQUIRE(rr.natural_unit_.n_bpu() == 1); + REQUIRE(rr.natural_unit_[0].native_dim() == dim::distance); + REQUIRE(rr.natural_unit_.abbrev() == flatstring("m")); + } + + { + constexpr natural_unit w + = (bpu_array_maker::make_bpu_array + (bpu(dim::time, scalefactor_ratio_type(1, 1), power_ratio_type(1, 1)))); + + static_assert(w.n_bpu() == 1); + + log && log(xtag("w.abbrev", w.abbrev())); + + constexpr auto rr = nu_ratio(v, w); + + log && log(xtag("rr", rr)); + + REQUIRE(rr.natural_unit_.n_bpu() == 3); + REQUIRE(rr.natural_unit_[0].native_dim() == dim::mass); + REQUIRE(rr.natural_unit_[1].native_dim() == dim::distance); + REQUIRE(rr.natural_unit_[2].native_dim() == dim::time); + REQUIRE(rr.natural_unit_.abbrev() == flatstring("kg.m.s^-1")); + } + + { + natural_unit w = v; + + REQUIRE(w.n_bpu() == 2); + REQUIRE(w[0].native_dim() == dim::mass); + + nu_ratio_inplace(&w, + bpu(dim::time, scalefactor_ratio_type(1, 1), power_ratio_type(2, 1))); + + REQUIRE(w.n_bpu() == 3); + REQUIRE(w[0].native_dim() == dim::mass); + REQUIRE(w[1].native_dim() == dim::distance); + REQUIRE(w[2].native_dim() == dim::time); + + log && log(xtag("w.abbrev", w.abbrev())); + + REQUIRE(w.n_bpu() == 3); + REQUIRE(w.abbrev() == flatstring("kg.m.s^-2")); + } + } + } /*TEST_CASE(natural_unit3)*/ + TEST_CASE("scaled_unit0", "[scaled_unit0]") { constexpr bool c_debug_flag = true; @@ -855,6 +956,29 @@ namespace xo { static_assert(ng2.natural_unit_.n_bpu() == 1); } /*TEST_CASE(scaled_unit0)*/ + TEST_CASE("scaled_unit1", "[scaled_unit1]") { + constexpr bool c_debug_flag = true; + + // can get bits from /dev/random by uncommenting the 2nd line below + //uint64_t seed = xxx; + //rng::Seed seed; + + //auto rng = xo::rng::xoshiro256ss(seed); + + scope log(XO_DEBUG2(c_debug_flag, "TEST_CASE.scaled_unit1")); + //log && log("(A)", xtag("foo", foo)); + + constexpr auto ng = su2::nanogram; + constexpr auto ng2 = ng / ng; + + log && log(xtag("ng", ng)); + log && log(xtag("ng/ng", ng2)); + //log && log(xtag("ng/ng", + + static_assert(ng.natural_unit_.n_bpu() == 1); + static_assert(ng2.natural_unit_.n_bpu() == 0); + } /*TEST_CASE(scaled_unit1)*/ + TEST_CASE("Quantity", "[Quantity]") { constexpr bool c_debug_flag = true; @@ -908,6 +1032,27 @@ namespace xo { scope log(XO_DEBUG2(c_debug_flag, "TEST_CASE.Quantity3")); //log && log("(A)", xtag("foo", foo)); + /* not constexpr until c++26 */ + Quantity ng = unit_qty(su2::nanogram); + auto ng0 = ng / ng; + + log && log(xtag("ng/ng", ng0)); + + REQUIRE(ng0.scale() == 1); + } /*TEST_CASE(Quantity3)*/ + + TEST_CASE("Quantity4", "[Quantity]") { + constexpr bool c_debug_flag = true; + + // can get bits from /dev/random by uncommenting the 2nd line below + //uint64_t seed = xxx; + //rng::Seed seed; + + //auto rng = xo::rng::xoshiro256ss(seed); + + scope log(XO_DEBUG2(c_debug_flag, "TEST_CASE.Quantity4")); + //log && log("(A)", xtag("foo", foo)); + /* not constexpr until c++26 */ Quantity ng = unit_qty(su2::nanogram); Quantity ug = unit_qty(su2::microgram); @@ -937,7 +1082,43 @@ namespace xo { } //REQUIRE(ng2.scale() == 1); - } /*TEST_CASE(Quantity3)*/ + } /*TEST_CASE(Quantity4)*/ + + TEST_CASE("Quantity5", "[Quantity]") { + constexpr bool c_debug_flag = true; + + // can get bits from /dev/random by uncommenting the 2nd line below + //uint64_t seed = xxx; + //rng::Seed seed; + + //auto rng = xo::rng::xoshiro256ss(seed); + + scope log(XO_DEBUG2(c_debug_flag, "TEST_CASE.Quantity5")); + //log && log("(A)", xtag("foo", foo)); + + /* not constexpr until c++26 */ + Quantity ng = unit_qty(su2::nanogram); + Quantity ug = unit_qty(su2::microgram); + + { + auto ratio1 = ng / ug; + log && log(xtag("ng/ug", ratio1)); + + /* units will be nanograms, since that's on lhs */ + REQUIRE(ratio1.unit().n_bpu() == 0); + REQUIRE(ratio1.scale() == 0.001); + } + + { + auto ratio2 = ug / ng; + log && log(xtag("ug/ng", ratio2)); + + REQUIRE(ratio2.unit().n_bpu() == 0); + REQUIRE(ratio2.scale() == 1000.0); + } + + //REQUIRE(ng2.scale() == 1); + } /*TEST_CASE(Quantity5)*/ } /*namespace ut*/ } /*namespace xo*/