/** * \file * \copyright * Copyright (c) 2012-2020, OpenGeoSys Community (http://www.opengeosys.org) * Distributed under a Modified BSD License. * See accompanying file LICENSE.txt or * http://www.opengeosys.org/project/license * */ #pragma once #ifdef OGS_USE_PYTHON #include <pybind11/pybind11.h> #endif // OGS_USE_PYTHON #include <algorithm> #include <iostream> #include <logog/include/logog.hpp> #include <vector> #include "BHEInflowPythonBoundaryConditionPythonSideInterface.h" #include "BaseLib/Error.h" #include "NumLib/IndexValueVector.h" #include "ProcessLib/BoundaryCondition/BoundaryCondition.h" #include "ProcessLib/BoundaryCondition/GenericNaturalBoundaryConditionLocalAssembler.h" #include "ProcessLib/HeatTransportBHE/BHE/BHETypes.h" #include "ProcessLib/Utils/CreateLocalAssemblers.h" #include "ProcessLib/Utils/ProcessUtils.h" #include "PythonBoundaryConditionLocalAssembler.h" namespace ProcessLib { //! A boundary condition whose values are computed by a Python script. template <typename BHEType> class BHEInflowPythonBoundaryCondition final : public BoundaryCondition { public: BHEInflowPythonBoundaryCondition( std::pair<GlobalIndexType, GlobalIndexType>&& in_out_global_indices, BHEType& bhe, BHEInflowPythonBoundaryConditionPythonSideInterface& py_bc_object) : _in_out_global_indices(std::move(in_out_global_indices)), _bhe(bhe), _py_bc_object(py_bc_object) { const auto g_idx_T_out = in_out_global_indices.second; // store the bc node ids to BHE network dataframe std::get<3>(_py_bc_object.dataframe_network).emplace_back(g_idx_T_out); } void getEssentialBCValues( const double t, const GlobalVector& /* x */, NumLib::IndexValueVector<GlobalIndexType>& bc_values) const override { bc_values.ids.resize(1); bc_values.values.resize(1); auto const& data_exchange = _py_bc_object.dataframe_network; // get the number of all boundary nodes const std::size_t n_bc_nodes = std::get<3>(data_exchange).size(); // get T_in bc_id bc_values.ids[0] = _in_out_global_indices.first; // get T_out bc_id auto const boundary_node_id = _in_out_global_indices.second; // return T_in from currently BHE dataframe column 2, // update flowrate and HeatTransferCoefficients for each BHE for (std::size_t i = 0; i < n_bc_nodes; i++) { // auto pair_flag_value = // _bc_data.bc_object->getDirichletBCValue(boundary_node_id); auto const dataframe_node_id = std::get<3>(data_exchange); auto const dataframe_Tin_val = std::get<1>(data_exchange); auto const dataframe_BHE_flowrate = std::get<4>(data_exchange); if (dataframe_node_id[i] == boundary_node_id) { bc_values.values[0] = dataframe_Tin_val[i]; _bhe.updateHeatTransferCoefficients(dataframe_BHE_flowrate[i]); break; } } // store the current time to network dataframe std::get<0>(_py_bc_object.dataframe_network) = t; } private: std::pair<GlobalIndexType, GlobalIndexType> const _in_out_global_indices; BHEType& _bhe; BHEInflowPythonBoundaryConditionPythonSideInterface& _py_bc_object; }; template <typename BHEType> std::unique_ptr<BHEInflowPythonBoundaryCondition<BHEType>> createBHEInflowPythonBoundaryCondition( std::pair<GlobalIndexType, GlobalIndexType>&& in_out_global_indices, BHEType& bhe, BHEInflowPythonBoundaryConditionPythonSideInterface& py_bc_object) { DBUG("Constructing BHEInflowPythonBoundaryCondition."); // In case of partitioned mesh the boundary could be empty, i.e. there is no // boundary condition. #ifdef USE_PETSC // For this special boundary condition the boundary condition is not empty // if the global indices are non-negative. if (in_out_global_indices.first < 0 && in_out_global_indices.second < 0) { return nullptr; } // If only one of the global indices (in or out) is negative the // implementation is not valid. if (in_out_global_indices.first < 0 || in_out_global_indices.second < 0) { OGS_FATAL( "The partition cuts the BHE into two independent parts. This " "behaviour is not implemented."); } #endif // USE_PETSC return std::make_unique<BHEInflowPythonBoundaryCondition<BHEType>>( std::move(in_out_global_indices), bhe, py_bc_object); } } // namespace ProcessLib