Skip to content

Commit 91e192b

Browse files
committed
added logMath aliases. replaced log and exp calls everywhere with aliases. a couple of derivatives wrt log fixed....need to finish this part
1 parent 2b2aeaf commit 91e192b

18 files changed

+269
-398
lines changed

src/common/logMath.hpp

Lines changed: 124 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,124 @@
1+
2+
#pragma once
3+
4+
#include "macros.hpp"
5+
#include <cmath>
6+
7+
8+
9+
namespace hpcReact
10+
{
11+
12+
#if defined(__GLIBC__) // glibc provides exp10/exp10f as an extension
13+
#define HPCREACT_HAS_EXP10 1
14+
#else
15+
#define HPCREACT_HAS_EXP10 0
16+
#endif
17+
18+
enum class LogBase { e, ten };
19+
enum class MathMode { accurate, fast };
20+
21+
template< LogBase Base, MathMode Mode = MathMode::accurate >
22+
struct LogExp
23+
{
24+
25+
template<class T>
26+
HPCREACT_HOST_DEVICE static constexpr
27+
T ln10()
28+
{
29+
return T(2.3025850929940456840179914546843642L);
30+
}
31+
32+
template< typename T >
33+
HPCREACT_HOST_DEVICE
34+
static constexpr inline
35+
T dLogConst()
36+
{
37+
if constexpr ( Base == LogBase::e ) { return T(1.0); }
38+
else { return T(1.0/ln10()); }
39+
}
40+
41+
42+
43+
// ***** log function *********************************************************************************************
44+
HPCREACT_HOST_DEVICE
45+
static inline
46+
double log( double const x )
47+
{
48+
if constexpr ( Base == LogBase::e ) { return ::log(x); }
49+
else { return ::log10(x); }
50+
}
51+
52+
HPCREACT_HOST_DEVICE
53+
static inline
54+
float log( float const x )
55+
{
56+
#if defined(__CUDA_ARCH__)
57+
if constexpr ( Mode == MathMode::fast )
58+
{
59+
if constexpr ( Base == LogBase::e) { return __logf(x); }
60+
else { return __log10f(x); }
61+
}
62+
#endif
63+
if constexpr ( Base == LogBase::e ) { return ::logf(x); }
64+
else { return ::log10f(x); }
65+
}
66+
67+
// ***** exp function *********************************************************************************************
68+
HPCREACT_HOST_DEVICE
69+
static inline
70+
double exp( double const x )
71+
{
72+
if constexpr ( Base == LogBase::e)
73+
{
74+
return ::exp(x);
75+
}
76+
else
77+
{
78+
#if HPCREACT_HAS_EXP10
79+
return ::exp10(x);
80+
#else
81+
return ::exp( x * ln10<double>() ); ;
82+
#endif
83+
}
84+
}
85+
86+
87+
HPCREACT_HOST_DEVICE
88+
static inline
89+
float exp( float const x )
90+
{
91+
#if defined(__CUDA_ARCH__)
92+
if constexpr ( Mode == MathMode::fast )
93+
{
94+
if constexpr ( Base == LogBase::e ) { return __expf(x); }
95+
else { return __exp10f(x); }
96+
}
97+
#endif
98+
if constexpr ( Base == LogBase::e)
99+
{
100+
return ::expf(x);
101+
}
102+
else
103+
{
104+
#if HPCREACT_HAS_EXP10
105+
return ::exp10f(x);
106+
#else
107+
return ::expf( x * ln10<float>() );
108+
#endif
109+
}
110+
}
111+
}; // struct LogExp
112+
113+
114+
#if !defined(HPC_REACT_LOG_TYPE)
115+
#define HPC_REACT_LOG_TYPE 1
116+
#endif
117+
118+
#if HPC_REACT_LOG_TYPE == 0
119+
using logmath = LogExp< LogBase::e, MathMode::fast >;
120+
#elif HPC_REACT_LOG_TYPE == 1
121+
using logmath = LogExp< LogBase::ten, MathMode::fast >;
122+
#endif
123+
124+
} // namespace hpcReact

src/reactions/exampleSystems/unitTests/testGenericChainReactions.cpp

Lines changed: 1 addition & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -44,12 +44,7 @@ TEST( testChainGenericKineticReactions, computeReactionRatesTest_chainReactionPa
4444
{ { 0.05, 0.0, 0.0 },
4545
{ 0.0, 0.03, 0.0 },
4646
{ 0.0, 0.0, 0.02 } };
47-
computeReactionRatesTest< double, false >( serialAllKineticParams.kineticReactionsParameters(),
48-
initialSpeciesConcentration,
49-
surfaceArea, // No use. Just to pass something here
50-
expectedReactionRates,
51-
expectedReactionRatesDerivatives );
52-
computeReactionRatesTest< double, true >( serialAllKineticParams.kineticReactionsParameters(),
47+
computeReactionRatesTest< double >( serialAllKineticParams.kineticReactionsParameters(),
5348
initialSpeciesConcentration,
5449
surfaceArea, // No use. Just to pass something here
5550
expectedReactionRates,

src/reactions/exampleSystems/unitTests/testKineticReactions.cpp

Lines changed: 3 additions & 18 deletions
Original file line numberDiff line numberDiff line change
@@ -29,12 +29,7 @@ TEST( testKineticReactions, computeReactionRatesTest_simpleKineticTestRateParams
2929
double const expectedReactionRatesDerivatives[][5] =
3030
{ { 2.0, -0.5, 0.0, 0.0, 0.0 },
3131
{ 0.0, 0.0, 0.5, 0.25, 0.0 } };
32-
computeReactionRatesTest< double, false >( bulkGeneric::simpleKineticTestRateParams.kineticReactionsParameters(),
33-
initialSpeciesConcentration,
34-
surfaceArea, // No use. Just to pass something here
35-
expectedReactionRates,
36-
expectedReactionRatesDerivatives );
37-
computeReactionRatesTest< double, true >( bulkGeneric::simpleKineticTestRateParams.kineticReactionsParameters(),
32+
computeReactionRatesTest< double >( bulkGeneric::simpleKineticTestRateParams.kineticReactionsParameters(),
3833
initialSpeciesConcentration,
3934
surfaceArea, // No use. Just to pass something here
4035
expectedReactionRates,
@@ -52,12 +47,8 @@ TEST( testKineticReactions, computeSpeciesRatesTest_simpleKineticTestRateParams
5247
{ 0.0, 0.0, -0.5, -0.25, 0.0 },
5348
{ 0.0, 0.0, 1.0, 0.5, 0.0 } };
5449

55-
computeSpeciesRatesTest< double, false >( bulkGeneric::simpleKineticTestRateParams.kineticReactionsParameters(),
56-
initialSpeciesConcentration,
57-
expectedSpeciesRates,
58-
expectedSpeciesRatesDerivatives );
5950

60-
computeSpeciesRatesTest< double, true >( bulkGeneric::simpleKineticTestRateParams.kineticReactionsParameters(),
51+
computeSpeciesRatesTest< double >( bulkGeneric::simpleKineticTestRateParams.kineticReactionsParameters(),
6152
initialSpeciesConcentration,
6253
expectedSpeciesRates,
6354
expectedSpeciesRatesDerivatives );
@@ -70,18 +61,12 @@ TEST( testKineticReactions, testTimeStep )
7061
double const initialSpeciesConcentration[5] = { 1.0, 1.0e-16, 0.5, 1.0, 1.0e-16 };
7162
double const expectedSpeciesConcentrations[5] = { 3.92138293924124e-01, 3.03930853037938e-01, 5.05945480771998e-01, 7.02014627734060e-01, 5.95970744531880e-01 };
7263

73-
timeStepTest< double, false >( bulkGeneric::simpleKineticTestRateParams.kineticReactionsParameters(),
64+
timeStepTest< double >( bulkGeneric::simpleKineticTestRateParams.kineticReactionsParameters(),
7465
2.0,
7566
10,
7667
initialSpeciesConcentration,
7768
expectedSpeciesConcentrations );
7869

79-
// ln(c) as the primary variable results in a singular system.
80-
// timeStepTest< double, true >( simpleKineticTestRateParams,
81-
// 2.0,
82-
// 10,
83-
// initialSpeciesConcentration,
84-
// expectedSpeciesConcentrations );
8570
}
8671

8772
int main( int argc, char * * argv )

src/reactions/exampleSystems/unitTests/testMomasEasyCase.cpp

Lines changed: 6 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -53,11 +53,11 @@ void testMoMasAllEquilibriumHelper()
5353

5454
double const logInitialPrimarySpeciesConcentration[numPrimarySpecies] =
5555
{
56-
log( initialPrimarySpeciesConcentration[0] ),
57-
log( initialPrimarySpeciesConcentration[1] ),
58-
log( initialPrimarySpeciesConcentration[2] ),
59-
log( initialPrimarySpeciesConcentration[3] ),
60-
log( initialPrimarySpeciesConcentration[4] )
56+
logmath::log( initialPrimarySpeciesConcentration[0] ),
57+
logmath::log( initialPrimarySpeciesConcentration[1] ),
58+
logmath::log( initialPrimarySpeciesConcentration[2] ),
59+
logmath::log( initialPrimarySpeciesConcentration[3] ),
60+
logmath::log( initialPrimarySpeciesConcentration[4] )
6161
};
6262

6363
EquilibriumReactionsType::enforceEquilibrium_Aggregate( 0,
@@ -78,7 +78,7 @@ void testMoMasAllEquilibriumHelper()
7878

7979
for( int r=0; r<numPrimarySpecies; ++r )
8080
{
81-
EXPECT_NEAR( exp( logPrimarySpeciesConcentration[r] ), expectedPrimarySpeciesConcentrations[r], 1.0e-8 * expectedPrimarySpeciesConcentrations[r] );
81+
EXPECT_NEAR( logmath::exp( logPrimarySpeciesConcentration[r] ), expectedPrimarySpeciesConcentrations[r], 1.0e-8 * expectedPrimarySpeciesConcentrations[r] );
8282
}
8383

8484
}

src/reactions/exampleSystems/unitTests/testMomasMediumCase.cpp

Lines changed: 6 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -54,11 +54,11 @@ void testMoMasMediumEquilibriumHelper()
5454

5555
double const logInitialPrimarySpeciesConcentration[numPrimarySpecies] =
5656
{
57-
log( initialPrimarySpeciesConcentration[0] ),
58-
log( initialPrimarySpeciesConcentration[1] ),
59-
log( initialPrimarySpeciesConcentration[2] ),
60-
log( initialPrimarySpeciesConcentration[3] ),
61-
log( initialPrimarySpeciesConcentration[4] )
57+
logmath::log( initialPrimarySpeciesConcentration[0] ),
58+
logmath::log( initialPrimarySpeciesConcentration[1] ),
59+
logmath::log( initialPrimarySpeciesConcentration[2] ),
60+
logmath::log( initialPrimarySpeciesConcentration[3] ),
61+
logmath::log( initialPrimarySpeciesConcentration[4] )
6262
};
6363

6464
EquilibriumReactionsType::enforceEquilibrium_Aggregate( 0,
@@ -79,7 +79,7 @@ void testMoMasMediumEquilibriumHelper()
7979

8080
for( int r=0; r<numPrimarySpecies; ++r )
8181
{
82-
EXPECT_NEAR( exp( logPrimarySpeciesConcentration[r] ), expectedPrimarySpeciesConcentrations[r], 1.0e-8 * expectedPrimarySpeciesConcentrations[r] );
82+
EXPECT_NEAR( logmath::exp( logPrimarySpeciesConcentration[r] ), expectedPrimarySpeciesConcentrations[r], 1.0e-8 * expectedPrimarySpeciesConcentrations[r] );
8383
}
8484

8585

src/reactions/geochemistry/unitTests/testGeochemicalEquilibriumReactions.cpp

Lines changed: 8 additions & 8 deletions
Original file line numberDiff line numberDiff line change
@@ -122,13 +122,13 @@ TEST( testEquilibriumReactions, testcarbonateSystemAllEquilibrium2 )
122122

123123
double const logInitialPrimarySpeciesConcentration[numPrimarySpecies] =
124124
{
125-
log( initialPrimarySpeciesConcentration[0] ),
126-
log( initialPrimarySpeciesConcentration[1] ),
127-
log( initialPrimarySpeciesConcentration[2] ),
128-
log( initialPrimarySpeciesConcentration[3] ),
129-
log( initialPrimarySpeciesConcentration[4] ),
130-
log( initialPrimarySpeciesConcentration[5] ),
131-
log( initialPrimarySpeciesConcentration[6] )
125+
logmath::log( initialPrimarySpeciesConcentration[0] ),
126+
logmath::log( initialPrimarySpeciesConcentration[1] ),
127+
logmath::log( initialPrimarySpeciesConcentration[2] ),
128+
logmath::log( initialPrimarySpeciesConcentration[3] ),
129+
logmath::log( initialPrimarySpeciesConcentration[4] ),
130+
logmath::log( initialPrimarySpeciesConcentration[5] ),
131+
logmath::log( initialPrimarySpeciesConcentration[6] )
132132
};
133133

134134
double logPrimarySpeciesConcentration[numPrimarySpecies];
@@ -150,7 +150,7 @@ TEST( testEquilibriumReactions, testcarbonateSystemAllEquilibrium2 )
150150

151151
for( int r=0; r<numPrimarySpecies; ++r )
152152
{
153-
EXPECT_NEAR( exp( logPrimarySpeciesConcentration[r] ), expectedPrimarySpeciesConcentrations[r], 1.0e-8 * expectedPrimarySpeciesConcentrations[r] );
153+
EXPECT_NEAR( logmath::exp( logPrimarySpeciesConcentration[r] ), expectedPrimarySpeciesConcentrations[r], 1.0e-8 * expectedPrimarySpeciesConcentrations[r] );
154154
}
155155

156156

src/reactions/geochemistry/unitTests/testGeochemicalKineticReactions.cpp

Lines changed: 8 additions & 24 deletions
Original file line numberDiff line numberDiff line change
@@ -81,12 +81,7 @@ TEST( testKineticReactions, computeReactionRatesTest_carbonateSystemAllKinetic )
8181

8282
};
8383

84-
computeReactionRatesTest< double, false >( carbonateSystemAllKinetic.kineticReactionsParameters(),
85-
initialSpeciesConcentration,
86-
surfaceArea, // No use. Just to pass something here
87-
expectedReactionRates,
88-
expectedReactionRatesDerivatives );
89-
computeReactionRatesTest< double, true >( carbonateSystemAllKinetic.kineticReactionsParameters(),
84+
computeReactionRatesTest< double >( carbonateSystemAllKinetic.kineticReactionsParameters(),
9085
initialSpeciesConcentration,
9186
surfaceArea, // No use. Just to pass something here
9287
expectedReactionRates,
@@ -124,12 +119,7 @@ TEST( testKineticReactions, computeReactionRatesQuotientTest_carbonateSystem )
124119
{ 0, 0, 0, 0, 0, 0, 0, 0, 0, 3.0877659574468075e-03, -3.0877659574468075e-03, -2.9999999999999997e-02, 0, 0, 0, 0 }
125120
};
126121

127-
computeReactionRatesTest< double, false >( carbonateSystem.kineticReactionsParameters(),
128-
initialSpeciesConcentration,
129-
surfaceArea,
130-
expectedReactionRates,
131-
expectedReactionRatesDerivatives );
132-
computeReactionRatesTest< double, true >( carbonateSystem.kineticReactionsParameters(),
122+
computeReactionRatesTest< double >( carbonateSystem.kineticReactionsParameters(),
133123
initialSpeciesConcentration,
134124
surfaceArea,
135125
expectedReactionRates,
@@ -218,18 +208,12 @@ TEST( testKineticReactions, testTimeStep_carbonateSystemAllKinetic )
218208
1.072307827865370e+00 // Na+1
219209
};
220210

221-
timeStepTest< double, false >( carbonateSystemAllKinetic.kineticReactionsParameters(),
222-
10.0,
223-
10000,
224-
initialSpeciesConcentration,
225-
expectedSpeciesConcentrations );
226-
227-
// ln(c) as the primary variable results in a singular system.
228-
// timeStepTest< double, true >( simpleKineticTestRateParams,
229-
// 2.0,
230-
// 10,
231-
// initialSpeciesConcentration,
232-
// expectedSpeciesConcentrations );
211+
timeStepTest< double >( carbonateSystemAllKinetic.kineticReactionsParameters(),
212+
10.0,
213+
10000,
214+
initialSpeciesConcentration,
215+
expectedSpeciesConcentrations );
216+
233217
}
234218

235219
int main( int argc, char * * argv )

src/reactions/geochemistry/unitTests/testGeochemicalMixedReactions.cpp

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -50,7 +50,7 @@ TEST( testMixedReactions, testTimeStep_carbonateSystem )
5050
1.070434904554991 // Na+1
5151
};
5252

53-
timeStepTest< double, true >( carbonateSystem,
53+
timeStepTest< double >( carbonateSystem,
5454
1.0,
5555
10,
5656
initialAggregateSpeciesConcentration,

0 commit comments

Comments
 (0)