This commit is contained in:
Jonathan Rampersad
2023-10-06 12:18:58 -04:00
parent 6b99a8d143
commit f10e84bf30
2 changed files with 319 additions and 188 deletions

View File

@ -16,7 +16,6 @@ namespace JRAMPERSAD
{
/**
* \brief Structure for options to be used when running one of the two genetic algorithms in a Function object
*
*/
struct GA_Options
{
@ -113,17 +112,17 @@ namespace JRAMPERSAD
});
}
template <int lrgst_expo> // Genetic Algorithm helper struct
// Genetic Algorithm helper struct
struct GA_Solution
{
unsigned short lrgst_expo;
double rank, x, y_val;
bool ranked;
GA_Solution() : rank(0), x(0), y_val(0), ranked(false) {}
GA_Solution(double Rank, double x_val, double y = 0) : rank(Rank), x(x_val), y_val(y), ranked(false) {}
GA_Solution() : lrgst_expo(0), rank(0), x(0), y_val(0) {}
GA_Solution(unsigned short Lrgst_expo, double Rank, double x_val, double y = 0) : lrgst_expo(Lrgst_expo), rank(Rank), x(x_val), y_val(y) {}
virtual ~GA_Solution() = default;
void fitness(const std::vector<int>& constants)
void fitness(const std::vector<int64_t>& constants)
{
double ans = 0;
for (int i = lrgst_expo; i >= 0; i--)
@ -137,153 +136,118 @@ namespace JRAMPERSAD
using namespace detail;
/**
* \brief A class representing an Exponential Function (e.g 2x^2 + 4x - 1),
* \tparam lrgst_expo The largest exponent in the function (e.g 2 means largest exponent is x^2)
* \brief class representing an Exponential Function (e.g 2x^2 + 4x - 1)
*/
template <int lrgst_expo>
class Function
{
private:
std::vector<int> constants;
const unsigned short lrgst_expo; /**< lrgst_expo The largest exponent in the function (e.g 2 means largest exponent is x^2) */
std::vector<int64_t> constants;
bool bInitialized;
void CanPerform() const { if (!bInitialized) throw std::logic_error("Function object not initialized fully! Please call .SetConstants() to initialize"); }
public:
// Speicialty function to get the real roots of a Quadratic Function without relying on a Genetic Algorithm to approximate
friend std::vector<double> QuadraticSolve(const Function<2>& f);
friend std::vector<double> QuadraticSolve(const Function& f);
public:
/**
* \brief Constructor for Function class
* \param constnts An array with the constants for the function (e.g 2, 1, 3 = 2x^2 + 1x - 3) size of array MUST be lrgst_expo + 1
* \param Lrgst_expo The largest exponent in the function (e.g 2 means largest exponent is x^2)
*/
Function(const std::vector<int>& constnts);
/**
* \brief Constructor for Function class
* \param constnts An array with the constants for the function (e.g 2, 1, 3 = 2x^2 + 1x - 3) size of array MUST be lrgst_expo + 1
*/
Function(std::vector<int>&& constnts);
Function(const Function& other) = default;
Function(Function&& other) noexcept = default;
Function(const unsigned short& Lrgst_expo) : lrgst_expo(Lrgst_expo), bInitialized(false)
{
if (lrgst_expo < 0)
throw std::logic_error("Function template argument must not be less than 0");
constants.reserve(Lrgst_expo);
}
/** \brief Destructor */
virtual ~Function();
/** \brief Copy Constructor */
Function(const Function& other) = default;
/** \brief Move Constructor */
Function(Function&& other) noexcept = default;
/** \brief Copy Assignment operator */
Function& operator=(const Function& other) = default;
/** \brief Move Assignment operator */
Function& operator=(Function&& other) noexcept = default;
// Operator function to display function object in a human readable format
friend std::ostream& operator<<(std::ostream& os, const Function<lrgst_expo> func)
{
if (lrgst_expo == 0)
{
os << func.constants[0];
return os;
}
/**
* \brief Sets the constants of the function
* \param constnts An array with the constants for the function (e.g 2, 1, 3 = 2x^2 + 1x - 3) size of array MUST be lrgst_expo + 1
*/
void SetConstants(const std::vector<int64_t>& constnts);
/**
* \brief Sets the constants of the function
* \param constnts An array with the constants for the function (e.g 2, 1, 3 = 2x^2 + 1x - 3) size of array MUST be lrgst_expo + 1
*/
void SetConstants(std::vector<int64_t>&& constnts);
if (func.constants[0] == 1)
os << "x";
else if (func.constants[0] == -1)
os << "-x";
else
os << func.constants[0] << "x";
friend std::ostream& operator<<(std::ostream& os, const Function func);
if (lrgst_expo != 1)
os << "^" << lrgst_expo;
friend Function operator+(const Function& f1, const Function& f2);
friend Function operator-(const Function& f1, const Function& f2);
for (int i = lrgst_expo - 1; i > 0; i--)
{
int n = func.constants[lrgst_expo - i];
if (n == 0) continue;
auto s = n > 0 ? " + " : " - ";
if (n != 1)
os << s << ABS(n) << "x";
else
os << s << "x";
if (i != 1)
os << "^" << i;
}
int n = func.constants[lrgst_expo];
if (n == 0) return os;
auto s = n > 0 ? " + " : " - ";
os << s;
os << ABS(n);
return os;
}
template<int e1, int e2, int r>
friend Function<r> operator+(const Function<e1>& f1, const Function<e2>& f2); // Operator to add two functions
template<int e1, int e2, int r>
friend Function<r> operator-(const Function<e1>& f1, const Function<e2>& f2); // Operator to subtract two functions
// Operators to multiply a function by a constant (Scaling it)
friend Function<lrgst_expo> operator*(const Function<lrgst_expo>& f, const int& c)
{
if (c == 1) return f;
if (c == 0) throw std::logic_error("Cannot multiply a function by 0");
std::vector<int> res;
for (auto& val : f.constants)
res.push_back(c * val);
return Function<lrgst_expo>(res);
}
Function<lrgst_expo>& operator*=(const int& c)
{
if (c == 1) return *this;
if (c == 0) throw std::logic_error("Cannot multiply a function by 0");
for (auto& val : this->constants)
val *= c;
return *this;
}
friend Function operator*(const Function& f, const int64_t& c);
Function& operator*=(const int64_t& c);
/**
* \brief Calculates the differential (dy/dx) of the function
* \returns a function representing the differential (dy/dx) of the calling function object
* \brief Calculates the differential (dy/dx) of the Function
* \returns a Function representing the differential (dy/dx) of the calling function object
*/
[[nodiscard("MATH::EXP::Function::differential() returns the differential, the calling object is not changed")]]
Function<lrgst_expo - 1> differential() const;
Function differential() const;
/**
* \brief Function that uses a genetic algorithm to find the approximate roots of the function
* \brief Uses a genetic algorithm to find the approximate roots of the function
* \param options GA_Options object specifying the options to run the algorithm
* \returns A vector containing a n number of approximate root values (n = sample_size as defined in options)
*/
[[nodiscard]] std::vector<double> get_real_roots(const GA_Options& options = GA_Options()) const;
/**
* \brief Function that solves for y when x = user value
* \brief Solves for y when x = user value
* \param x_val the X Value you'd like the function to use
* \returns the Y value the function returns based on the entered X value
*/
[[nodiscard]] double solve_y(const double& x_val) const noexcept;
[[nodiscard]] double solve_y(const double& x_val) const;
/**
* \brief Function that uses a genetic algorithm to find the values of x where y = user value
* \brief Uses a genetic algorithm to find the values of x where y = user value
* \param y_val The return value that you would like to find the approximate x values needed to solve when entered into the function
* \param options GA_Options object specifying the options to run the algorithm
* \returns A vector containing a n number of x values that cause the function to approximately equal the y_val (n = sample_size as defined in options)
*/
[[nodiscard]] std::vector<double> solve_x(const double& y_val, const GA_Options& options = GA_Options()) const;
/** \returns lrgst_expo */
[[nodiscard]] auto GetWhatIsTheLargestExponent() const { return lrgst_expo; }
};
/**
* \brief Uses the quadratic function to solve the roots of an entered quadratic equation
* \param f Quadratic function you'd like to find the roots of (Quadratic Function object is a Function<2> object
* \param f Quadratic function you'd like to find the roots of (Quadratic Function object is a Function object who's lrgst_expo value = 2
* \returns a vector containing the roots
*/
std::vector<double> QuadraticSolve(const Function<2>& f)
std::vector<double> QuadraticSolve(const Function& f)
{
try
{
if (f.lrgst_expo != 2) throw std::logic_error("Function f is not a quadratic function");
f.CanPerform();
}
catch (const std::exception& e)
{
throw e;
}
std::vector<double> res;
const int& a = f.constants[0];
const int& b = f.constants[1];
const int& c = f.constants[2];
const auto& a = f.constants[0];
const auto& b = f.constants[1];
const auto& c = f.constants[2];
const double sqr_val = static_cast<double>(POW(b, 2) - (4 * a * c));
@ -297,16 +261,114 @@ namespace JRAMPERSAD
return res;
}
template<int e1, int e2, int r = (e1 > e2 ? e1 : e2)>
Function<r> operator+(const Function<e1>& f1, const Function<e2>& f2)
Function::~Function()
{
std::vector<int> res;
constants.clear();
}
void Function::SetConstants(const std::vector<int64_t>& constnts)
{
if (constnts.size() != lrgst_expo + 1)
throw std::logic_error("Function<n> must be created with (n+1) integers in vector object");
if (constnts[0] == 0)
throw std::logic_error("First value should not be 0");
constants = constnts;
bInitialized = true;
}
void Function::SetConstants(std::vector<int64_t>&& constnts)
{
if (constnts.size() != lrgst_expo + 1)
throw std::logic_error("Function<n> must be created with (n+1) integers in vector object");
if (constnts[0] == 0)
throw std::logic_error("First value should not be 0");
constants = std::move(constnts);
bInitialized = true;
}
/** Operator function to display function object in a human readable format */
std::ostream& operator<<(std::ostream& os, const Function func)
{
try
{
func.CanPerform();
}
catch (const std::exception& e)
{
throw e;
}
if (func.lrgst_expo == 0)
{
os << func.constants[0];
return os;
}
if (func.constants[0] == 1)
os << "x";
else if (func.constants[0] == -1)
os << "-x";
else
os << func.constants[0] << "x";
if (func.lrgst_expo != 1)
os << "^" << func.lrgst_expo;
for (auto i = func.lrgst_expo - 1; i > 0; i--)
{
auto n = func.constants[func.lrgst_expo - i];
if (n == 0) continue;
auto s = n > 0 ? " + " : " - ";
if (n != 1)
os << s << ABS(n) << "x";
else
os << s << "x";
if (i != 1)
os << "^" << i;
}
auto n = func.constants[func.lrgst_expo];
if (n == 0) return os;
auto s = n > 0 ? " + " : " - ";
os << s;
os << ABS(n);
return os;
}
/** Operator to add two functions */
Function operator+(const Function& f1, const Function& f2)
{
try
{
f1.CanPerform();
f2.CanPerform();
}
catch (const std::exception& e)
{
throw e;
}
auto e1 = f1.lrgst_expo;
auto e2 = f2.lrgst_expo;
auto r = e1 > e2 ? e1 : e2;
std::vector<int64_t> res;
if (e1 > e2)
{
for (auto& val : f1.constants)
res.push_back(val);
int i = e1 - e2;
auto i = e1 - e2;
for (auto& val : f2.constants)
{
res[i] += val;
@ -326,19 +388,35 @@ namespace JRAMPERSAD
}
}
return Function<r>{res};
Function f(r);
f.SetConstants(res);
return f;
}
template<int e1, int e2, int r = (e1 > e2 ? e1 : e2)>
Function<r> operator-(const Function<e1>& f1, const Function<e2>& f2)
/** Operator to subtract two functions */
Function operator-(const Function& f1, const Function& f2)
{
std::vector<int> res;
try
{
f1.CanPerform();
f2.CanPerform();
}
catch (const std::exception& e)
{
throw e;
}
auto e1 = f1.lrgst_expo;
auto e2 = f2.lrgst_expo;
auto r = e1 > e2 ? e1 : e2;
std::vector<int64_t> res;
if (e1 > e2)
{
for (auto& val : f1.constants)
res.push_back(val);
int i = e1 - e2;
auto i = e1 - e2;
for (auto& val : f2.constants)
{
res[i] -= val;
@ -362,78 +440,109 @@ namespace JRAMPERSAD
}
}
return Function<r>{res};
Function f(r);
f.SetConstants(res);
return f;
}
template <int lrgst_expo>
Function<lrgst_expo>::Function(const std::vector<int>& constnts)
/** Operator to multiply a function by a constant (Scaling it) */
Function operator*(const Function& f, const int64_t& c)
{
if (lrgst_expo < 0)
throw std::logic_error("Function template argument must not be less than 0");
if (constnts.size() != lrgst_expo + 1)
throw std::logic_error("Function<n> must be created with (n+1) integers in vector object");
if (constnts[0] == 0)
throw std::logic_error("First value should not be 0");
constants = constnts;
try
{
f.CanPerform();
}
catch (const std::exception& e)
{
throw e;
}
template<int lrgst_expo>
Function<lrgst_expo>::Function(std::vector<int>&& constnts)
{
if (lrgst_expo < 0)
throw std::logic_error("Function template argument must not be less than 0");
if (c == 1) return f;
if (c == 0) throw std::logic_error("Cannot multiply a function by 0");
if (constnts.size() != lrgst_expo + 1)
throw std::logic_error("Function<n> must be created with (n+1) integers in vector object");
std::vector<int64_t> res;
for (auto& val : f.constants)
res.push_back(c * val);
if (constnts[0] == 0)
throw std::logic_error("First value should not be 0");
Function f_res(f.lrgst_expo);
f_res.SetConstants(res);
constants = std::move(constnts);
return f_res;
}
template <int lrgst_expo>
Function<lrgst_expo>::~Function()
/** Operator to multiply a function by a constant (Scaling it) */
Function& Function::operator*=(const int64_t& c)
{
constants.clear();
try
{
this->CanPerform();
}
catch (const std::exception& e)
{
throw e;
}
template <int lrgst_expo>
Function<lrgst_expo - 1> Function<lrgst_expo>::differential() const
if (c == 1) return *this;
if (c == 0) throw std::logic_error("Cannot multiply a function by 0");
for (auto& val : this->constants)
val *= c;
return *this;
}
Function Function::differential() const
{
try
{
this->CanPerform();
}
catch (const std::exception& e)
{
throw e;
}
if (lrgst_expo == 0)
throw std::logic_error("Cannot differentiate a number (Function<0>)");
std::vector<int> result;
std::vector<int64_t> result;
for (int i = 0; i < lrgst_expo; i++)
{
result.push_back(constants[i] * (lrgst_expo - i));
}
return Function<lrgst_expo - 1>{result};
Function f{ (unsigned short)(lrgst_expo - 1) };
f.SetConstants(result);
return f;
}
template<int lrgst_expo>
std::vector<double> Function<lrgst_expo>::get_real_roots(const GA_Options& options) const
std::vector<double> Function::get_real_roots(const GA_Options& options) const
{
try
{
this->CanPerform();
}
catch (const std::exception& e)
{
throw e;
}
// Create initial random solutions
std::random_device device;
std::uniform_real_distribution<double> unif(options.min_range, options.max_range);
std::vector<GA_Solution<lrgst_expo>> solutions;
std::vector<GA_Solution> solutions;
solutions.resize(options.data_size);
for (unsigned int i = 0; i < options.sample_size; i++)
solutions[i] = (GA_Solution<lrgst_expo>{0, unif(device)});
solutions[i] = (GA_Solution{lrgst_expo, 0, unif(device)});
float timer{ 0 };
for (unsigned int count = 0; count < options.num_of_generations; count++)
{
std::generate(std::execution::par, solutions.begin() + options.sample_size, solutions.end(), [&unif, &device]() {
return GA_Solution<lrgst_expo>{0, unif(device)};
std::generate(std::execution::par, solutions.begin() + options.sample_size, solutions.end(), [this, &unif, &device]() {
return GA_Solution{lrgst_expo, 0, unif(device)};
});
// Run our fitness function
@ -446,7 +555,7 @@ namespace JRAMPERSAD
});
// Take top solutions
std::vector<GA_Solution<lrgst_expo>> sample;
std::vector<GA_Solution> sample;
std::copy(
solutions.begin(),
solutions.begin() + options.sample_size,
@ -495,40 +604,49 @@ namespace JRAMPERSAD
return ans;
}
template<int lrgst_expo>
double Function<lrgst_expo>::solve_y(const double& x_val) const noexcept
double Function::solve_y(const double& x_val) const
{
std::vector<bool> exceptions;
for (int i : constants)
exceptions.push_back(i != 0);
try
{
this->CanPerform();
}
catch (const std::exception& e)
{
throw e;
}
double ans{ 0 };
for (int i = lrgst_expo; i >= 0; i--)
{
if (exceptions[i])
ans += constants[i] * POW(x_val, (lrgst_expo - i));
}
return ans;
}
template<int lrgst_expo>
inline std::vector<double> Function<lrgst_expo>::solve_x(const double& y_val, const GA_Options& options) const
inline std::vector<double> Function::solve_x(const double& y_val, const GA_Options& options) const
{
try
{
this->CanPerform();
}
catch (const std::exception& e)
{
throw e;
}
// Create initial random solutions
std::random_device device;
std::uniform_real_distribution<double> unif(options.min_range, options.max_range);
std::vector<GA_Solution<lrgst_expo>> solutions;
std::vector<GA_Solution> solutions;
solutions.resize(options.data_size);
for (unsigned int i = 0; i < options.sample_size; i++)
solutions[i] = (GA_Solution<lrgst_expo>{0, unif(device), y_val});
solutions[i] = (GA_Solution{lrgst_expo, 0, unif(device), y_val});
for (unsigned int count = 0; count < options.num_of_generations; count++)
{
std::generate(std::execution::par, solutions.begin() + options.sample_size, solutions.end(), [&unif, &device, &y_val]() {
return GA_Solution<lrgst_expo>{0, unif(device), y_val};
std::generate(std::execution::par, solutions.begin() + options.sample_size, solutions.end(), [this, &unif, &device, &y_val]() {
return GA_Solution{lrgst_expo, 0, unif(device), y_val};
});
@ -542,7 +660,7 @@ namespace JRAMPERSAD
});
// Take top solutions
std::vector<GA_Solution<lrgst_expo>> sample;
std::vector<GA_Solution> sample;
std::copy(
solutions.begin(),
solutions.begin() + options.sample_size,
@ -593,4 +711,7 @@ namespace JRAMPERSAD
}
}
#define INITIALIZE_EXPO_FUNCTION(func, ...) \
func.SetConstants(__VA_ARGS__)
#endif // !JONATHAN_RAMPERSAD_EXPONENTIAL_H_

View File

@ -7,13 +7,11 @@
using namespace JRAMPERSAD;
template <int n>
using Function = EXPONENTIAL::Function<n>;
using EXPONENTIAL::Function;
typedef TIMER::Timer timer;
template<int exp>
void CalcRoots(std::mutex& m, const Function<exp>& func, EXPONENTIAL::GA_Options options)
void CalcRoots(std::mutex& m, const Function& func, EXPONENTIAL::GA_Options options)
{
m.lock();
std::cout << "Starting calculation...\n";
@ -33,8 +31,7 @@ void CalcRoots(std::mutex& m, const Function<exp>& func, EXPONENTIAL::GA_Options
m.unlock();
}
template<int exp>
void SolveX(std::mutex& m, const Function<exp>& func, EXPONENTIAL::GA_Options options, const double& y)
void SolveX(std::mutex& m, const Function& func, EXPONENTIAL::GA_Options options, const double& y)
{
timer t;
auto res = func.solve_x(y, options);
@ -52,27 +49,37 @@ void SolveX(std::mutex& m, const Function<exp>& func, EXPONENTIAL::GA_Options op
int main()
{
std::vector<int> vec{ 1, 5, 4 };
Function<2> f{ vec };
Function<3> g{ { 1, -6, 11, -6 } };
std::vector<int64_t> vec{ 1, 5, 4 };
Function f{2};
INITIALIZE_EXPO_FUNCTION(f, vec);
Function g{3};
INITIALIZE_EXPO_FUNCTION(g, { 1, -6, 11, -6 });
EXPONENTIAL::GA_Options options;
options.mutation_percentage = 0.005;
options.num_of_generations = 10;
options.sample_size = 1000;
options.data_size = 100000;
options.min_range = 4.9;
options.max_range = 5;
options.num_of_generations = 1;
options.sample_size = 1;
options.data_size = 2;
options.min_range = 0.13;
options.max_range = 0.14;
auto res = (f + g).get_real_roots(options);
std::for_each(res.begin(), res.end(),
[](const auto& val) {
std::cout << "x:" << val << '\n';
});
std::cout << (f + g) << " when x = 0.13056\n" << (f + g).solve_y(0.13056);
std::mutex m;
std::thread th(CalcRoots<3>, std::ref(m), std::cref(g), options);
//std::thread th1(SolveX<3>, std::ref(m), std::cref(g), options, 5);
//std::thread th2(SolveX<3>, std::ref(m), std::cref(g), options, 23);
//std::thread th(CalcRoots, std::ref(m), std::cref(g), options);
//std::thread th1(SolveX, std::ref(m), std::cref(g), options, 5);
//std::thread th2(SolveX, std::ref(m), std::cref(g), options, 23);
//CalcRoots<3>(m, g);
//CalcRoots(m, g);
m.lock();
std::cout << g << " when x = 4.961015\n" << "y = " << g.solve_y(4.961015) << "\n\n";
//std::cout << g << " when x = 4.961015\n" << "y = " << g.solve_y(4.961015) << "\n\n";
//std::cout << g << " when x = 4.30891\n" << "y = " << g.solve_y(4.30891) << "\n\n";
//std::cout << g << " when x = 2\n" << "y = " << g.solve_y(2) << "\n\n";
//std::cout << g << " when x = 3\n" << "y = " << g.solve_y(3) << "\n\n";
@ -82,10 +89,13 @@ int main()
//std::cout << "Calculating Roots for function f(x) = " << g << '\n';
//std::cout << "The y-intercept of the function f(x) is " << g.solve_y(0) << '\n';
std::cout << "dy/dx of f(x) is " << g.differential() << '\n';
//std::cout << "dy/dx of f(x) is " << g.differential() << '\n';
//std::cout << "f(x) = " << f << std::endl;
//std::cout << "g(x) = " << g << std::endl;
//std::cout << "f(x) + g(x) = " << f + g << std::endl;
m.unlock();
th.join();
//th.join();
//th1.join();
//th2.join();
return 0;