emath.cpp

Go to the documentation of this file.

 
#include "emath.h"
void emth::Arithmetic::get_divisors(int n, std::vector<int>& v) noexcept {
    int limit(std::sqrt(n));
    for (int i = 1; i <= limit; i++) {
        if (std::fmod(n, i) == 0) {
            v.push_back(i);
            if (i != n / i) {
                v.push_back(n / i);
            }
        }
    }
    return;
}
emth::Monomial::Monomial() noexcept : coeff(0), degree(0) {}
emth::Monomial::Monomial(const double cf, const int dgr) noexcept : coeff(cf), degree(dgr) {}
emth::Monomial::Monomial(const Monomial& ogn) noexcept : coeff(ogn.coeff), degree(ogn.degree) {}
emth::Monomial::Monomial(Monomial&& ogn) noexcept : coeff(std::move(ogn.coeff)), degree(std::move(ogn.degree)) {}
emth::Monomial& emth::Monomial::operator =(const emth::Monomial& ogn) noexcept {Monomial copy(ogn); this->swap(copy); return *this;}
emth::Monomial& emth::Monomial::operator =(emth::Monomial&& ogn) noexcept {this->swap(ogn); return *this;}
bool emth::Monomial::operator ==(const emth::Monomial& ogn) const noexcept {return ((this->coeff == ogn.coeff) && (this->degree == ogn.degree));}
bool emth::Monomial::operator !=(const emth::Monomial& ogn) const noexcept {return (!(*this == ogn));}
bool emth::Monomial::operator >(const emth::Monomial& ogn) const noexcept {return (this->degree > ogn.degree);};
bool emth::Monomial::operator <(const emth::Monomial& ogn) const noexcept {return (this->degree < ogn.degree);};
bool emth::Monomial::operator >=(const emth::Monomial& ogn) const noexcept {return (this->degree >= ogn.degree);};
bool emth::Monomial::operator <=(const emth::Monomial& ogn) const noexcept {return (this->degree <= ogn.degree);};
emth::Monomial emth::Monomial::operator +(const emth::Monomial& ogn) const noexcept {return emth::Monomial(*this) += ogn;}
emth::Monomial emth::Monomial::operator -(const emth::Monomial& ogn) const noexcept {return emth::Monomial(*this) -= ogn;}
emth::Monomial emth::Monomial::operator *(const emth::Monomial& ogn) const noexcept {return emth::Monomial(*this) *= ogn;}
emth::Monomial emth::Monomial::operator /(const emth::Monomial& ogn) const noexcept {return emth::Monomial(*this) /= ogn;}
emth::Monomial emth::Monomial::operator %(const emth::Monomial& ogn) const noexcept {return emth::Monomial(*this) %= ogn;}
emth::Monomial& emth::Monomial::operator +=(const emth::Monomial& ogn) noexcept {
    if ((this->degree == ogn.degree) || this->coeff == 0) {
        if (this->coeff == 0) {this->degree = ogn.degree;}
        this->coeff += ogn.coeff;
    }
    return *this;
}
emth::Monomial& emth::Monomial::operator -=(const emth::Monomial& ogn) noexcept {
    if ((this->degree == ogn.degree) || this->coeff == 0) {
        if (this->coeff == 0) {this->degree = ogn.degree;}
        this->coeff -= ogn.coeff;
    }
    return *this;
}
emth::Monomial& emth::Monomial::operator *=(const emth::Monomial& ogn) noexcept {
    this->coeff *= ogn.coeff;
    this->degree += ogn.degree;
    return *this;
}
emth::Monomial& emth::Monomial::operator /=(const emth::Monomial& ogn) noexcept {
    if (ogn.coeff != 0) {
        this->coeff /= ogn.coeff;
        this->degree -= ogn.degree;
    }
    return *this;
}
emth::Monomial& emth::Monomial::operator %=(const emth::Monomial& ogn) noexcept {
    if (ogn.coeff != 0) {
        this->coeff = std::fmod(this->coeff, ogn.coeff);
        this->degree -= ogn.degree;
    }
    return *this;
}
std::stringstream& emth::operator<<(std::stringstream& ss, const emth::Monomial& m) {
    if (m.coeff != 1.0 || m.degree == 0) {ss << m.coeff;} //if (m.coeff > 0) {ss << "+ ";} else {ss << "- ";}
    switch (m.degree) {
        case 0: break;
        case 1: ss << "x"; break;
        default: ss << "x^" << m.degree; break;
    } return ss;
}
std::ostream& emth::operator<<(std::ostream& out, const emth::Monomial& m) {
    if (m.coeff != 1.0 || m.degree == 0) {out << m.coeff;} //if (m.coeff > 0) {out << "+ ";} else {out << "- ";}
    switch (m.degree) {
        case 0: break;
        case 1: out << "x"; break;
        default: out << "x^" << m.degree; break;
    } return out;
}
double emth::Monomial::get_coeff() const noexcept {return this->coeff;}
int emth::Monomial::get_degree() const noexcept {return this->degree;}
void emth::Monomial::set_coeff(double cf) noexcept {this->coeff = cf; return;}
void emth::Monomial::set_degree(int dgr) noexcept {this->degree = dgr; return;}
std::string emth::Monomial::get_expression() const noexcept {return (std::stringstream() << *this).str();}
double emth::Monomial::get_value(const double& x) const noexcept {return (this->coeff * std::pow(x, this->degree));};
void emth::Monomial::swap(emth::Monomial& ogn) noexcept {
    std::swap(this->coeff, ogn.coeff);
    std::swap(this->degree, ogn.degree);
    return;
}
emth::Monomial::~Monomial() = default;
emth::Polynomial::Polynomial() noexcept {this->monomials[0] = emth::Monomial(0, 0);}
emth::Polynomial::Polynomial(const std::map<int, emth::Monomial> mns) noexcept : monomials(mns) {}
emth::Polynomial::Polynomial(const std::initializer_list<emth::Monomial> mns) noexcept {for(const emth::Monomial& mn: mns) {this->monomials.emplace(mn.get_degree(), mn);}}
emth::Polynomial::Polynomial(const std::vector<emth::Monomial> mns) noexcept {for(const emth::Monomial& mn: mns) {this->monomials.emplace(mn.get_degree(), mn);}}
emth::Polynomial::Polynomial(const emth::Polynomial& ogn) noexcept : monomials(ogn.monomials) {} 
emth::Polynomial::Polynomial(emth::Polynomial&& ogn) noexcept : monomials(std::move(ogn.monomials)) {}
emth::Polynomial& emth::Polynomial::operator =(const emth::Polynomial& ogn) noexcept {emth::Polynomial copy(ogn); this->swap(copy); return *this;}
emth::Polynomial& emth::Polynomial::operator =(emth::Polynomial&& ogn) noexcept {this->swap(ogn); return *this;}
bool emth::Polynomial::operator ==(const emth::Polynomial& ogn) const noexcept {return ((this->get_degree() == ogn.get_degree()) && (this->monomials == ogn.monomials));}
bool emth::Polynomial::operator !=(const emth::Polynomial& ogn) const noexcept {return (!(*this == ogn));}
bool emth::Polynomial::operator >(const emth::Polynomial& ogn) const noexcept {return (this->get_degree() > ogn.get_degree());};
bool emth::Polynomial::operator <(const emth::Polynomial& ogn) const noexcept {return (this->get_degree() < ogn.get_degree());};
bool emth::Polynomial::operator >=(const emth::Polynomial& ogn) const noexcept {return (this->get_degree() >= ogn.get_degree());};
bool emth::Polynomial::operator <=(const emth::Polynomial& ogn) const noexcept {return (this->get_degree() <= ogn.get_degree());};
emth::Polynomial emth::Polynomial::operator +(const emth::Polynomial& ogn) const noexcept {return emth::Polynomial(*this) += ogn;}
emth::Polynomial emth::Polynomial::operator -(const emth::Polynomial& ogn) const noexcept {return emth::Polynomial(*this) -= ogn;}
emth::Polynomial emth::Polynomial::operator *(const emth::Polynomial& ogn) const noexcept {return emth::Polynomial(*this) *= ogn;}
emth::Polynomial emth::Polynomial::operator /(const emth::Polynomial& ogn) const noexcept {return emth::Polynomial(*this) /= ogn;}
emth::Polynomial emth::Polynomial::operator %(const emth::Polynomial& ogn) const noexcept {return emth::Polynomial(*this) %= ogn;}
emth::Polynomial& emth::Polynomial::operator +=(const emth::Polynomial& ogn) noexcept {
    for(const std::pair<const int, emth::Monomial>& kv : ogn.monomials) {this->monomials[kv.first] += kv.second;}
    return this->redux();
}
emth::Polynomial& emth::Polynomial::operator -=(const emth::Polynomial& ogn) noexcept {
    for(const std::pair<const int, emth::Monomial>& kv : ogn.monomials) {this->monomials[kv.first] -= kv.second;}
    return this->redux();
}
emth::Polynomial& emth::Polynomial::operator *=(const emth::Polynomial& ogn) noexcept {
    std::map<int, emth::Monomial> rs;
    for(const std::pair<const int, emth::Monomial>& kv : this->monomials) {
        for(const std::pair<const int, emth::Monomial>& kkvv : ogn.monomials) {
            rs[kv.first + kkvv.first] += kv.second * kkvv.second;
        }
    }
    this->monomials = std::move(rs);
    return this->redux();
}
emth::Polynomial& emth::Polynomial::operator /=(const emth::Polynomial& ogn) noexcept {
    if (!ogn.monomials.empty()) {
        emth::Polynomial quotient;
        while (!this->monomials.empty() && this->monomials.crbegin()->second.get_degree() >= ogn.monomials.crbegin()->second.get_degree()) {
            emth::Monomial mon = this->monomials.crbegin()->second / ogn.monomials.crbegin()->second;
            quotient.monomials[mon.get_degree()] = mon;
            *this -= (ogn * emth::Polynomial({mon}));
        }
        this->monomials = quotient.monomials;
    }
    return this->redux();
}
emth::Polynomial& emth::Polynomial::operator %=(const emth::Polynomial& ogn) noexcept {
    if (!ogn.monomials.empty()) {
        while (!this->monomials.empty() && this->monomials.crbegin()->second.get_degree() >= ogn.monomials.crbegin()->second.get_degree()) {
            emth::Monomial mon = this->monomials.crbegin()->second / ogn.monomials.crbegin()->second;
            *this -= (ogn * emth::Polynomial({mon}));
        }
    }
    return this->redux();
}
std::stringstream& emth::operator<<(std::stringstream& ss, const emth::Polynomial& p) noexcept {
    std::reverse_iterator<std::map<const int, emth::Monomial>::const_iterator> it = p.monomials.crbegin(); ss << it++->second;
    for(std::reverse_iterator<std::map<const int, emth::Monomial>::const_iterator> _it = it; _it != p.monomials.crend(); _it++) {
        if (_it->second.get_coeff() > 0) {ss << " + "; ss << _it->second;} else {ss << " - "; ss << _it->second.get_expression().erase(0, 1);}
    }
    return ss;
}
std::ostream& emth::operator<<(std::ostream& out, const emth::Polynomial& p) noexcept {
    std::reverse_iterator<std::map<const int, emth::Monomial>::const_iterator> it = p.monomials.crbegin(); out << it++->second;
    for(std::reverse_iterator<std::map<const int, emth::Monomial>::const_iterator> _it = it; _it != p.monomials.crend(); _it++) {
        if (_it->second.get_coeff() > 0) {out << " + "; out << _it->second;} else {out << " - "; out << _it->second.get_expression().erase(0, 1);}
    }
    return out;
}
int emth::Polynomial::get_degree() const noexcept {return this->monomials.crbegin()->first;}
std::map<int, emth::Monomial> emth::Polynomial::get_monomials() const noexcept {return this->monomials;}
bool emth::Polynomial::is_empty() const noexcept {return this->monomials.empty();}
bool emth::Polynomial::push_monomial(const emth::Monomial& m) noexcept {return this->push_monomial(emth::Monomial(m));};
bool emth::Polynomial::push_monomial(emth::Monomial&& m) noexcept {
    if (m.get_coeff() == 0) {return false;}
    this->monomials[m.get_degree()] += m;
    this->redux();
    return true;
};
std::string emth::Polynomial::get_expression() const noexcept {return (std::stringstream() << *this).str();}
double emth::Polynomial::get_value(const double& x) const noexcept {
    double y = 0;
    for(const std::pair<const int, emth::Monomial>& kv: this->monomials) {y += kv.second.get_value(x);}
    return y;
};
emth::Polynomial emth::Polynomial::get_derivative() const noexcept {
    std::map<int, emth::Monomial> fdx;
    for(const std::pair<const int, emth::Monomial>& kv: this->monomials) {
        if (kv.second.get_degree() != 0) {fdx[kv.second.get_degree() - 1] = emth::Monomial(kv.second.get_coeff() * kv.second.get_degree(), kv.second.get_degree() - 1);}
    }
    return emth::Polynomial(std::move(fdx)).redux();
};
emth::Polynomial emth::Polynomial::get_integral() const noexcept {
    std::map<int, emth::Monomial> sdx;
    for(const std::pair<const int, emth::Monomial>& kv: this->monomials) {
        sdx[kv.second.get_degree() + 1] = emth::Monomial(kv.second.get_coeff() / (kv.second.get_degree() + 1), kv.second.get_degree() + 1);
    }
    return emth::Polynomial(std::move(sdx)).redux();
};
std::vector<std::complex<double>> emth::Polynomial::get_roots() const noexcept {
    Eigen::VectorXd cfs(this->get_degree()+1); for(const std::pair<const int, emth::Monomial>& kv: this->monomials) {cfs[kv.first] = kv.second.get_coeff();}
    const Eigen::PolynomialSolver<double, Eigen::Dynamic>::RootsType rts =  Eigen::PolynomialSolver<double, Eigen::Dynamic>(cfs).roots();
    std::vector<std::complex<double>> roots; for(const auto& rt: rts) {roots.push_back(std::complex(rt.real(), rt.imag()));};
    return roots;
}
void emth::Polynomial::swap(emth::Polynomial& ogn) noexcept {
    std::swap(this->monomials, ogn.monomials);
    return;
}
emth::Polynomial& emth::Polynomial::redux() noexcept {
    std::vector<int> ks;
    for(const std::pair<const int, emth::Monomial>& kv: this->monomials) {
        if (kv.second.get_coeff() == 0) {ks.push_back(kv.first);}
    }
    for(int& k: ks) {this->monomials.erase(k);}
    return *this;
}
emth::Polynomial::~Polynomial() = default;