ikarus-project / ikarus

Ikarus is a C++-based library that uses DUNE modules to solve partial differential equations with the finite element method and more.
https://ikarus-project.github.io/
Other
5 stars 3 forks source link

add a wrapper to evaluate the required results #143

Closed tarun-mitruka closed 1 year ago

tarun-mitruka commented 1 year ago
 ResultTypeMap<double> result;
  auto gridView        = fe.localView().globalBasis().gridView();
  auto scalarBasis     = makeConstSharedBasis(gridView, lagrangeDG<1>());
  auto localScalarView = scalarBasis->localView();
  std::vector<Dune::FieldVector<double, 3>> stressVector(scalarBasis->size());
  auto ele = elements(gridView).begin();

  localScalarView.bind(*ele);
  const auto& fe2              = localScalarView.tree().finiteElement();
  const auto& referenceElement = Dune::ReferenceElements<double, gridDim>::general(ele->type());
  for (auto c = 0UL; c < fe2.size(); ++c) {
    const auto fineKey                        = fe2.localCoefficients().localKey(c);
    const auto nodalPositionInChildCoordinate = referenceElement.position(fineKey.subEntity(), fineKey.codim());
    auto coord                                = toEigen(nodalPositionInChildCoordinate);
    fe.calculateAt(resultRequirements, coord, result);
    Eigen::Vector3d computedResult = result.getResult(ResultType::cauchyStress).eval();
    const auto nodeIndex           = localScalarView.index(localScalarView.tree().localIndex(c))[0];
    stressVector[nodeIndex]        = toDune(computedResult);
    for (auto voigtIndex = 0UL; voigtIndex < 3; ++voigtIndex) {
      const auto FEStressComponent       = stressVector[nodeIndex][voigtIndex];
      const auto expectedStressComponent = expectedStress[c][voigtIndex];
      const bool isStressCorrect
          = Dune::FloatCmp::eq<double, Dune::FloatCmp::CmpStyle::absolute>(FEStressComponent, expectedStressComponent);
      t.check(isStressCorrect) << "Stress component " << voigtIndex << " at node " << c << " is not correctly computed"
                               << messageIfFailed;
    }
  }
rath3t commented 1 year ago

Closed with #165