b2super_element.H

//------------------------------------------------------------------------
// b2super_element.H --
//
// written by doreille
//            Thomas Blome <thomas.blome@dlr.de>
//
// Copyright (c) 2010 SMR Engineering & Development SA
//                         2502 Bienne, Switzerland
//
// (c) 2023 Deutsches Zentrum für Luft- und Raumfahrt (DLR) e.V.
//          Linder Höhe, 51147 Köln
//
//
// All Rights Reserved.  Proprietary source code.  The contents of
// this file may not be disclosed to third parties, copied or
// duplicated in any form, in whole or in part, without the prior
// written permission of SMR.
//------------------------------------------------------------------------
 
#ifndef B2SUPER_ELEMENT_H_
#define B2SUPER_ELEMENT_H_
 
#include <vector>
 
#include "elements/properties/b2super_element_properties.H"
#include "model/b2element.H"
#include "model/b2node.H"
#include "utils/b2linear_algebra.H"
 
namespace b2000 {
 
template <typename T, typename MATRIX_TYPE>
class SuperElement : public TypedElement<T> {
public:
    using property_type = SuperElementProperty;
 
    using type_t = ObjectTypeComplete<SuperElement, typename TypedElement<T>::type_t>;
 
    int get_number_of_dof() const { return property->get_number_of_internal_dof(); }
 
    size_t set_global_dof_numbering(size_t index) {
        internal_dof_index = index;
        return index + property->get_number_of_internal_dof();
    }
 
    std::pair<size_t, size_t> get_global_dof_numbering() const {
        return std::pair<size_t, size_t>(
              internal_dof_index, internal_dof_index + property->get_number_of_internal_dof());
    }
 
    std::pair<size_t, size_t> get_interface_reduction_internal_dof() {
        std::pair<size_t, size_t> r = property->get_interface_reduction_internal_dof();
        r.first += internal_dof_index;
        r.second += internal_dof_index;
        return r;
    }
 
    void set_nodes(std::pair<int, Node* const*> nodes_) {
        nodes.assign(nodes_.second, nodes_.second + nodes_.first);
    }
 
    std::pair<int, Node* const*> get_nodes() const {
        return std::pair<int, Node* const*>(nodes.size(), &nodes[0]);
    }
 
    void set_property(ElementProperty* property_) {
        property = dynamic_cast<property_type*>(property_);
        if (property == nullptr) {
            TypeError() << "Super element has wrong property type " << typeid(*property_) << THROW;
        }
        // maybe print  property->get_object_name()  << std::endl;
    }
 
    const ElementProperty* get_property() const { return property; }
 
    const std::vector<Element::VariableInfo> get_value_info() const {
        assert(property != nullptr);
        return std::vector<Element::VariableInfo>(property->get_order(), Element::linear);
    }
 
    void get_value(
          Model& model, const bool linear, const EquilibriumSolution equilibrium_solution,
          const double time, const double delta_time,
          const b2linalg::Matrix<T, b2linalg::Mrectangle_constref>& dof,
          GradientContainer* gradient_container, SolverHints* solver_hints,
          b2linalg::Index& dof_numbering, b2linalg::Vector<T, b2linalg::Vdense>& value,
          const std::vector<bool>& d_value_d_dof_flags,
          std::vector<b2linalg::Matrix<T, MATRIX_TYPE> >& d_value_d_dof,
          b2linalg::Vector<T, b2linalg::Vdense>& d_value_d_time) {
        if (!dof_numbering.is_null()) {
            const size_t id = property->get_number_of_internal_dof();
            const size_t ed = property->get_number_of_external_dof();
            dof_numbering.resize(id + ed);
            size_t* ptr = &dof_numbering[0];
            for (size_t k = 0; k != nodes.size(); ++k) {
                std::pair<size_t, size_t> dn = nodes[k]->get_global_dof_numbering();
                while (dn.first != dn.second) {
                    if (ptr == &dof_numbering[0] + ed) {
                        Exception() << "Number of degree of freedom mismatch in super element "
                                    << Object::get_object_name() << THROW;
                    }
                    *ptr++ = dn.first++;
                }
            }
            if (ptr != &dof_numbering[0] + ed) {
                Exception() << "Number of degree of freedom mismatch in super element "
                            << Object::get_object_name() << THROW;
            }
            size_t is = internal_dof_index;
            size_t ie = is + id;
            for (; is != ie; ++is, ++ptr) { *ptr = is; }
        }
        if (dof.is_null() || value.is_null()) {
            property->get_value(value, d_value_d_dof_flags, d_value_d_dof);
        } else {
            // computation of the linear contribution of the dof to the value
            std::vector<b2linalg::Matrix<T, MATRIX_TYPE> > d_value_d_dof_tmp;
            std::vector<bool> d_value_d_dof_flags(dof.size2(), true);
            property->get_value(
                  value, d_value_d_dof_flags,
                  d_value_d_dof.empty() ? d_value_d_dof_tmp : d_value_d_dof);
            for (size_t k = 0; k != d_value_d_dof.size(); ++k) {
                value += (d_value_d_dof.empty() ? d_value_d_dof_tmp : d_value_d_dof)[k]
                         * dof[k](dof_numbering);
            }
        }
    }
 
    int face_field_order(const int face_id, const std::string& field_name) {
        return property->face_field_order(face_id, field_name);
    }
 
    bool face_field_linear_on_dof(const int face_id, const std::string& field_name) {
        return property->face_field_linear_on_dof(face_id, field_name);
    }
 
    int face_field_polynomial_sub_face(
          const int face_id, const std::string& field_name,
          b2linalg::Matrix<double, b2linalg::Mrectangle>& sub_nodes,
          std::vector<Element::Triangle>& sub_faces) {
        return property->face_field_polynomial_sub_face(face_id, field_name, sub_nodes, sub_faces);
    }
 
    void face_geom(
          const int face_id,
          const b2linalg::Vector<double, b2linalg::Vdense_constref>& internal_coor,
          b2linalg::Vector<double>& geom,
          b2linalg::Matrix<double, b2linalg::Mrectangle>& d_geom_d_icoor) {
        property->face_geom(face_id, internal_coor, geom, d_geom_d_icoor);
    }
 
    void face_field_value(
          const int face_id, const std::string& field_name,
          const b2linalg::Vector<double, b2linalg::Vdense_constref>& internal_coor,
          const b2linalg::Matrix<double, b2linalg::Mrectangle_constref>& dof, const double time,
          b2linalg::Vector<double, b2linalg::Vdense>& value,
          b2linalg::Matrix<double, b2linalg::Mrectangle>& d_value_d_icoor,
          b2linalg::Index& dof_numbering,
          b2linalg::Matrix<double, b2linalg::Mrectangle>& d_value_d_dof,
          b2linalg::Index& d_value_d_dof_dep_col) {
        property->face_field_value(
              face_id, field_name, internal_coor, dof, time, value, d_value_d_icoor, dof_numbering,
              d_value_d_dof, d_value_d_dof_dep_col);
        if (!dof_numbering.is_null()) {
            for (size_t i = 0; i != dof_numbering.size(); ++i) {
                dof_numbering[i] += internal_dof_index;
            }
        }
    }
 
    static type_t type;
 
private:
    size_t internal_dof_index;
    std::vector<Node*> nodes;
    property_type* property;
};
 
}  // namespace b2000
 
#endif /* B2SUPER_ELEMENT_H_ */