]> https://gitweb.dealii.org/ - code-gallery.git/commitdiff
Technical API code changes
authorDavid Schneider <david.schneider@ipvs.uni-stuttgart.de>
Tue, 6 Feb 2024 08:23:47 +0000 (09:23 +0100)
committerDavid Schneider <david.schneider@ipvs.uni-stuttgart.de>
Tue, 6 Feb 2024 08:23:47 +0000 (09:23 +0100)
coupled_laplace_problem/CMakeLists.txt
coupled_laplace_problem/coupled_laplace_problem.cc
coupled_laplace_problem/fancy_boundary_condition.cc

index 92a3e15ab20f36d229571c73a8f0e1b86ae80371..7995817d2f25fc27062783c82a86124921f25768 100644 (file)
@@ -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}")
index 1330bec9653a47e4f138b8e48500fc41332d624b..b760c41fbd2ccbf69b3ef01eed77727893c31e8a 100644 (file)
@@ -28,7 +28,7 @@
 #include <deal.II/numerics/vector_tools.h>
 // In addition to the deal.II header files, we include the preCICE API in order
 // to obtain access to preCICE specific functionality
-#include <precice/SolverInterface.hpp>
+#include <precice/precice.hpp>
 
 #include <fstream>
 #include <iostream>
@@ -62,20 +62,20 @@ template <int dim, typename ParameterClass>
 class Adapter
 {
 public:
-  Adapter(const ParameterClass &   parameters,
+  Adapter(const ParameterClass    &parameters,
           const types::boundary_id dealii_boundary_interface_id);
 
   double
-  initialize(const DoFHandler<dim> &                    dof_handler,
+  initialize(const DoFHandler<dim>                     &dof_handler,
              std::map<types::global_dof_index, double> &boundary_data,
-             const MappingQGeneric<dim> &               mapping);
+             const MappingQGeneric<dim>                &mapping);
 
   double
   advance(std::map<types::global_dof_index, double> &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 <int dim, typename ParameterClass>
 Adapter<dim, ParameterClass>::Adapter(
-  const ParameterClass &   parameters,
+  const ParameterClass    &parameters,
   const types::boundary_id deal_boundary_interface_id)
   : precice(parameters.participant_name,
             parameters.config_file,
@@ -156,18 +154,12 @@ Adapter<dim, ParameterClass>::Adapter(
 template <int dim, typename ParameterClass>
 double
 Adapter<dim, ParameterClass>::initialize(
-  const DoFHandler<dim> &                    dof_handler,
+  const DoFHandler<dim>                     &dof_handler,
   std::map<types::global_dof_index, double> &boundary_data,
-  const MappingQGeneric<dim> &               mapping)
+  const MappingQGeneric<dim>                &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<dim, ParameterClass>::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<dim, ParameterClass>::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);
 
index 46fc24a29406d3144aab8f6721307cabc1b8ef38..60fab5689d8a0d8808368881c43a4c19f6f3f990 100644 (file)
@@ -1,9 +1,10 @@
 // This program does not use any deal.II functionality and depends only on
 // preCICE and the standard libraries.
-#include <precice/SolverInterface.hpp>
+#include <precice/precice.hpp>
 
 #include <iostream>
 #include <sstream>
+#include <vector>
 
 // 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<double> writeData(numberOfVertices);
   std::vector<int>    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;

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.