From b1b6f011ee8060dc3178decad02b67e4da766bdc Mon Sep 17 00:00:00 2001 From: David Schneider Date: Tue, 6 Feb 2024 09:23:47 +0100 Subject: [PATCH] Technical API code changes --- coupled_laplace_problem/CMakeLists.txt | 2 +- .../coupled_laplace_problem.cc | 65 +++++++------------ .../fancy_boundary_condition.cc | 28 +++----- 3 files changed, 35 insertions(+), 60 deletions(-) diff --git a/coupled_laplace_problem/CMakeLists.txt b/coupled_laplace_problem/CMakeLists.txt index 92a3e15..7995817 100644 --- a/coupled_laplace_problem/CMakeLists.txt +++ b/coupled_laplace_problem/CMakeLists.txt @@ -38,7 +38,7 @@ DEAL_II_INITIALIZE_CACHED_VARIABLES() PROJECT(${TARGET}) DEAL_II_INVOKE_AUTOPILOT() -FIND_PACKAGE(precice REQUIRED +FIND_PACKAGE(precice 3.0 REQUIRED HINTS ${PRECICE_DIR} ) MESSAGE(STATUS "Using the preCICE version found at ${precice_CONFIG}") diff --git a/coupled_laplace_problem/coupled_laplace_problem.cc b/coupled_laplace_problem/coupled_laplace_problem.cc index 1330bec..b760c41 100644 --- a/coupled_laplace_problem/coupled_laplace_problem.cc +++ b/coupled_laplace_problem/coupled_laplace_problem.cc @@ -28,7 +28,7 @@ #include // In addition to the deal.II header files, we include the preCICE API in order // to obtain access to preCICE specific functionality -#include +#include #include #include @@ -62,20 +62,20 @@ template class Adapter { public: - Adapter(const ParameterClass & parameters, + Adapter(const ParameterClass ¶meters, const types::boundary_id dealii_boundary_interface_id); double - initialize(const DoFHandler & dof_handler, + initialize(const DoFHandler &dof_handler, std::map &boundary_data, - const MappingQGeneric & mapping); + const MappingQGeneric &mapping); double advance(std::map &boundary_data, const double computed_timestep_length); // public precCICE solver interface - precice::SolverInterface precice; + precice::Participant precice; // Boundary ID of the deal.II triangulation, associated with the coupling // interface. The variable is defined in the constructor of this class and @@ -94,8 +94,6 @@ private: // These IDs are filled by preCICE during the initialization. We set a default // value of -1 in order to detect potential errors more easily. - int mesh_id; - int read_data_id; int n_interface_nodes; // DoF IndexSet, containing relevant coupling DoF indices at the coupling @@ -132,7 +130,7 @@ private: // which is associated with the coupling interface. template Adapter::Adapter( - const ParameterClass & parameters, + const ParameterClass ¶meters, const types::boundary_id deal_boundary_interface_id) : precice(parameters.participant_name, parameters.config_file, @@ -156,18 +154,12 @@ Adapter::Adapter( template double Adapter::initialize( - const DoFHandler & dof_handler, + const DoFHandler &dof_handler, std::map &boundary_data, - const MappingQGeneric & mapping) + const MappingQGeneric &mapping) { Assert(dim > 1, ExcNotImplemented()); - AssertDimension(dim, precice.getDimensions()); - - // In a first step, we get preCICE specific IDs from preCICE and store them in - // the respective variables. Later, they are used for data transfer. - mesh_id = precice.getMeshID(mesh_name); - read_data_id = precice.getDataID(read_data_name, mesh_id); - + AssertDimension(dim, precice.getMeshDimensions(mesh_name)); // Afterwards, we extract the number of interface nodes and the coupling DoFs // at the coupling interface from our deal.II solver via @@ -219,30 +211,24 @@ Adapter::initialize( // Now we have all information to define the coupling mesh and pass the // information to preCICE. - precice.setMeshVertices(mesh_id, - n_interface_nodes, - interface_nodes_positions.data(), - interface_nodes_ids.data()); + precice.setMeshVertices(mesh_name, + interface_nodes_positions, + interface_nodes_ids); // Then, we initialize preCICE internally calling the API function // `initialize()` - const double max_delta_t = precice.initialize(); - + precice.initialize(); + const double max_delta_t = precice.getMaxTimeStepSize(); // read first coupling data from preCICE if available (i.e. deal.II is // the second participant in a serial coupling scheme) - if (precice.isReadDataAvailable()) - { - precice.readBlockScalarData(read_data_id, - n_interface_nodes, - interface_nodes_ids.data(), - read_data.data()); - - // After receiving the coupling data in `read_data`, we convert it to - // the std::map `boundary_data` which is later needed in order to apply - // Dirichlet boundary conditions - format_precice_to_dealii(boundary_data); - } + precice.readData( + mesh_name, read_data_name, interface_nodes_ids, max_delta_t, read_data); + + // After receiving the coupling data in `read_data`, we convert it to + // the std::map `boundary_data` which is later needed in order to apply + // Dirichlet boundary conditions + format_precice_to_dealii(boundary_data); return max_delta_t; } @@ -261,16 +247,15 @@ Adapter::advance( // return, preCICE tells us the maximum admissible time-step size our // participant is allowed to compute in order to not exceed the next coupling // time step. - const double max_delta_t = precice.advance(computed_timestep_length); + precice.advance(computed_timestep_length); + const double max_delta_t = precice.getMaxTimeStepSize(); // As a next step, we obtain data, i.e. the boundary condition, from another // participant. We have already all IDs and just need to convert our obtained // data to the deal.II compatible 'boundary map' , which is done in the // format_deal_to_precice function. - precice.readBlockScalarData(read_data_id, - n_interface_nodes, - interface_nodes_ids.data(), - read_data.data()); + precice.readData( + mesh_name, read_data_name, interface_nodes_ids, max_delta_t, read_data); format_precice_to_dealii(boundary_data); diff --git a/coupled_laplace_problem/fancy_boundary_condition.cc b/coupled_laplace_problem/fancy_boundary_condition.cc index 46fc24a..60fab56 100644 --- a/coupled_laplace_problem/fancy_boundary_condition.cc +++ b/coupled_laplace_problem/fancy_boundary_condition.cc @@ -1,9 +1,10 @@ // This program does not use any deal.II functionality and depends only on // preCICE and the standard libraries. -#include +#include #include #include +#include // The program computes a time-varying parabolic boundary condition, which is // passed to preCICE and serves as Dirichlet boundary condition for the other @@ -47,17 +48,11 @@ main() const int commRank = 0; const int commSize = 1; - precice::SolverInterface precice(solverName, - configFileName, - commRank, - commSize); + precice::Participant precice(solverName, configFileName, commRank, commSize); - const int meshID = precice.getMeshID(meshName); - const int dimensions = precice.getDimensions(); + const int dimensions = precice.getMeshDimensions(meshName); const int numberOfVertices = 6; - const int dataID = precice.getDataID(dataWriteName, meshID); - // Set up data structures std::vector writeData(numberOfVertices); std::vector vertexIDs(numberOfVertices); @@ -81,13 +76,11 @@ main() } // Pass the vertices to preCICE - precice.setMeshVertices(meshID, - numberOfVertices, - vertices.data(), - vertexIDs.data()); + precice.setMeshVertices(meshName, vertices, vertexIDs); // initialize the Solverinterface - double dt = precice.initialize(); + precice.initialize(); + double dt = precice.getMaxTimeStepSize(); // Start time loop const double end_time = 1; @@ -99,13 +92,10 @@ main() { std::cout << "Boundary participant: writing coupling data \n"; - precice.writeBlockScalarData(dataID, - numberOfVertices, - vertexIDs.data(), - writeData.data()); + precice.writeData(meshName, dataWriteName, vertexIDs, writeData); } - dt = precice.advance(dt); + precice.advance(dt); std::cout << "Boundary participant: advancing in time\n"; time += dt; -- 2.39.5