Chombo + EB  3.2
RelaxSolver.H
Go to the documentation of this file.
1 #ifdef CH_LANG_CC
2 /*
3  * _______ __
4  * / ___/ / ___ __ _ / / ___
5  * / /__/ _ \/ _ \/ V \/ _ \/ _ \
6  * \___/_//_/\___/_/_/_/_.__/\___/
7  * Please refer to Copyright.txt, in Chombo's root directory.
8  */
9 #endif
10
11 // BVS, Feb 7, 2005
12
13 #ifndef _RELAXSOLVER_H_
14 #define _RELAXSOLVER_H_
15
16 #include "LinearSolver.H"
17 #include "parstream.H"
18 #include "CH_Timer.H"
20
21 ///
22 /**
23  Iterative solver which only uses the operator's preconditioner.
24  Probably only useful as a bottom solver.
25  */
26 template <class T>
27 class RelaxSolver : public LinearSolver<T>
28 {
29 public:
30  RelaxSolver();
31  virtual ~RelaxSolver();
32
33  ///
34  /**
35  Set whether the solver uses only homogeneous boundary conditions
36  */
37  virtual void setHomogeneous(bool a_homogeneous)
38  {
39  m_homogeneous = a_homogeneous;
40  }
41
42  ///
43  /**
44  Define the solver. a_op is the linear operator.
45  */
46  virtual void define(LinearOp<T>* a_op, bool a_homogeneous);
47
48  ///
49  /**
50  Solve L(phi) = rho
51  */
52  virtual void solve(T& a_phi, const T& a_rhs);
53
54  ///
55  /**
56  public member data: whether or not to use inhomogeneous boundary conditions.
57  */
59
60  ///
61  /**
62  public member data: linear operator.
63  */
65
66  ///
67  /**
68  public member data: maximum number of iterations
69  */
70  int m_imax;
71
72  ///
73  /**
74  public member data: which type of norm to use (default is 2-norm)
75  */
77
78  ///
79  /**
80  public member data: how much screen output user wants
81  */
83
84  ///
85  /**
86  public member data: solver tolerance
87  */
89
90  ///
91  /**
92  public member data: hang when min(previous residual norms) < m_hang * current residual norm (default is zero -- don't check for hanging)
93  */
95  private:
97 };
98
99 // *******************************************************
100 // RelaxSolver Implementation
101 // *******************************************************
102
103 template <class T>
105  :m_homogeneous(false),
106  m_op(NULL),
107  m_imax(40),
108  m_normType(2),
109  m_verbosity(0),
110  m_eps(1.0E-6),
111  m_hang(0.0)
112 {
113 }
114
115 template <class T>
117 {
118 }
119
120 template <class T>
121 void RelaxSolver<T>::define(LinearOp<T>* a_operator, bool a_homogeneous)
122 {
123  m_homogeneous = a_homogeneous;
124  m_op = a_operator;
125 }
126
127 template <class T>
128 void RelaxSolver<T>::solve(T& a_phi, const T& a_rhs)
129 {
130  CH_assert(m_op != NULL);
131  CH_TIME("RelaxSolver::solve");
132  T r, e;
133
134  m_op->create(r, a_rhs);
135  m_op->create(e, a_phi);
136
137  m_op->setToZero(r); // added by petermc, 27 Nov 2013, to zero out ghosts
138  m_op->residual(r, a_phi, a_rhs, m_homogeneous);
139
140  Real norm = m_op->norm(r, m_normType);
141  m_minNorm = norm;
142
143  if (m_verbosity >= 3)
144  {
145  pout() << " RelaxSolver:: initial Residual norm = "
146  << norm << "\n";
147  }
148
149  if (norm > 0.) // Do not iterate if norm == 0.
150  {
151  Real initialNorm = norm;
152  int iter =0;
153  while (iter < m_imax)
154  {
155  m_op->setToZero(e);
156  m_op->preCond(e, r);
157  m_op->incr(a_phi, e, 1.0);
158  m_op->residual(r, a_phi, a_rhs, m_homogeneous);
159  Real oldNorm = norm;
160  norm = m_op->norm(r, m_normType);
161  iter++;
162
163  if (m_verbosity > 3)
164  {
165  pout() << " RelaxSolver:: iteration = " << iter
166  << " residual norm = " << norm
167  << ", rate = " << oldNorm/norm << endl;
168  }
169  if (norm < m_eps*initialNorm) break;
170
171  if (norm < m_minNorm)
172  {
173  m_minNorm = norm;
174  }
175  else if (norm * m_hang > m_minNorm)
176  {
177  if (m_verbosity > 3)
178  {
179  pout() << " RelaxSolver:: exit due to solver hang " << endl;
180  }
181  break;
182  }
183  }
184  }
185
186  m_op->clear(r);
187  m_op->clear(e);
188 }
189
190 #include "NamespaceFooter.H"
191 #endif /*_RELAXSOLVER_H_*/
std::ostream & pout()
Use this in place of std::cout for program output.
Definition: LinearSolver.H:28
#define CH_assert(cond)
Definition: CHArray.H:37
virtual void define(LinearOp< T > *a_op, bool a_homogeneous)
Definition: RelaxSolver.H:121
Real norm(const BoxLayoutData< FArrayBox > &A, const Interval &interval, const int &p=2)
virtual void setHomogeneous(bool a_homogeneous)
Definition: RelaxSolver.H:37
int m_normType
Definition: RelaxSolver.H:76
#define CH_TIME(name)
Definition: CH_Timer.H:82
Real m_eps
Definition: RelaxSolver.H:88
double Real
Definition: REAL.H:33
Real m_hang
Definition: RelaxSolver.H:94
virtual void solve(T &a_phi, const T &a_rhs)
Definition: RelaxSolver.H:128
Real m_minNorm
Definition: RelaxSolver.H:96
bool m_homogeneous
Definition: RelaxSolver.H:58
Definition: LinearSolver.H:156
int m_verbosity
Definition: RelaxSolver.H:82
Definition: RelaxSolver.H:27
virtual ~RelaxSolver()
Definition: RelaxSolver.H:116
int m_imax
Definition: RelaxSolver.H:70
LinearOp< T > * m_op
Definition: RelaxSolver.H:64
RelaxSolver()
Definition: RelaxSolver.H:104