www.mooseframework.org
EigenExecutionerBase.C
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://www.mooseframework.org
3 //*
4 //* All rights reserved, see COPYRIGHT for full restrictions
5 //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 //*
7 //* Licensed under LGPL 2.1, please see LICENSE for details
8 //* https://www.gnu.org/licenses/lgpl-2.1.html
9 
10 #include "EigenExecutionerBase.h"
11 
12 // MOOSE includes
13 #include "AuxiliarySystem.h"
14 #include "DisplacedProblem.h"
15 #include "FEProblem.h"
16 #include "MooseApp.h"
17 #include "MooseEigenSystem.h"
18 #include "UserObject.h"
19 
20 template <>
23 {
25  params.addClassDescription("Executioner for Eigen value problems.");
26 
27  params.addRequiredParam<PostprocessorName>("bx_norm", "To evaluate |Bx| for the eigenvalue");
28  params.addParam<PostprocessorName>("normalization", "To evaluate |x| for normalization");
29  params.addParam<Real>("normal_factor", "Normalize x to make |x| equal to this factor");
30  params.addParam<bool>(
31  "output_before_normalization", true, "True to output a step before normalization");
32  params.addParam<bool>("auto_initialization", true, "True to ask the solver to set initial");
33  params.addParam<Real>("time", 0.0, "System time");
34 
35  params.addPrivateParam<bool>("_eigen", true);
36 
37  params.addParamNamesToGroup("normalization normal_factor output_before_normalization",
38  "Normalization");
39  params.addParamNamesToGroup("auto_initialization time", "Advanced");
40 
41  params.addPrivateParam<bool>("_eigen", true);
42 
43  return params;
44 }
45 
46 const Real &
48 {
49  return _source_integral_old;
50 }
51 
53  : Executioner(parameters),
54  _problem(_fe_problem),
55  _eigen_sys(static_cast<MooseEigenSystem &>(_problem.getNonlinearSystemBase())),
56  _eigenvalue(declareRestartableData("eigenvalue", 1.0)),
57  _source_integral(getPostprocessorValue("bx_norm")),
58  _source_integral_old(1),
59  _normalization(isParamValid("normalization")
60  ? getPostprocessorValue("normalization")
61  : getPostprocessorValue("bx_norm")), // use |Bx| for normalization by default
62  _final_timer(registerTimedSection("final", 1))
63 {
64  // FIXME: currently we have to use old and older solution vectors for power iteration.
65  // We will need 'step' in the future.
66  _problem.transient(true);
67 
68  // we want to tell the App about what our system time is (in case anyone else is interested).
69  Real system_time = getParam<Real>("time");
70  _app.setStartTime(system_time);
71 
72  // set the system time
73  _problem.time() = system_time;
74  _problem.timeOld() = system_time;
75 
76  // used for controlling screen print-out
77  _problem.timeStep() = 0;
78  _problem.dt() = 1.0;
79 }
80 
81 void
83 {
86 
87  if (getParam<bool>("auto_initialization"))
88  {
89  // Initialize the solution of the eigen variables
90  // Note: initial conditions will override this if there is any by _problem.initialSetup()
92  }
95 
96  // check when the postprocessors are evaluated
97  const ExecFlagEnum & bx_exec =
98  _problem.getUserObject<UserObject>(getParam<PostprocessorName>("bx_norm")).getExecuteOnEnum();
99  if (!bx_exec.contains(EXEC_LINEAR))
100  mooseError("Postprocessor " + getParam<PostprocessorName>("bx_norm") +
101  " requires execute_on = 'linear'");
102 
103  if (isParamValid("normalization"))
104  _norm_exec = _problem.getUserObject<UserObject>(getParam<PostprocessorName>("normalization"))
105  .getExecuteOnEnum();
106  else
107  _norm_exec = bx_exec;
108 
109  // check if _source_integral has been evaluated during initialSetup()
110  if (!bx_exec.contains(EXEC_INITIAL))
112 
113  if (_source_integral == 0.0)
114  mooseError("|Bx| = 0!");
115 
116  // normalize solution to make |Bx|=_eigenvalue, _eigenvalue at this point has the initialized
117  // value
119 
120  if (_problem.getDisplacedProblem() != NULL)
121  _problem.getDisplacedProblem()->syncSolutions();
122 
123  /* a time step check point */
125 }
126 
127 void
129 {
130  Real consistency_tolerance = 1e-10;
131 
132  // Scale the solution so that the postprocessor is equal to k.
133  // Note: all dependent objects of k must be evaluated on linear!
134  // We have a fix point loop here, in case the postprocessor is a nonlinear function of the scaling
135  // factor.
136  // FIXME: We have assumed this loop always converges.
137  while (std::fabs(k - _source_integral) > consistency_tolerance * std::fabs(k))
138  {
139  // On the first time entering, the _source_integral has been updated properly in
140  // FEProblemBase::initialSetup()
143  std::stringstream ss;
144  ss << std::fixed << std::setprecision(10) << _source_integral;
145  _console << "\n|Bx| = " << ss.str() << std::endl;
146  }
147 }
148 
149 void
151 {
152  // check to make sure that we don't have any time kernels in this simulation
154  mooseError("You have specified time kernels in your steady state eigenvalue simulation");
156  mooseError("You have not specified any eigen kernels in your eigenvalue simulation");
157 }
158 
159 bool
161  unsigned int max_iter,
162  Real pfactor,
163  bool cheb_on,
164  Real tol_eig,
165  bool echo,
166  PostprocessorName xdiff,
167  Real tol_x,
168  Real & k,
169  Real & initial_res)
170 {
171  mooseAssert(max_iter >= min_iter,
172  "Maximum number of power iterations must be greater than or equal to its minimum");
173  mooseAssert(pfactor > 0.0, "Invaid linear convergence tolerance");
174  mooseAssert(tol_eig > 0.0, "Invalid eigenvalue tolerance");
175  mooseAssert(tol_x > 0.0, "Invalid solution norm tolerance");
176 
177  // obtain the solution diff
178  const PostprocessorValue * solution_diff = NULL;
179  if (!xdiff.empty())
180  {
181  solution_diff = &getPostprocessorValueByName(xdiff);
182  const ExecFlagEnum & xdiff_exec = _problem.getUserObject<UserObject>(xdiff).getExecuteOnEnum();
183  if (!xdiff_exec.contains(EXEC_LINEAR))
184  mooseError("Postprocessor " + xdiff + " requires execute_on = 'linear'");
185  }
186 
187  // not perform any iteration when max_iter==0
188  if (max_iter == 0)
189  return true;
190 
191  // turn off nonlinear flag so that RHS kernels opterate on previous solutions
193 
194  // FIXME: currently power iteration use old and older solutions,
195  // so save old and older solutions before they are changed by the power iteration
197  if (_problem.getDisplacedProblem() != NULL)
198  _problem.getDisplacedProblem()->saveOldSolutions();
199 
200  // save solver control parameters to be modified by the power iteration
201  Real tol1 = _problem.es().parameters.get<Real>("linear solver tolerance");
202  unsigned int num1 =
203  _problem.es().parameters.get<unsigned int>("nonlinear solver maximum iterations");
204  Real tol2 = _problem.es().parameters.get<Real>("nonlinear solver relative residual tolerance");
205 
206  // every power iteration is a linear solve, so set nonlinear iteration number to one
207  _problem.es().parameters.set<Real>("linear solver tolerance") = pfactor;
208  // disable nonlinear convergence check
209  _problem.es().parameters.set<unsigned int>("nonlinear solver maximum iterations") = 1;
210  _problem.es().parameters.set<Real>("nonlinear solver relative residual tolerance") = 1 - 1e-8;
211 
212  if (echo)
213  {
214  _console << '\n';
215  _console << " Power iterations starts\n";
216  _console << " ________________________________________________________________________________ "
217  << std::endl;
218  }
219 
220  // some iteration variables
221  Chebyshev_Parameters chebyshev_parameters;
222 
223  std::vector<Real> keff_history;
224  std::vector<Real> diff_history;
225 
226  bool converged;
227 
228  unsigned int iter = 0;
229 
230  // power iteration loop...
231  // Note: |Bx|/k will stay constant one!
232  makeBXConsistent(k);
233  while (true)
234  {
235  if (echo)
236  _console << " Power iteration= " << iter << std::endl;
237 
238  // Important: we do not call _problem.advanceState() because we do not
239  // want to overwrite the old postprocessor values and old material
240  // properties in stateful materials.
243  if (_problem.getDisplacedProblem() != NULL)
244  {
245  _problem.getDisplacedProblem()->nlSys().copyOldSolutions();
246  _problem.getDisplacedProblem()->auxSys().copyOldSolutions();
247  }
248 
249  Real k_old = k;
251 
252  preIteration();
253  _problem.solve();
254  converged = _problem.converged();
255  if (!converged)
256  break;
257  postIteration();
258 
259  // save the initial residual
260  if (iter == 0)
262 
263  // update eigenvalue
265  _eigenvalue = k;
266 
267  if (echo)
268  {
269  // output on screen the convergence history only when we want to and MOOSE output system is
270  // not used
271  keff_history.push_back(k);
272  if (solution_diff)
273  diff_history.push_back(*solution_diff);
274 
275  std::stringstream ss;
276  if (solution_diff)
277  {
278  ss << '\n';
279  ss << " +================+=====================+=====================+\n";
280  ss << " | iteration | eigenvalue | solution_difference |\n";
281  ss << " +================+=====================+=====================+\n";
282  unsigned int j = 0;
283  if (keff_history.size() > 10)
284  {
285  ss << " : : : :\n";
286  j = keff_history.size() - 10;
287  }
288  for (; j < keff_history.size(); j++)
289  ss << " | " << std::setw(14) << j << " | " << std::setw(19) << std::scientific
290  << std::setprecision(8) << keff_history[j] << " | " << std::setw(19) << std::scientific
291  << std::setprecision(8) << diff_history[j] << " |\n";
292  ss << " +================+=====================+=====================+\n" << std::flush;
293  }
294  else
295  {
296  ss << '\n';
297  ss << " +================+=====================+\n";
298  ss << " | iteration | eigenvalue |\n";
299  ss << " +================+=====================+\n";
300  unsigned int j = 0;
301  if (keff_history.size() > 10)
302  {
303  ss << " : : :\n";
304  j = keff_history.size() - 10;
305  }
306  for (; j < keff_history.size(); j++)
307  ss << " | " << std::setw(14) << j << " | " << std::setw(19) << std::scientific
308  << std::setprecision(8) << keff_history[j] << " |\n";
309  ss << " +================+=====================+\n" << std::flush;
310  ss << std::endl;
311  }
312  _console << ss.str();
313  }
314 
315  // increment iteration number here
316  iter++;
317 
318  if (cheb_on)
319  {
320  chebyshev(chebyshev_parameters, iter, solution_diff);
321  if (echo)
322  _console << " Chebyshev step: " << chebyshev_parameters.icheb << std::endl;
323  }
324 
325  if (echo)
326  _console
327  << " ________________________________________________________________________________ "
328  << std::endl;
329 
330  // not perform any convergence check when number of iterations is less than min_iter
331  if (iter >= min_iter)
332  {
333  // no need to check convergence of the last iteration
334  if (iter != max_iter)
335  {
336  Real keff_error = fabs(k_old - k) / k;
337  if (keff_error > tol_eig)
338  converged = false;
339  if (solution_diff)
340  if (*solution_diff > tol_x)
341  converged = false;
342  if (converged)
343  break;
344  }
345  else
346  {
347  converged = false;
348  break;
349  }
350  }
351  }
352 
353  // restore parameters changed by the executioner
354  _problem.es().parameters.set<Real>("linear solver tolerance") = tol1;
355  _problem.es().parameters.set<unsigned int>("nonlinear solver maximum iterations") = num1;
356  _problem.es().parameters.set<Real>("nonlinear solver relative residual tolerance") = tol2;
357 
358  // FIXME: currently power iteration use old and older solutions, so restore them
360  if (_problem.getDisplacedProblem() != NULL)
361  _problem.getDisplacedProblem()->restoreOldSolutions();
362 
363  return converged;
364 }
365 
366 void
368 {
369 }
370 
371 void
373 {
374 }
375 
376 void
378 {
379  if (getParam<bool>("output_before_normalization"))
380  {
381  _problem.timeStep()++;
382  Real t = _problem.time();
385  _problem.time() = t;
386  }
387 
388  Real s = 1.0;
390  {
391  _console << " Cannot let the normalization postprocessor on custom.\n";
392  _console << " Normalization is abandoned!" << std::endl;
393  }
394  else
395  {
397  s = normalizeSolution(force);
398  if (!MooseUtils::absoluteFuzzyEqual(s, 1.0))
399  _console << " Solution is rescaled with factor " << s << " for normalization!" << std::endl;
400  }
401 
402  if ((!getParam<bool>("output_before_normalization")) || !MooseUtils::absoluteFuzzyEqual(s, 1.0))
403  {
404  _problem.timeStep()++;
405  Real t = _problem.time();
408  _problem.time() = t;
409  }
410 
411  {
412  TIME_SECTION(_final_timer)
415  }
416 }
417 
418 Real
420 {
421  if (force)
423 
424  Real factor;
425  if (isParamValid("normal_factor"))
426  factor = getParam<Real>("normal_factor");
427  else
428  factor = _eigenvalue;
429  Real scaling = factor / _normalization;
430 
431  if (!MooseUtils::absoluteFuzzyEqual(scaling, 1.0))
432  {
433  // FIXME: we assume linear scaling here!
435  // update all aux variables and user objects
436 
437  for (const ExecFlagType & flag : _app.getExecuteOnEnum().items())
438  _problem.execute(flag);
439  }
440  return scaling;
441 }
442 
443 void
445 {
446  std::ostringstream ss;
447  ss << '\n';
448  ss << "*******************************************************\n";
449  ss << " Eigenvalue = " << std::fixed << std::setprecision(10) << _eigenvalue << '\n';
450  ss << "*******************************************************";
451 
452  _console << ss.str() << std::endl;
453 }
454 
456  : n_iter(50), fsmooth(2), finit(6), lgac(0), icheb(0), flux_error_norm_old(1), icho(0)
457 {
458 }
459 
460 void
462 {
463  finit = 6;
464  lgac = 0;
465  icheb = 0;
466  flux_error_norm_old = 1;
467  icho = 0;
468 }
469 
470 void
472  unsigned int iter,
473  const PostprocessorValue * solution_diff)
474 {
475  if (!solution_diff)
476  mooseError("solution diff is required for Chebyshev acceleration");
477 
478  if (chebyshev_parameters.lgac == 0)
479  {
480  if (chebyshev_parameters.icho == 0)
481  chebyshev_parameters.ratio = *solution_diff / chebyshev_parameters.flux_error_norm_old;
482  else
483  {
484  chebyshev_parameters.ratio = chebyshev_parameters.ratio_new;
485  chebyshev_parameters.icho = 0;
486  }
487 
488  if (iter > chebyshev_parameters.finit && chebyshev_parameters.ratio >= 0.4 &&
489  chebyshev_parameters.ratio <= 1)
490  {
491  chebyshev_parameters.lgac = 1;
492  chebyshev_parameters.icheb = 1;
493  chebyshev_parameters.error_begin = *solution_diff;
494  chebyshev_parameters.iter_begin = iter;
495  double alp = 2 / (2 - chebyshev_parameters.ratio);
496  std::vector<double> coef(2);
497  coef[0] = alp;
498  coef[1] = 1 - alp;
502  }
503  }
504  else
505  {
506  chebyshev_parameters.icheb++;
507  double gamma = acosh(2 / chebyshev_parameters.ratio - 1);
508  double alp = 4 / chebyshev_parameters.ratio *
509  std::cosh((chebyshev_parameters.icheb - 1) * gamma) /
510  std::cosh(chebyshev_parameters.icheb * gamma);
511  double beta = (1 - chebyshev_parameters.ratio / 2) * alp - 1;
512  /* if (iter<int(chebyshev_parameters.iter_begin+chebyshev_parameters.n_iter))
513  {
514  std::vector<double> coef(3);
515  coef[0] = alp;
516  coef[1] = 1-alp+beta;
517  coef[2] = -beta;
518  _eigen_sys.combineSystemSolution(NonlinearSystem::EIGEN, coef);
519  }
520  else
521  {*/
522  double gamma_new =
523  (*solution_diff / chebyshev_parameters.error_begin) *
524  (std::cosh((chebyshev_parameters.icheb - 1) * acosh(2 / chebyshev_parameters.ratio - 1)));
525  if (gamma_new < 1.0)
526  gamma_new = 1.0;
527 
528  chebyshev_parameters.ratio_new =
529  chebyshev_parameters.ratio / 2 *
530  (std::cosh(acosh(gamma_new) / (chebyshev_parameters.icheb - 1)) + 1);
531  if (gamma_new > 1.01)
532  {
533  chebyshev_parameters.lgac = 0;
534  // chebyshev_parameters.icheb = 0;
535  // if (chebyshev_parameters.icheb>30)
536  // {
537  if (chebyshev_parameters.icheb > 0)
538  {
539  chebyshev_parameters.icho = 1;
540  chebyshev_parameters.finit = iter;
541  }
542  else
543  {
544  chebyshev_parameters.icho = 0;
545  chebyshev_parameters.finit = iter + chebyshev_parameters.fsmooth;
546  }
547  }
548  else
549  {
550  std::vector<double> coef(3);
551  coef[0] = alp;
552  coef[1] = 1 - alp + beta;
553  coef[2] = -beta;
557  }
558  // }
559  }
560  chebyshev_parameters.flux_error_norm_old = *solution_diff;
561 }
562 
563 bool
564 EigenExecutionerBase::nonlinearSolve(Real rel_tol, Real abs_tol, Real pfactor, Real & k)
565 {
566  makeBXConsistent(k);
567 
568  // turn on nonlinear flag so that eigen kernels opterate on the current solutions
570 
571  // set nonlinear solver controls
572  Real tol1 = _problem.es().parameters.get<Real>("nonlinear solver absolute residual tolerance");
573  Real tol2 = _problem.es().parameters.get<Real>("linear solver tolerance");
574  Real tol3 = _problem.es().parameters.get<Real>("nonlinear solver relative residual tolerance");
575 
576  _problem.es().parameters.set<Real>("nonlinear solver absolute residual tolerance") = abs_tol;
577  _problem.es().parameters.set<Real>("nonlinear solver relative residual tolerance") = rel_tol;
578  _problem.es().parameters.set<Real>("linear solver tolerance") = pfactor;
579 
580  // call nonlinear solve
581  _problem.solve();
582 
583  k = _source_integral;
584  _eigenvalue = k;
585 
586  _problem.es().parameters.set<Real>("nonlinear solver absolute residual tolerance") = tol1;
587  _problem.es().parameters.set<Real>("linear solver tolerance") = tol2;
588  _problem.es().parameters.set<Real>("nonlinear solver relative residual tolerance") = tol3;
589 
590  return _problem.converged();
591 }
void scaleSystemSolution(SYSTEMTAG tag, Real scaling_factor)
Scale the solution vector.
A MultiMooseEnum object to hold "execute_on" flags.
Definition: ExecFlagEnum.h:24
virtual Real & time() const
bool absoluteFuzzyEqual(const T &var1, const T2 &var2, const T3 &tol=libMesh::TOLERANCE *libMesh::TOLERANCE)
Function to check whether two variables are equal within an absolute tolerance.
Definition: MooseUtils.h:231
T & getUserObject(const std::string &name, unsigned int tid=0) const
Get the user object by its name.
const ExecFlagType EXEC_CUSTOM
EigenExecutionerBase(const InputParameters &parameters)
Constructor.
virtual void copyOldSolutions()
Shifts the solutions backwards in time.
Definition: SystemBase.C:1052
NonlinearSystemBase & getNonlinearSystemBase()
void addPrivateParam(const std::string &name, const T &value)
These method add a parameter to the InputParameters object which can be retrieved like any other para...
const ExecFlagEnum & getExecuteOnEnum() const
Return the app level ExecFlagEnum, this contains all the available flags for the app.
Definition: MooseApp.h:667
virtual void makeBXConsistent(Real k)
Normalize solution so that |Bx| = k.
const PostprocessorValue & getPostprocessorValueByName(const PostprocessorName &name)
Retrieve the value of the Postprocessor.
const Real & eigenvalueOld()
The old eigenvalue used by inverse power iterations.
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
virtual void onTimestepEnd() override
const ExecFlagType EXEC_TIMESTEP_END
virtual void postIteration()
Override this for actions that should take place after linear solve of each inverse power iteration...
void mooseError(Args &&... args) const
Definition: MooseObject.h:147
virtual void checkIntegrity()
Make sure time kernel is not presented.
void addRequiredParam(const std::string &name, const std::string &doc_string)
This method adds a parameter and documentation string to the InputParameters object that will be extr...
virtual EquationSystems & es() override
void buildSystemDoFIndices(SYSTEMTAG tag=ALL)
Build DoF indices for a system.
void initSystemSolution(SYSTEMTAG tag, Real v)
Initialize the solution vector with a constant value.
bool contains(const std::string &value) const
Contains methods for seeing if a value is in the MultiMooseEnum.
virtual bool converged() override
virtual void execute(const ExecFlagType &exec_type)
Convenience function for performing execution of MOOSE systems.
void eigenKernelOnOld()
Ask eigenkernels to operate on old or current solution vectors.
void combineSystemSolution(SYSTEMTAG tag, const std::vector< Real > &coefficients)
Linear combination of the solution vectors.
const std::set< ExecFlagType > & items() const
Reference the all the available items.
Definition: ExecFlagEnum.h:73
virtual std::shared_ptr< DisplacedProblem > getDisplacedProblem()
void setStartTime(const Real time)
Set the starting time for the simulation.
Definition: MooseApp.C:1042
virtual bool nonlinearSolve(Real rel_tol, Real abs_tol, Real pfactor, Real &k)
Perform nonlinear solve with the initial guess of the solution.
Real PostprocessorValue
MOOSE typedefs.
Definition: MooseTypes.h:154
const ExecFlagEnum & getExecuteOnEnum() const
Return the execute on MultiMooseEnum for this object.
InputParameters validParams< Executioner >()
Definition: Executioner.C:25
Executioners are objects that do the actual work of solving your problem.
Definition: Executioner.h:32
const ExecFlagType EXEC_LINEAR
AuxiliarySystem & getAuxiliarySystem()
virtual int & timeStep() const
const Real & _normalization
Postprocessor for normalization.
void chebyshev(Chebyshev_Parameters &params, unsigned int iter, const PostprocessorValue *solution_diff)
Real & _eigenvalue
Storage for the eigenvalue computed by the executioner.
Class for containing MooseEnum item information.
Definition: MooseEnumItem.h:21
virtual void restoreOldSolutions()
Restore old solutions from the backup vectors and deallocate them.
bool containsEigenKernel() const
Weather or not the system contains eigen kernels.
MooseEigenSystem & _eigen_sys
MooseApp & _app
The MooseApp this object is associated with.
Definition: MooseObject.h:177
InputParameters validParams< EigenExecutionerBase >()
virtual void printEigenvalue()
Print eigenvalue.
virtual void transient(bool trans)
virtual void saveOldSolutions()
Allocate vectors and save old solutions into them.
void addClassDescription(const std::string &doc_string)
This method adds a description of the class that will be displayed in the input file syntax dump...
virtual Real & timeOld() const
void addParam(const std::string &name, const S &value, const std::string &doc_string)
These methods add an option parameter and a documentation string to the InputParameters object...
virtual Real normalizeSolution(bool force=true)
Normalize the solution vector based on the postprocessor value for normalization. ...
virtual void initialSetup()
const ConsoleStream _console
An instance of helper class to write streams to the Console objects.
void initSystemSolutionOld(SYSTEMTAG tag, Real v)
virtual void solve() override
bool isParamValid(const std::string &name) const
Test if the supplied parameter is valid.
Definition: MooseObject.h:89
virtual void postExecute() override
Override this for actions that should take place after the main solve.
virtual void init() override
Initialize the executioner.
virtual Real & dt() const
const ExecFlagType EXEC_FINAL
Base class for user-specific data.
Definition: UserObject.h:37
virtual void preIteration()
Override this for actions that should take place before linear solve of each inverse power iteration...
virtual void outputStep(ExecFlagType type)
Output the current step.
virtual bool inversePowerIteration(unsigned int min_iter, unsigned int max_iter, Real pfactor, bool cheb_on, Real tol_eig, bool echo, PostprocessorName xdiff, Real tol_x, Real &k, Real &initial_res)
Perform inverse power iterations with the initial guess of the solution.
void addParamNamesToGroup(const std::string &space_delim_names, const std::string group_name)
This method takes a space delimited list of parameter names and adds them to the specified group name...
const ExecFlagType EXEC_INITIAL