tsgCoreOneDimensional.hpp

Go to the documentation of this file.

/*
 * Copyright (c) 2017, Miroslav Stoyanov
 *
 * This file is part of
 * Toolkit for Adaptive Stochastic Modeling And Non-Intrusive ApproximatioN: TASMANIAN
 *
 * Redistribution and use in source and binary forms, with or without modification, are permitted provided that the following conditions are met:
 *
 * 1. Redistributions of source code must retain the above copyright notice, this list of conditions and the following disclaimer.
 *
 * 2. Redistributions in binary form must reproduce the above copyright notice, this list of conditions
 *    and the following disclaimer in the documentation and/or other materials provided with the distribution.
 *
 * 3. Neither the name of the copyright holder nor the names of its contributors may be used to endorse
 *    or promote products derived from this software without specific prior written permission.
 *
 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES,
 * INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED.
 * IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY,
 * OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA,
 * OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY,
 * OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
 *
 * UT-BATTELLE, LLC AND THE UNITED STATES GOVERNMENT MAKE NO REPRESENTATIONS AND DISCLAIM ALL WARRANTIES, BOTH EXPRESSED AND IMPLIED.
 * THERE ARE NO EXPRESS OR IMPLIED WARRANTIES OF MERCHANTABILITY OR FITNESS FOR A PARTICULAR PURPOSE, OR THAT THE USE OF THE SOFTWARE WILL NOT INFRINGE ANY PATENT,
 * COPYRIGHT, TRADEMARK, OR OTHER PROPRIETARY RIGHTS, OR THAT THE SOFTWARE WILL ACCOMPLISH THE INTENDED RESULTS OR THAT THE SOFTWARE OR ITS USE WILL NOT RESULT IN INJURY OR DAMAGE.
 * THE USER ASSUMES RESPONSIBILITY FOR ALL LIABILITIES, PENALTIES, FINES, CLAIMS, CAUSES OF ACTION, AND COSTS AND EXPENSES, CAUSED BY, RESULTING FROM OR ARISING OUT OF,
 * IN WHOLE OR IN PART THE USE, STORAGE OR DISPOSAL OF THE SOFTWARE.
 */
 
#ifndef __TSG_CORE_ONE_DIMENSIONAL_HPP
#define __TSG_CORE_ONE_DIMENSIONAL_HPP
 
#include "tsgIOHelpers.hpp"
#include "tsgLinearSolvers.hpp"
#include "tsgSequenceOptimizer.hpp"
 
namespace TasGrid{
 
class CustomTabulated{
public:
    CustomTabulated() : num_levels(0){}
    template<typename iomode> CustomTabulated(std::istream &is, iomode) : CustomTabulated(){
        if (std::is_same<iomode, IO::mode_ascii_type>::value)
            read<mode_ascii>(is); else read<mode_binary>(is);
    }
    CustomTabulated(std::vector<int> &&cnum_nodes, std::vector<int> &&cprecision,
                    std::vector<std::vector<double>> &&cnodes, std::vector<std::vector<double>> &&cweights,
                    std::string &&cdescription) :
        num_levels(cnum_nodes.size()), num_nodes(std::move(cnum_nodes)), precision(std::move(cprecision)),
        nodes(std::move(cnodes)), weights(std::move(cweights)), description(std::move(cdescription))
    {}
 
    template<bool useAscii> void write(std::ostream &os) const;
    template<bool useAscii> void read(std::istream &is);
 
    void read(const char* filename);
 
    int getNumLevels() const{ return num_levels; }
    int getNumPoints(int level) const{ checkLevel(level, "number of points"); return num_nodes[level]; }
    int getIExact(int level) const{ checkLevel(level, "i-exactness"); return num_nodes[level] -1; }
    int getQExact(int level) const{ checkLevel(level, "q-exactness"); return precision[level]; }
 
    void getWeightsNodes(int level, std::vector<double> &w, std::vector<double> &x) const;
    void getWeightsNodes(int level, double w[], double x[]) const;
 
    const char* getDescription() const;
 
protected:
    void checkLevel(int level, std::string const &op) const{
        if (level >= num_levels)
            throw std::runtime_error(std::string("ERROR: needed custom rule ") + op + " with level " + std::to_string(level) + " but the table ends at " + std::to_string(num_levels - 1));
    }
 
private:
    int num_levels;
    std::vector<int> num_nodes;
    std::vector<int> precision;
 
    std::vector<std::vector<double>> nodes;
    std::vector<std::vector<double>> weights;
    std::string description;
};
 
CustomTabulated getSubrules(CustomTabulated &ct, int start_index, int stride, std::string description);
 
namespace OneDimensionalMeta{
    int getNumPoints(int level, TypeOneDRule rule);
    int getIExact(int level, TypeOneDRule rule);
    int getQExact(int level, TypeOneDRule rule);
 
    bool isNonNested(TypeOneDRule rule);
    bool isSequence(TypeOneDRule rule);
    bool isGlobal(TypeOneDRule rule);
    bool isSingleNodeGrowth(TypeOneDRule rule);
    bool isLocalPolynomial(TypeOneDRule rule);
    bool isWavelet(TypeOneDRule rule);
    bool isFourier(TypeOneDRule rule);
 
    const char* getHumanString(TypeOneDRule rule);
 
 
    inline TypeDepth getControurType(TypeDepth type){
        if ((type == type_level) || (type == type_iptotal) || (type == type_qptotal)){
            return type_level;
        }else if ((type == type_curved) || (type == type_ipcurved) || (type == type_qpcurved)){
            return type_curved;
        }else{
            return type_hyperbolic;
        }
    }
 
    inline bool isExactLevel(TypeDepth type){
        return (type == type_level) || (type == type_curved) || (type == type_hyperbolic) || (type == type_tensor);
    }
 
    inline bool isExactInterpolation(TypeDepth type){
        return (type == type_iptotal) || (type == type_ipcurved) || (type == type_iphyperbolic) || (type == type_iptensor);
    }
 
    inline bool isExactQuadrature(TypeDepth type){
        return (type == type_qptotal) || (type == type_qpcurved) || (type == type_qphyperbolic) || (type == type_qptensor);
    }
 
 
    inline TypeDepth getSelectionType(TypeDepth type){
        if ((type == type_level) || (type == type_curved) || (type == type_hyperbolic)){
            return type_level;
        }else if ((type == type_iptotal) || (type == type_ipcurved) || (type == type_iphyperbolic)){
            return type_iptotal;
        }else{
            return type_qptotal;
        }
    }
 
    inline bool isTypeCurved(TypeDepth type){ return (getControurType(type) == type_curved); }
}
 
namespace OneDimensionalNodes{
    // non-nested rules
    void getGaussLegendre(int m, std::vector<double> &w, std::vector<double> &x);
    void getChebyshev(int m, std::vector<double> &w, std::vector<double> &x);
    void getGaussChebyshev1(int m, std::vector<double> &w, std::vector<double> &x);
    void getGaussChebyshev2(int m, std::vector<double> &w, std::vector<double> &x);
    void getGaussJacobi(int m, std::vector<double> &w, std::vector<double> &x, double alpha, double beta);
    void getGaussHermite(int m, std::vector<double> &w, std::vector<double> &x, double alpha);
    void getGaussLaguerre(int m, std::vector<double> &w, std::vector<double> &x, double alpha);
 
    // nested rules
    std::vector<double> getClenshawCurtisNodes(int level);
    double getClenshawCurtisWeight(int level, int point);
 
    std::vector<double> getClenshawCurtisNodesZero(int level); // assuming zero boundary
    double getClenshawCurtisWeightZero(int level, int point); // assuming zero boundary
 
    std::vector<double> getFejer2Nodes(int level);
    double getFejer2Weight(int level, int point);
 
    std::vector<double> getRLeja(int n);
    std::vector<double> getRLejaCentered(int n);
    std::vector<double> getRLejaShifted(int n);
 
    std::vector<double> getFourierNodes(int level);
}
 
}
 
#endif