HL_SignedDistTraverse.h 5.49 KB
Newer Older
1
2
3
4
5
6
7
/******************************************************************************
 *
 * AMDiS - Adaptive multidimensional simulations
 *
 * Copyright (C) 2013 Dresden University of Technology. All Rights Reserved.
 * Web: https://fusionforge.zih.tu-dresden.de/projects/amdis
 *
8
 * Authors:
9
10
11
12
13
14
15
16
17
 * Simon Vey, Thomas Witkowski, Andreas Naumann, Simon Praetorius, et al.
 *
 * This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE
 * WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
 *
 *
 * This file is part of AMDiS
 *
 * See also license.opensource.txt in the distribution.
18
 *
19
 ******************************************************************************/
20
21
22
23




24
25
26
27
28
29
30
31
32
33
34
35
36
37
#ifndef HL_SIGNEDDISTTRAVERSE
#define HL_SIGNEDDISTTRAVERSE

#include "ElInfo.h"
#include "FixVec.h"
#include "Traverse.h"
#include "ElementLevelSet.h"
#include "BoundaryElementDist.h"
#include "ElementUpdate.h"
#include "ElementUpdate_2d.h"
#include "ElementUpdate_3d.h"
#include "HL_SignedDist.h"
#include "VelocityExt.h"

38
39
namespace reinit {

40
41
using namespace AMDiS;

42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
#ifdef HAVE_PARALLEL_DOMAIN_AMDIS
  /**
   * Functor which assigns a the maximum of a and b.
   * (Used for synchronizing DOFVectors on interior boundaries in parallel AMDiS)
   */
  struct max_assigner: FunctorBase
  {
    int getDegree(int d0) const { return d0; }
    static void eval(double& v0, const double& v1) { v0 = std::max(v0, v1); }
    void operator()(double &v0, const double& v1) { eval(v0, v1); }
  };
#endif

#ifdef HAVE_PARALLEL_DOMAIN_AMDIS
  /**
   * Functor which assigns a the minimum of a and b with respect to 0.
   * (Used for synchronizing DOFVectors on interior boundaries in parallel AMDiS)
   */
  struct min_to_zero_assigner: FunctorBase
  {
    int getDegree(int d0) const { return d0; }
    static void eval(double& v0, const double& v1) { v0 = (std::abs(v0) < std::abs(v1) ? (v0) : (v1)); }
    void operator()(double &v0, const double& v1) { eval(v0, v1); }
  };
#endif

68
69
class HL_SignedDistTraverse : public HL_SignedDist
{
70
public:
71

72
  HL_SignedDistTraverse(const char *name_,
73
74
75
76
			int dim_,
			bool doVelocityExt = false,
			Flag velExtType_ = VEL_EXT)
    : HL_SignedDist(name_, dim_, doVelocityExt, velExtType_),
77
78
      sDOld_DOF(NULL),
      update_DOF(NULL),
79
80
81
      tol_reached(false),
      elVert(dim_, NO_INIT),
      uhVal(dim_, NO_INIT)
82
83
84
85
  {
    FUNCNAME("HL_SignedDistTraverse::HL_SignedDistTraverse");

    // ===== Read parameters from init file. =====
86
87
88
    Parameters::get(name + "->tolerance", tol);
    Parameters::get(name + "->maximal number of iteration steps", maxIt);
    Parameters::get(name + "->Gauss-Seidel iteration", GaussSeidelFlag);
89
    Parameters::get(name + "->periodic boundaries", boundaryTypes);
90
91
92
93
94
95
96

    TEST_EXIT(tol > 0)("illegal tolerance !\n");

    // ---> for test purposes: initialization of counter variables
    calcUpdate_Cntr = 0;
    setUpdate_Cntr = 0;
    // ---> end: for test purposes
97
  }
98

99

100
101
102
  ~HL_SignedDistTraverse()
  {
    if (sDOld_DOF)
103
      delete sDOld_DOF;
104
105
106
107

    // ---> for test purposes: print result of update counting
    printUpdateCntr();
    // ---> end: for test purposes
108
  }
109
110
111

 protected:
  /**
112
113
114
   * Initializes the boundary: calculation of the distance of boundary
   * vertices to the interface.
   * Interface is given by lS_DOF and result is stored in
115
116
117
118
119
120
121
122
123
124
   * sD_DOF.
   */
  void initializeBoundary();

  /**
   * Calculates the distance function and stores result in sD_DOF.
   * Requirement: The boundary values are already set in sD_DOF.
   */
  void HL_updateIteration();

125
  /// Hopf-Lax element update on element elInfo.
126
127
  void HL_elementUpdate(ElInfo *elInfo);

128
  /// Calculate the update for the vertex xh of element elInfo.
129
  double calcElementUpdate(ElInfo *elInfo,
130
131
132
133
134
135
136
137
138
139
140
141
142
			   int nXh,
			   const DegreeOfFreedom *locInd);

  /**
   * In iteration loop: checks whether tolerance is reached.
   *    true  -  tolerance is reached
   *    false -  tolerance is not reached
   */
  bool checkTol();

  // ---> for test purposes: print result of update counting
  void printUpdateCntr()
  {
143
144
145
146
147
148
    /* cout << "\n"; */
    /* cout << "\tUpdate statistic: \n"; */
    /* cout << "\t-----------------\n"; */
    /* cout << "\tcalculated updates: " << calcUpdate_Cntr << "\n"; */
    /* cout << "\tset updates: " << setUpdate_Cntr << "\n"; */
    /* cout << "\n"; */
149
  }
150
151
152
153
  // ---> end: for test purposes

 protected:
  /**
154
   * DOF vector for the last iteration step.
155
156
157
158
159
160
161
162
163
164
165
166
   * Used during Jacobi iteration and to check whether tolerance is reached.
   */
  DOFVector<double> *sDOld_DOF;

  /**
   * Pointer to the DOF vector holding the values of the last iteration
   * step.
   *      Gauss-Seidel iteration  -  update_DOF == sD_DOF
   *      Jacobi iteration        -  update_DOF == sDOld_DOF
   */
  DOFVector<double> *update_DOF;

167
  /// Tolerance for Hopf-Lax update iteration loop.
168
169
  double tol;

170
  /// Flag showing whether tolerance tol is reached.
171
172
  bool tol_reached;

173
  /// Maximal number of mesh iterations for Hopf-Lax update.
174
175
176
177
178
179
180
181
182
  int maxIt;

  /**
   * Indicates whether Gauss-Seidel or Jacobi iteration is used.
   *     0 - Jacobi
   *   !=0 - Gauss-Seidel
   */
  int GaussSeidelFlag;

183
184
185
186
187
188
  std::vector<DegreeOfFreedom> locInd;

  FixVec<WorldVector<double> *, VERTEX> elVert;

  FixVec<double, VERTEX> uhVal;

189
190
191
192
193
  // ---> for test purposes: variables to count calculated and set updates
  int calcUpdate_Cntr;
  int setUpdate_Cntr;
  // ---> end: for test purposes

194
195
  // boundary types for periodic mesh boundaries
  std::vector<signed int> boundaryTypes;
196
197
};

198
199
200
201
} // end namespace reinit

using reinit::HL_SignedDistTraverse;

202
#endif  // HL_SIGNEDDISTTRAVERSE