threed_beam_fea.h

Go to the documentation of this file.

// Copyright 2015. All rights reserved.
//
// Redistribution and use in source and binary forms, with or without
// modification, are permitted provided that the following conditions are met:
//
// * Redistributions of source code must retain the above copyright notice,
//   this list of conditions and the following disclaimer.
// * 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.
//
// 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.
//
// Author: ryan.latture@gmail.com (Ryan Latture)

#ifndef THREED_BEAM_FEA_H
#define THREED_BEAM_FEA_H

#ifdef EIGEN_USE_MKL_ALL
#include <Eigen/PardisoSupport>
#else
#include <Eigen/SparseLU>
#endif

#include <Eigen/Core>
#include <Eigen/Geometry>
#include <Eigen/SparseCore>

#include "containers.h"
#include "options.h"
#include "summary.h"
#include "csv_parser.h"

namespace fea {

    typedef Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic> GlobalStiffMatrix;

    typedef Eigen::Matrix<double, 12, 12, Eigen::RowMajor> LocalMatrix;

    typedef Eigen::Matrix<double, Eigen::Dynamic, 1> ForceVector;

    typedef Eigen::SparseMatrix<double> SparseMat;

    inline double norm(const Node &n1, const Node &n2);

    class GlobalStiffAssembler {

    public:

        GlobalStiffAssembler() {
            Kelem.setZero();
            Klocal.setZero();
            Aelem.setZero();
            AelemT.setZero();
            SparseKelem.resize(12, 12);
            SparseKelem.reserve(40);
        };

        void operator()(SparseMat &Kg, const Job &job, const std::vector<Tie> &ties);

        void calcKelem(unsigned int i, const Job &job);

        void calcAelem(const Eigen::Vector3d &nx, const Eigen::Vector3d &nz);

        LocalMatrix getKelem() {
            return Kelem;
        }

        LocalMatrix getAelem() {
            return Aelem;
        }

    private:
        LocalMatrix Kelem;
        LocalMatrix Klocal;
        LocalMatrix Aelem;
        LocalMatrix AelemT;
        SparseMat SparseKelem;
    };

    void loadBCs(SparseMat &Kg, ForceVector &force_vec, const std::vector<BC> &BCs, unsigned int num_nodes);

    void loadEquations(SparseMat &Kg, const std::vector<Equation> &equations, unsigned int num_nodes, unsigned int num_bcs);

    void loadTies(std::vector<Eigen::Triplet<double> > &triplets, const std::vector<Tie> &ties);


    std::vector<std::vector<double> > computeTieForces(const std::vector<Tie> &ties,
                                                       const std::vector<std::vector<double> > &nodal_displacements);

    void loadForces(SparseMat &force_vec, const std::vector<Force> &forces);

    Summary solve(const Job &job,
                  const std::vector<BC> &BCs,
                  const std::vector<Force> &forces,
                  const std::vector<Tie> &ties,
                  const std::vector<Equation> &equations,
                  const Options &options);
} // namespace fea

#endif // THREED_BEAM_FEA_H