/*************************************************************************** * Copyright (c) 2023 Ondsel, Inc. * * * * This file is part of OndselSolver. * * * * See LICENSE file for details about copyright. * ***************************************************************************/ #pragma once #include #include "Array.h" #include "SimulationStoppingError.h" namespace MbD { template class FullVector : public Array { public: FullVector() : Array() {} FullVector(std::vector vec) : Array(vec) {} FullVector(size_t count) : Array(count) {} FullVector(size_t count, const T& value) : Array(count, value) {} FullVector(typename std::vector::iterator begin, typename std::vector::iterator end) : Array(begin, end) {} FullVector(std::initializer_list list) : Array{ list } {} double dot(std::shared_ptr> vec); void atiplusNumber(size_t i, T value); void atiminusNumber(size_t i, T value); double sumOfSquares() override; size_t numberOfElements() override; void zeroSelf() override; void atiplusFullVector(size_t i, std::shared_ptr> fullVec); void atiplusFullVectortimes(size_t i, std::shared_ptr> fullVec, T factor); void equalSelfPlusFullVectortimes(std::shared_ptr> fullVec, T factor); double maxMagnitude() override; void normalizeSelf(); double length(); virtual void conditionSelf(); virtual void conditionSelfWithTol(double tol); std::shared_ptr> clonesptr(); bool isIncreasing(); bool isIncreasingIfExceptionsAreLessThan(double tol); bool isDecreasingIfExceptionsAreLessThan(double tol); std::ostream& printOn(std::ostream& s) const override; }; template inline double FullVector::dot(std::shared_ptr> vec) { auto n = this->size(); double answer = 0.0; for (size_t i = 0; i < n; i++) { answer += this->at(i) * vec->at(i); } return answer; } template inline void FullVector::atiplusNumber(size_t i, T value) { this->at(i) += value; } template inline void FullVector::atiminusNumber(size_t i, T value) { this->at(i) -= value; } template<> inline double FullVector::sumOfSquares() { double sum = 0.0; for (size_t i = 0; i < this->size(); i++) { double element = this->at(i); sum += element * element; } return sum; } template inline double FullVector::sumOfSquares() { throw SimulationStoppingError("To be implemented."); return 0.0; } template inline size_t FullVector::numberOfElements() { return this->size(); } template<> inline void FullVector::zeroSelf() { for (size_t i = 0; i < this->size(); i++) { this->at(i) = 0.0; } } template inline void FullVector::zeroSelf() { throw SimulationStoppingError("To be implemented."); } template inline void FullVector::atiplusFullVector(size_t i1, std::shared_ptr> fullVec) { for (size_t ii = 0; ii < fullVec->size(); ii++) { auto i = i1 + ii; this->at(i) += fullVec->at(ii); } } template inline void FullVector::atiplusFullVectortimes(size_t i1, std::shared_ptr> fullVec, T factor) { for (size_t ii = 0; ii < fullVec->size(); ii++) { auto i = i1 + ii; this->at(i) += fullVec->at(ii) * factor; } } template inline void FullVector::equalSelfPlusFullVectortimes(std::shared_ptr> fullVec, T factor) { for (size_t i = 0; i < this->size(); i++) { this->atiplusNumber(i, fullVec->at(i) * factor); } } template<> inline double FullVector::maxMagnitude() { double max = 0.0; for (size_t i = 0; i < this->size(); i++) { double element = this->at(i); if (element < 0.0) element = -element; if (max < element) max = element; } return max; } template inline double FullVector::maxMagnitude() { throw SimulationStoppingError("To be implemented."); return 0.0; } template<> inline void FullVector::normalizeSelf() { double length = this->length(); if (length == 0.0) throw std::runtime_error("Cannot normalize a null vector."); this->magnifySelf(1.0 / length); } template inline double FullVector::length() { double ssq = 0.0; for (size_t i = 0; i < this->size(); i++) { double elem = this->at(i); ssq += elem * elem; } return std::sqrt(ssq); } template inline void FullVector::conditionSelf() { constexpr double epsilon = std::numeric_limits::epsilon(); double tol = this->maxMagnitude() * epsilon; this->conditionSelfWithTol(tol); } template<> inline void FullVector::conditionSelfWithTol(double tol) { for (size_t i = 0; i < this->size(); i++) { double element = this->at(i); if (element < 0.0) element = -element; if (element < tol) this->atiput(i, 0.0); } } template inline void FullVector::conditionSelfWithTol([[maybe_unused]] double tol) { assert(false && tol != tol); // clang++ flips out with warnings if you don't use 'tol' // but suppressing that warning breaks Visual Studio. return; // Visual Studio demands the unused return } template inline std::shared_ptr> FullVector::clonesptr() { //Return shallow copy of *this wrapped in shared_ptr throw SimulationStoppingError("To be implemented."); return std::make_shared>(*this); } template inline bool FullVector::isIncreasing() { return isIncreasingIfExceptionsAreLessThan(0.0); } template inline bool FullVector::isIncreasingIfExceptionsAreLessThan(double tol) { //"Test if elements are increasing." //"Ok if spoilers are less than tol." auto next = this->at(0); for (size_t i = 1; i < this->size(); i++) { auto previous = next; next = this->at(i); if (previous > next && (previous - tol > next)) return false; } return true; } template inline bool FullVector::isDecreasingIfExceptionsAreLessThan(double tol) { //"Test if elements are increasing." //"Ok if spoilers are less than tol." auto next = this->at(0); for (size_t i = 1; i < this->size(); i++) { auto previous = next; next = this->at(i); if (previous < next && (previous + tol < next)) return false; } return true; } template inline std::ostream& FullVector::printOn(std::ostream& s) const { s << "FullVec{"; s << this->at(0); for (size_t i = 1; i < this->size(); i++) { s << ", " << this->at(i); } s << "}"; return s; } }