PolynomialPatch.cpp

Go to the documentation of this file.

/*
 *  Copyright (c) 2015, Chris Rogers
 *  All rights reserved.
 *  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 STFC 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 OWNER 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.
 */
 
#include "Utilities/GeneralOpalException.h"
 
#include "Fields/Interpolation/PolynomialPatch.h"
#include "Fields/Interpolation/SolveFactory.h"
 
namespace interpolation {
    PolynomialPatch::PolynomialPatch(
            Mesh* grid_points, Mesh* validity_region,
            std::vector<SquarePolynomialVector*> polynomials) {
        grid_points_     = grid_points;
        validity_region_ = validity_region;
        points_          = polynomials;
        // validate input data
        if (grid_points_ == nullptr)
            throw GeneralOpalException(
                    "PolynomialPatch::PolynomialPatch", "PolynomialPatch grid_points_ was nullptr");
        if (validity_region_ == nullptr)
            throw GeneralOpalException(
                    "PolynomialPatch::PolynomialPatch",
                    "PolynomialPatch validity_region_ was nullptr");
        if (points_.size() == 0) {
            throw GeneralOpalException(
                    "PolynomialPatch::PolynomialPatch",
                    "Could not make PolynomialPatch with 0 length polynomials vector");
        }
        Mesh::Iterator end = grid_points_->end();
        if (int(points_.size()) != end.toInteger()) {
            throw GeneralOpalException(
                    "PolynomialPatch::PolynomialPatch",
                    "Could not make PolynomialPatch with bad length polynomials vector");
        }
        point_dimension_ = grid_points_->getPositionDimension();
        if (validity_region_->getPositionDimension() != int(point_dimension_)) {
            throw GeneralOpalException(
                    "PolynomialPatch::PolynomialPatch",
                    "PolynomialPatch validity_region_ has bad dimension");
        }
        value_dimension_ = points_[0]->ValueDimension();
        for (size_t i = 0; i < points_.size(); ++i) {
            if (points_[i] == nullptr)
                throw GeneralOpalException(
                        "PolynomialPatch::PolynomialPatch",
                        "PolynomialPatch points_ element was nullptr");
            if (points_[i]->PointDimension() != point_dimension_) {
                throw GeneralOpalException(
                        "PolynomialPatch::PolynomialPatch",
                        "Polynomial with mismatched PointDimension in PolynomialPatch");
            }
            if (points_[i]->ValueDimension() != value_dimension_)
                throw GeneralOpalException(
                        "PolynomialPatch::PolynomialPatch",
                        "Polynomial with mismatched ValueDimension in PolynomialPatch");
        }
    }
 
    VectorMap* PolynomialPatch::clone() const {
        Mesh* new_mesh = nullptr;
        if (grid_points_ != nullptr) {
            new_mesh = grid_points_->clone();
        }
        Mesh* new_validity = nullptr;
        if (validity_region_ != nullptr) {
            new_validity = validity_region_->clone();
        }
        std::vector<SquarePolynomialVector*> new_points(points_.size());
        for (size_t i = 0; i < points_.size(); ++i) {
            if (points_[i] == nullptr) {
                new_points[i] = nullptr;
            } else {
                new_points[i] = new SquarePolynomialVector(*points_[i]);
            }
        }
        return new PolynomialPatch(new_mesh, new_validity, new_points);
    }
 
    PolynomialPatch::PolynomialPatch()
        : validity_region_(nullptr),
          grid_points_(nullptr),
          points_(),
          point_dimension_(0),
          value_dimension_(0) {}
 
    PolynomialPatch::~PolynomialPatch() {
        delete grid_points_;
        delete validity_region_;
        for (size_t i = 0; i < points_.size(); ++i)
            delete points_[i];
    }
 
    void PolynomialPatch::function(const double* point, double* value) const {
        Mesh::Iterator nearest = grid_points_->getNearest(point);
        std::vector<double> point_temp(point_dimension_);
        std::vector<double> nearest_pos = nearest.getPosition();
        for (size_t i = 0; i < point_dimension_; ++i)
            point_temp[i] = point[i] - nearest_pos[i];
        int points_index = nearest.toInteger();
        points_[points_index]->F(&point_temp[0], value);
    }
 
    SquarePolynomialVector* PolynomialPatch::getPolynomialVector(const double* point) const {
        Mesh::Iterator nearest = grid_points_->getNearest(point);
        int points_index       = nearest.toInteger();
        return points_[points_index];
    }
}  // namespace interpolation