File indexing completed on 2026-05-05 08:35:21
0001
0002
0003
0004
0005
0006
0007 #pragma once
0008
0009 #include <cmath>
0010
0011 #include "corecel/Assert.hh"
0012 #include "corecel/Macros.hh"
0013 #include "corecel/Types.hh"
0014
0015 #include "detail/SoftEqualTraits.hh"
0016
0017 namespace celeritas
0018 {
0019
0020
0021
0022
0023
0024
0025
0026
0027
0028
0029 CELER_CONSTEXPR_FUNCTION real_type sqrt_tol()
0030 {
0031 return detail::SoftEqualTraits<real_type>::sqrt_prec();
0032 }
0033
0034
0035
0036
0037
0038
0039
0040
0041
0042
0043
0044
0045
0046
0047
0048
0049
0050
0051
0052
0053 template<class RealType = ::celeritas::real_type>
0054 class SoftEqual
0055 {
0056 public:
0057
0058
0059 using value_type = RealType;
0060
0061
0062 public:
0063
0064 CELER_CONSTEXPR_FUNCTION SoftEqual();
0065
0066
0067 explicit CELER_FUNCTION SoftEqual(value_type rel);
0068
0069
0070 CELER_FUNCTION SoftEqual(value_type rel, value_type abs);
0071
0072
0073
0074
0075 bool CELER_FUNCTION operator()(value_type a, value_type b) const;
0076
0077
0078
0079
0080 CELER_CONSTEXPR_FUNCTION value_type rel() const { return rel_; }
0081
0082
0083 CELER_CONSTEXPR_FUNCTION value_type abs() const { return abs_; }
0084
0085 private:
0086 value_type rel_;
0087 value_type abs_;
0088
0089 using SETraits = detail::SoftEqualTraits<value_type>;
0090 };
0091
0092
0093
0094
0095
0096
0097
0098 template<class F>
0099 class EqualOr : public F
0100 {
0101 public:
0102
0103 template<class... C>
0104 CELER_FUNCTION EqualOr(C&&... args) : F{std::forward<C>(args)...}
0105 {
0106 }
0107
0108
0109 template<class T, class U>
0110 bool CELER_FUNCTION operator()(T a, U b) const
0111 {
0112 return a == b || static_cast<F const&>(*this)(a, b);
0113 }
0114 };
0115
0116
0117
0118
0119
0120 template<class RealType = ::celeritas::real_type>
0121 class SoftClose
0122 {
0123 public:
0124
0125
0126 using value_type = RealType;
0127
0128
0129 public:
0130
0131 CELER_CONSTEXPR_FUNCTION SoftClose();
0132
0133
0134 CELER_FUNCTION SoftClose(value_type abs);
0135
0136
0137
0138
0139 bool CELER_FUNCTION operator()(value_type a, value_type b) const;
0140
0141
0142
0143
0144 CELER_CONSTEXPR_FUNCTION value_type abs() const { return abs_; }
0145
0146 private:
0147 value_type abs_;
0148
0149 using SETraits = detail::SoftEqualTraits<value_type>;
0150 };
0151
0152
0153
0154
0155
0156 template<class RealType = ::celeritas::real_type>
0157 class SoftZero
0158 {
0159 public:
0160
0161
0162 using argument_type = RealType;
0163 using value_type = RealType;
0164
0165
0166 public:
0167
0168 CELER_CONSTEXPR_FUNCTION SoftZero();
0169
0170
0171 explicit CELER_FUNCTION SoftZero(value_type abs);
0172
0173
0174
0175
0176 inline CELER_FUNCTION bool operator()(value_type actual) const;
0177
0178
0179
0180
0181 CELER_CONSTEXPR_FUNCTION value_type abs() const { return abs_; }
0182
0183 private:
0184 value_type abs_;
0185
0186 using SETraits = detail::SoftEqualTraits<value_type>;
0187 };
0188
0189
0190
0191
0192 template<class T>
0193 CELER_FUNCTION SoftEqual(T) -> SoftEqual<T>;
0194 template<class T>
0195 CELER_FUNCTION SoftEqual(T, T) -> SoftEqual<T>;
0196 template<class F>
0197 CELER_FUNCTION EqualOr(F&&) -> EqualOr<F>;
0198
0199
0200
0201
0202
0203
0204
0205 template<class RealType>
0206 CELER_CONSTEXPR_FUNCTION SoftEqual<RealType>::SoftEqual()
0207 : rel_{SETraits::rel_prec()}, abs_{SETraits::abs_thresh()}
0208 {
0209 }
0210
0211
0212
0213
0214
0215 template<class RealType>
0216 CELER_FUNCTION SoftEqual<RealType>::SoftEqual(value_type rel)
0217 : SoftEqual(rel, rel * (SETraits::abs_thresh() / SETraits::rel_prec()))
0218 {
0219 CELER_EXPECT(rel > 0);
0220 }
0221
0222
0223
0224
0225
0226
0227
0228
0229
0230
0231 template<class RealType>
0232 CELER_FUNCTION SoftEqual<RealType>::SoftEqual(value_type rel, value_type abs)
0233 : rel_{rel}, abs_{abs}
0234 {
0235 CELER_EXPECT(rel > 0);
0236 CELER_EXPECT(abs > 0);
0237 }
0238
0239
0240
0241
0242
0243 template<class RealType>
0244 CELER_FUNCTION bool
0245 SoftEqual<RealType>::operator()(value_type a, value_type b) const
0246 {
0247 real_type rel = rel_ * std::fmax(std::fabs(a), std::fabs(b));
0248 return std::fabs(a - b) < std::fmax(abs_, rel);
0249 }
0250
0251
0252
0253
0254
0255 template<class RealType>
0256 CELER_CONSTEXPR_FUNCTION SoftClose<RealType>::SoftClose()
0257 : abs_(SETraits::abs_thresh())
0258 {
0259 }
0260
0261
0262
0263
0264
0265
0266
0267 template<class RealType>
0268 CELER_FUNCTION SoftClose<RealType>::SoftClose(value_type abs) : abs_(abs)
0269 {
0270 CELER_EXPECT(abs > 0);
0271 }
0272
0273
0274
0275
0276
0277 template<class RealType>
0278 CELER_FUNCTION bool
0279 SoftClose<RealType>::operator()(value_type a, value_type b) const
0280 {
0281 return std::fabs(a - b) < abs_;
0282 }
0283
0284
0285
0286
0287
0288 template<class RealType>
0289 CELER_CONSTEXPR_FUNCTION SoftZero<RealType>::SoftZero()
0290 : abs_(SETraits::abs_thresh())
0291 {
0292 }
0293
0294
0295
0296
0297
0298
0299
0300 template<class RealType>
0301 CELER_FUNCTION SoftZero<RealType>::SoftZero(value_type abs) : abs_(abs)
0302 {
0303 CELER_EXPECT(abs > 0);
0304 }
0305
0306
0307
0308
0309
0310
0311
0312 template<class RealType>
0313 CELER_FUNCTION bool SoftZero<RealType>::operator()(value_type actual) const
0314 {
0315 return std::fabs(actual) < abs_;
0316 }
0317
0318
0319
0320 template<class RealType>
0321 inline CELER_FUNCTION bool soft_equal(RealType expected, RealType actual)
0322 {
0323 return SoftEqual<RealType>()(expected, actual);
0324 }
0325
0326
0327
0328 template<class RealType>
0329 inline CELER_FUNCTION bool
0330 soft_near(RealType expected, RealType actual, RealType rel_error)
0331 {
0332 return SoftEqual<RealType>(rel_error)(expected, actual);
0333 }
0334
0335
0336
0337 template<class RealType>
0338 inline CELER_FUNCTION bool soft_zero(RealType actual)
0339 {
0340 return SoftZero<RealType>()(actual);
0341 }
0342
0343
0344
0345 template<class RealType>
0346 inline CELER_FUNCTION bool soft_mod(RealType dividend, RealType divisor)
0347 {
0348 auto remainder = std::fmod(dividend, divisor);
0349
0350 SoftEqual<RealType> seq(detail::SoftEqualTraits<RealType>::rel_prec());
0351
0352 return seq(0, remainder) || seq(divisor, remainder);
0353 }
0354
0355
0356 }