EshelbianRestart.cpp

Go to the documentation of this file.

/**
 * @file EshelbianRestart.cpp
 * @brief Checked sidecar metadata for native restart vectors.
 */
 
#include <MoFEM.hpp>
using namespace MoFEM;
 
#include <EshelbianPlasticity.hpp>
#include <EshelbianRestart.hpp>
 
#include <cerrno>
#include <cstdio>
#include <cstring>
#include <filesystem>
#include <fstream>
#include <limits>
#include <set>
#include <sstream>
 
namespace EshelbianPlasticity {
namespace {
 
using JsonObject = boost::json::object;
 
// Complete rank-local I/O and validation before any rank returns an error.
MoFEMErrorCode checkCollectiveLayoutError(MPI_Comm comm,
                                          const std::string &local_error) {
  MoFEMFunctionBegin;
  int rank, mpi_size;
  CHKERR MPI_Comm_rank(comm, &rank);
  CHKERR MPI_Comm_size(comm, &mpi_size);
  const int local_rank = local_error.empty() ? mpi_size : rank;
  int error_rank = mpi_size;
  CHKERR MPI_Allreduce(&local_rank, &error_rank, 1, MPI_INT, MPI_MIN, comm);
  if (error_rank != mpi_size) {
    std::string message;
    if (rank == error_rank)
      message = local_error.substr(0, 16384);
    int length = static_cast<int>(message.size());
    CHKERR MPI_Bcast(&length, 1, MPI_INT, error_rank, comm);
    message.resize(length);
    CHKERR MPI_Bcast(message.data(), length, MPI_CHAR, error_rank, comm);
    SETERRQ(comm, MOFEM_DATA_INCONSISTENCY, "%s", message.c_str());
  }
  MoFEMFunctionReturn(0);
}
 
MoFEMErrorCode getRestartLayout(const EshelbianCore &ep, Vec state,
                                JsonObject &layout) {
  MoFEMFunctionBegin;
  const auto comm = ep.mField.get_comm();
  CHKERR checkCollectiveLayoutError(
      comm, state && ep.dmElastic
                ? ""
                : "Restart layout requires an elastic DM and vector");
  int comparison = MPI_UNEQUAL;
  CHKERR MPI_Comm_compare(
      comm, PetscObjectComm(reinterpret_cast<PetscObject>(state)), &comparison);
  CHKERR checkCollectiveLayoutError(
      comm, comparison == MPI_IDENT || comparison == MPI_CONGRUENT
                ? ""
                : "Restart vector and elastic DM use different communicators");
 
  int mpi_size;
  CHKERR MPI_Comm_size(comm, &mpi_size);
  PetscInt global_size, local_size;
  CHKERR VecGetSize(state, &global_size);
  CHKERR VecGetLocalSize(state, &local_size);
  std::vector<PetscInt> partition_sizes(mpi_size);
  CHKERR MPI_Allgather(&local_size, 1, MPIU_INT, partition_sizes.data(), 1,
                       MPIU_INT, comm);
  boost::json::array partitions;
  for (const auto size : partition_sizes)
    partitions.emplace_back(static_cast<std::int64_t>(size));
 
  const Problem *problem_ptr = nullptr;
  CHKERR DMMoFEMGetProblemPtr(ep.dmElastic.get(), &problem_ptr);
  CHKERR checkCollectiveLayoutError(
      comm, problem_ptr && problem_ptr->getNumeredRowDofsPtr()
                ? ""
                : "Restart elastic DM has no numbered row DOFs");
  CHKERR checkCollectiveLayoutError(
      comm, problem_ptr->getNbDofsRow() == global_size
                ? ""
                : "Restart vector size differs from the elastic DM layout");
 
  const auto &features = ep.physicalEquations->getFeatures();
  layout["version"] = 1;
  layout["formulation"] =
      features.test(PhysicalEquations::AUXILIARY_LOGARITHMIC_STRESS)
          ? "auxiliary_logarithmic_stress"
          : ((features & PhysicalEquations::noStretchMask).any()
                 ? "no_stretch"
                 : "direct_stretch");
  // Native material IDs must remain independent of option-list positions.
  const auto material_model = ep.physicalEquations->getMaterialModel();
  layout["material_model"] = material_model == PhysicalEquations::Meta
                                 ? 5
                                 : static_cast<int>(material_model);
  // Keep legacy restart values independent of the feature bit positions.
  layout["stretch_handling"] =
      features.test(PhysicalEquations::NO_STRETCH_NONLINEAR)
          ? 2
          : (features.test(PhysicalEquations::NO_STRETCH_LINEAR) ? 1 : 0);
  layout["stretch_selector"] = static_cast<int>(ep.stretchSelector);
  layout["gradient_approximator"] = static_cast<int>(ep.gradApproximator);
  layout["rotation_selector"] = static_cast<int>(ep.rotSelector);
  layout["space_order"] = ep.spaceOrder;
  layout["broken_hdiv_basis"] = ApproximationBaseNames[ep.brokenHdivBase];
  layout["l2_user_base_scale"] = static_cast<bool>(ep.l2UserBaseScale);
  layout["vector_global_size"] = static_cast<std::int64_t>(global_size);
  layout["mpi_ranks"] = mpi_size;
  layout["vector_partition_sizes"] = std::move(partitions);
 
  std::set<std::string> field_names;
  for (const auto &dof : *problem_ptr->getNumeredRowDofsPtr())
    if (dof->getPetscGlobalDofIdx() >= 0)
      field_names.insert(dof->getName());
  std::string local_names;
  for (const auto &name : field_names)
    local_names += name + '\n';
  CHKERR checkCollectiveLayoutError(
      comm, local_names.size() <= std::numeric_limits<int>::max()
                ? ""
                : "Restart field names exceed the MPI metadata size limit");
  const int local_length = static_cast<int>(local_names.size());
  std::vector<int> lengths(mpi_size), offsets(mpi_size);
  CHKERR MPI_Allgather(&local_length, 1, MPI_INT, lengths.data(), 1, MPI_INT,
                       comm);
  std::size_t total_length = 0;
  for (int rank = 0; rank != mpi_size; ++rank) {
    offsets[rank] = static_cast<int>(total_length);
    total_length += lengths[rank];
  }
  CHKERR checkCollectiveLayoutError(
      comm, total_length <= std::numeric_limits<int>::max()
                ? ""
                : "Restart field layout exceeds the MPI metadata size limit");
  std::string all_names(total_length, '\0');
  CHKERR MPI_Allgatherv(local_names.data(), local_length, MPI_CHAR,
                        all_names.data(), lengths.data(), offsets.data(),
                        MPI_CHAR, comm);
  std::istringstream names_stream(all_names);
  for (std::string name; std::getline(names_stream, name);)
    if (!name.empty())
      field_names.insert(name);
 
  const auto material_fields =
      ep.physicalEquations->getMaterialFieldDefinitions(ep);
  boost::json::array material_names;
  for (const auto &field : material_fields) {
    if (field.order < 0)
      continue;
    CHKERR checkCollectiveLayoutError(
        comm, field_names.count(field.name)
                  ? ""
                  : "Restart material field has no active DOFs: " + field.name);
    material_names.emplace_back(field.name);
  }
  layout["material_fields"] = std::move(material_names);
 
  boost::json::array fields;
  for (const auto &field_name : field_names) {
    CHKERR checkCollectiveLayoutError(comm, ep.mField.check_field(field_name)
                                                ? ""
                                                : "Restart field is missing: " +
                                                      field_name);
    const auto field_ptr = ep.mField.get_field_structure(field_name);
    IS raw_is = nullptr;
    CHKERR DMMoFEMGetFieldIS(ep.dmElastic.get(), RowColData::ROW,
                             field_name.c_str(), &raw_is);
    SmartPetscObj<IS> field_is(raw_is);
    PetscInt field_size, local_field_size;
    CHKERR ISGetSize(field_is, &field_size);
    CHKERR ISGetLocalSize(field_is, &local_field_size);
    const auto &dofs =
        problem_ptr->getNumeredRowDofsPtr()->get<PetscGlobalIdx_mi_tag>();
    const PetscInt *indices = nullptr;
    CHKERR ISGetIndices(field_is, &indices);
    int minimum_order = std::numeric_limits<int>::max();
    int maximum_order = -1;
    std::string local_error;
    for (PetscInt n = 0; n != local_field_size; ++n) {
      const auto dof = dofs.find(indices[n]);
      if (dof == dofs.end()) {
        local_error = "Restart DM is missing a numbered DOF of " + field_name;
        break;
      }
      const int order = (*dof)->getDofOrder();
      minimum_order = std::min(minimum_order, order);
      maximum_order = std::max(maximum_order, order);
    }
    CHKERR ISRestoreIndices(field_is, &indices);
    CHKERR checkCollectiveLayoutError(comm, local_error);
    CHKERR MPI_Allreduce(MPI_IN_PLACE, &minimum_order, 1, MPI_INT, MPI_MIN,
                         comm);
    CHKERR MPI_Allreduce(MPI_IN_PLACE, &maximum_order, 1, MPI_INT, MPI_MAX,
                         comm);
    CHKERR checkCollectiveLayoutError(
        comm, field_size > 0 && minimum_order >= 0
                  ? ""
                  : "Restart field has no active DOFs: " + field_name);
 
    JsonObject field{
        {"name", field_name},
        {"components", static_cast<int>(field_ptr->getNbOfCoeffs())},
        {"space", field_ptr->getSpaceName()},
        {"basis", field_ptr->getApproxBaseName()},
        {"continuity", field_ptr->getContinuityName()},
        {"minimum_dof_order", minimum_order},
        {"maximum_dof_order", maximum_order},
        {"global_dofs", static_cast<std::int64_t>(field_size)}};
    const auto material_field = std::find_if(
        material_fields.begin(), material_fields.end(), [&](const auto &entry) {
          return entry.order >= 0 && entry.name == field_name;
        });
    if (material_field != material_fields.end()) {
      field["order_policy"] = "FieldOrders::stretch(space_order)";
      field["requested_order"] = material_field->order;
    }
    fields.emplace_back(std::move(field));
  }
  layout["fields"] = std::move(fields);
  MoFEMFunctionReturn(0);
}
 
std::string compareRestartLayout(const boost::json::value &stored,
                                 const JsonObject &expected,
                                 const std::string &file_name) {
  if (!stored.is_object())
    return "Restart metadata is not a JSON object: " + file_name;
  const auto &object = stored.as_object();
  for (const auto &entry : expected) {
    const auto stored_value = object.if_contains(entry.key());
    if (!stored_value || *stored_value != entry.value())
      return "Restart layout mismatch for '" +
             std::string(entry.key().data(), entry.key().size()) + "' in " +
             file_name +
             "; use the formulation, fields, orders and MPI partition of the "
             "saved native restart";
  }
  if (object.size() != expected.size())
    return "Unexpected entries in restart layout metadata: " + file_name;
  return "";
}
 
} // namespace
 
MoFEMErrorCode writeRestartLayout(const EshelbianCore &ep, Vec state,
                                  const std::string &binary_file) {
  MoFEMFunctionBegin;
  JsonObject layout;
  CHKERR getRestartLayout(ep, state, layout);
  const auto comm = ep.mField.get_comm();
  int rank;
  CHKERR MPI_Comm_rank(comm, &rank);
  const std::string file_name = binary_file + ".layout.json";
  std::string local_error;
  if (rank == 0) {
    const std::string temporary_file = file_name + ".tmp";
    try {
      std::ofstream output(temporary_file, std::ios::out | std::ios::trunc);
      output << boost::json::serialize(layout) << '\n';
      output.flush();
      const bool written = output.good();
      output.close();
      if (!written || output.fail())
        local_error = "Cannot write restart layout metadata: " + file_name;
      else if (std::rename(temporary_file.c_str(), file_name.c_str()))
        local_error = "Cannot install restart layout metadata " + file_name +
                      ": " + std::strerror(errno);
    } catch (const std::exception &error) {
      local_error = "Cannot write restart layout metadata " + file_name + ": " +
                    error.what();
    }
    if (!local_error.empty())
      std::remove(temporary_file.c_str());
  }
  CHKERR checkCollectiveLayoutError(comm, local_error);
  MoFEMFunctionReturn(0);
}
 
MoFEMErrorCode validateRestartLayout(const EshelbianCore &ep, Vec state,
                                     const std::string &binary_file) {
  MoFEMFunctionBegin;
  const auto comm = ep.mField.get_comm();
  int rank;
  CHKERR MPI_Comm_rank(comm, &rank);
  const std::string file_name = binary_file + ".layout.json";
  int exists = 0;
  std::string contents, local_error;
  if (rank == 0) {
    try {
      std::error_code error;
      exists = std::filesystem::exists(file_name, error);
      if (error) {
        local_error = "Cannot inspect restart layout metadata " + file_name +
                      ": " + error.message();
      } else if (exists) {
        std::ifstream input(file_name);
        std::ostringstream buffer;
        buffer << input.rdbuf();
        if (!input.is_open() || input.bad() || buffer.fail())
          local_error = "Cannot read restart layout metadata: " + file_name;
        else
          contents = buffer.str();
        if (contents.size() >
            static_cast<std::size_t>(std::numeric_limits<int>::max()))
          local_error = "Restart layout metadata is too large: " + file_name;
      }
    } catch (const std::exception &error) {
      local_error = "Cannot read restart layout metadata " + file_name + ": " +
                    error.what();
    }
  }
  CHKERR checkCollectiveLayoutError(comm, local_error);
  CHKERR MPI_Bcast(&exists, 1, MPI_INT, 0, comm);
  if (!exists) {
    CHKERR checkCollectiveLayoutError(
        comm, ep.physicalEquations->getFeatures().test(
                  PhysicalEquations::AUXILIARY_LOGARITHMIC_STRESS)
                  ? "Auxiliary logarithmic-stress restart requires " +
                        file_name +
                        "; an old raw stretch vector cannot be loaded into "
                        "the D/theta/Td fields"
                  : "");
    MoFEMFunctionReturnHot(0);
  }
 
  int length = static_cast<int>(contents.size());
  CHKERR MPI_Bcast(&length, 1, MPI_INT, 0, comm);
  contents.resize(length);
  CHKERR MPI_Bcast(contents.data(), length, MPI_CHAR, 0, comm);
  JsonObject expected;
  CHKERR getRestartLayout(ep, state, expected);
  try {
    boost::system::error_code parse_error;
    const auto stored = boost::json::parse(contents, parse_error);
    local_error = parse_error
                      ? "Invalid restart layout JSON in " + file_name + ": " +
                            parse_error.message()
                      : compareRestartLayout(stored, expected, file_name);
  } catch (const std::exception &error) {
    local_error = "Cannot validate restart layout metadata " + file_name +
                  ": " + error.what();
  }
  CHKERR checkCollectiveLayoutError(comm, local_error);
  MoFEMFunctionReturn(0);
}
 
} // namespace EshelbianPlasticity