BCTopSurfaceEnergyBalance.h

Go to the documentation of this file.

#ifndef BOUNDARY_CONDITION_SURFACE_ENERGY_BALANCE
#define BOUNDARY_CONDITION_SURFACE_ENERGY_BALANCE
 
#include <cmath>
#include <memory>
#include <string>
#include <valarray>
#include <vector>
 
#include "GenericBoundaryCondition.h"
#include "IniParser.h"
#include "PhysicalConstants.h"
#include "SimulationClassBase.h"
 
template <typename T>
class BCTopSurfaceEnergyBalance : public GenericBoundaryCondition<T>
{
    private:
        T m_area{};
        T m_stefan_boltzmann_const{};
        std::string m_temperature_type{"Temperature"};
 
    public:
        void writeBoundaryCondition() override;
 
        BCTopSurfaceEnergyBalance(SimulationClassBase<T>* sim, std::vector<int> boundary_points,
                                  std::vector<std::vector<int>> connected_points, int facet_nr,
                                  T area, std::string face_prefix);
        ~BCTopSurfaceEnergyBalance() override = default;
};
 
template <typename T>
BCTopSurfaceEnergyBalance<T>::BCTopSurfaceEnergyBalance(
    SimulationClassBase<T>* sim, std::vector<int> boundary_points,
    std::vector<std::vector<int>> connected_points, int facet_nr, T area, std::string face_prefix)
    : GenericBoundaryCondition<T>{sim, boundary_points, connected_points, facet_nr, face_prefix},
      m_area{area},
      m_stefan_boltzmann_const(PhysicalConstants::stefan_boltzmann_const)
{
    try
    {
        m_temperature_type =
            this->m_sim->m_simulation_config.getStringParameters("boundary_temperature");
    }
    catch (const BadInput&)
    {
        // This value is default set, therefor no error is propagated.
    }
}
 
template <typename T>
void BCTopSurfaceEnergyBalance<T>::writeBoundaryCondition()
{
    T q_sol;
    try
    {
        q_sol = this->m_sim->getFieldValue(
            "QSol", this->m_facet_number);  // It would be more efficient to not save them to a
                                            // variable, but I think code readability could suffer
                                            // too much, if I only used temporaries :(
    }
    catch (const BadInput& e)
    {
        std::cerr << "[BCTopSurfaceEnergyBalance]: SolarFlux module not "
                     "loaded or accessable!\n "
                     "[BCTopSurfaceEnergyBalance]:"
                  << e.what() << '\n';
        throw BadInput(static_cast<std::string>("[BCTopSurfaceEnergyBalance]: SolarFlux module not "
                                                "loaded or accessable!\n "
                                                "[BCTopSurfaceEnergyBalance]:")
                       + static_cast<std::string>(e.what()));
    }
    for (int i{0}; i < this->m_boundary_points.size(); i++)
    {
        int number_of_boundary_points{this->m_boundary_points[i]};
        T surface_temperature, emissivity;
        try
        {
            surface_temperature =
                this->m_sim->getFieldValue(m_temperature_type, this->m_facet_number);
            emissivity = this->m_sim->getFieldValue("Emissivity", this->m_facet_number);
        }
        catch (const BadInput& e)
        {
            std::cerr << "[BCTopSurfaceEnergyBalance]: Emissivity module not "
                         "loaded or accessable or wrong temperature_type name!\n "
                         "[BCTopSurfaceEnergyBalance]:"
                      << e.what() << '\n';
            throw BadInput(
                static_cast<std::string>("[BCTopSurfaceEnergyBalance]: Emissivity module not "
                                         "loaded or accessable or wrong temperature_type name!\n "
                                         "[BCTopSurfaceEnergyBalance]:")
                + static_cast<std::string>(e.what()));
        }
        for (auto& val : this->m_capacitance_matrix_elements[i])
        {
            val *= 0.0;
        }
        this->m_stiffness_matrix_elements[i][0] += 4 * m_stefan_boltzmann_const * emissivity
                                                   * std::pow(surface_temperature, 3)
                                                   * m_area;  // G_rad_fo * area
        this->m_forcing_vector_elements[i] +=
            (3 * m_stefan_boltzmann_const * emissivity * std::pow(surface_temperature, 4) + q_sol)
            * m_area;  // Q_sol * area - G_rad_const * area
    }
}
 
#endif