Liebe Gitlab-Nutzerin, lieber Gitlab-Nutzer,
es ist nun möglich sich mittels des ZIH-Logins/LDAP an unserem Dienst anzumelden. Die Konten der externen Nutzer:innen sind über den Reiter "Standard" erreichbar.
Die Administratoren


Dear Gitlab user,
it is now possible to log in to our service using the ZIH login/LDAP. The accounts of external users can be accessed via the "Standard" tab.
The administrators

periodic.cc 2.81 KB
Newer Older
Praetorius, Simon's avatar
Praetorius, Simon committed
1 2
#include <config.h>

3 4 5 6 7 8 9 10 11
#include <iostream>

#if HAVE_DUNE_SPGRID
#include <dune/grid/spgrid.hh>
#endif

#include <amdis/AMDiS.hpp>
#include <amdis/Integrate.hpp>
#include <amdis/ProblemStat.hpp>
12
#include <amdis/LocalOperators.hpp>
13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43

using namespace AMDiS;


template <class Grid>
void print(Grid const& grid)
{
  auto const& indexSet = grid.leafIndexSet();
  std::cout << "vertices:\n";
  for (auto const& v : vertices(grid.leafGridView())) {
    auto coord = v.geometry().corner(0);
    std::cout << indexSet.index(v) << ": [" << coord[0] << ", " << coord[1] << "]\n";
  }
  std::cout << "\n";

  std::cout << "elements:\n";
  for (auto const& e : elements(grid.leafGridView())) {
    std::cout << indexSet.index(e) << ": {" << indexSet.subIndex(e,0,2);
    for (unsigned int i = 1; i < e.subEntities(2); ++i) {
      std::cout << ", " << indexSet.subIndex(e,i,2);
     }
     std::cout << "}\n";
  }
  std::cout << "\n";

  std::cout << "boundary-intersections:\n";
  for (auto const& e : elements(grid.leafGridView())) {
    for (auto const& inter : intersections(grid.leafGridView(), e)) {
      if (!inter.boundary())
        continue;

44 45
      std::cout << inter.boundarySegmentIndex() << ": [" << inter.geometry().center()
                << "] {" << indexSet.index(inter.inside());
46 47 48 49 50 51 52 53 54 55
      if (inter.neighbor())
        std::cout << ", " << indexSet.index(inter.outside());
      std::cout << "}\n";
    }
  }
}

template <class Grid>
void run(Grid& grid)
{
56
  using Traits = LagrangeBasis<Grid, 1>;
57 58 59

  ProblemStat<Traits> prob("ellipt", grid);
  prob.initialize(INIT_ALL);
60
  prob.boundaryManager()->setBoxBoundary({-1,-1,1,1});
61 62 63

  print(grid);

64 65 66 67 68 69 70
  using M = typename ProblemStat<Traits>::SystemMatrix;
  using X = typename ProblemStat<Traits>::SolutionVector;
  using B = typename ProblemStat<Traits>::SystemVector;
  auto periodicBC = makePeriodicBC<M,X,B>(*prob.globalBasis(),
                                          {*prob.boundaryManager(), -1},
                                          {{{1.0,0.0}, {0.0,1.0}}, {1.0, 0.0}});
  periodicBC.init();
71 72 73 74 75 76 77 78 79 80 81 82 83
  std::cout << "periodicNodes:\n";
  for (std::size_t i = 0; i < periodicBC.periodicNodes().size(); ++i)
    std::cout << i << ": " << periodicBC.periodicNodes()[i] << "\n";
  std::cout << "\n";

  std::cout << "associations:\n";
  for (auto const& a : periodicBC.associations())
    std::cout << a.first << " -> " << a.second << "\n";
  std::cout << "\n";
}

int main(int argc, char** argv)
{
84
  Environment env(argc, argv);
85 86 87 88 89 90 91 92 93 94 95 96 97 98

#if HAVE_DUNE_SPGRID
  Dune::SPCube<double,2> cube({0.0,0.0},{1.0,1.0});
  Dune::SPDomain<double,2> domain({cube}, Dune::SPTopology<2>(0b01));
  Dune::SPGrid<double,2> grid1(domain, Dune::SPMultiIndex<2>({2,2}));
  run(grid1);
#endif

  // NOTE: 'periodic' YaspGrid not supported yet, but a simple YaspGrid can be made periodic
  Dune::YaspGrid<2> grid2({1.0,1.0}, {2,2});
  run(grid2);

  return 0;
}