www.mooseframework.org
QuadraturePointMarker.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 "QuadraturePointMarker.h"
11 #include "FEProblem.h"
12 #include "MooseEnum.h"
13 #include "Assembly.h"
14 
15 #include "libmesh/quadrature.h"
16 
19 {
22  MooseEnum third_state("DONT_MARK=-1 COARSEN DO_NOTHING REFINE", "DONT_MARK");
23  params.addParam<MooseEnum>(
24  "third_state",
25  third_state,
26  "The Marker state to apply to values falling in-between the coarsen and refine thresholds.");
27 
28  params.addParam<bool>("invert",
29  false,
30  "If this is true then values _below_ 'refine' will be "
31  "refined and _above_ 'coarsen' will be coarsened.");
32  params.addRequiredParam<VariableName>("variable",
33  "The values of this variable will be compared "
34  "to 'refine' and 'coarsen' to see what should "
35  "be done with the element");
36  return params;
37 }
38 
40  : Marker(parameters),
42  false,
43  "variable",
47  _u(mooseVariableField().sln()),
48  _qrule(_assembly.qRule()),
49  _q_point(_assembly.qPoints()),
50  _qp(0),
51  _third_state(getParam<MooseEnum>("third_state").getEnum<MarkerValue>())
52 {
54 }
55 
58 {
59  MarkerValue current_mark = DONT_MARK;
60 
61  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
62  {
63  MarkerValue new_mark = computeQpMarker();
64 
65  current_mark = std::max(current_mark, new_mark);
66  }
67 
68  return current_mark;
69 }
VarFieldType
Definition: MooseTypes.h:634
Definition: Marker.h:35
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
static InputParameters validParams()
MarkerValue
This mirrors the main refinement flag values in libMesh in Elem::RefinementState but adds "dont_mark"...
Definition: Marker.h:53
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...
auto max(const L &left, const R &right)
unsigned int _qp
The current quadrature point.
static InputParameters validParams()
const QBase *const & _qrule
The quadrature rule for the system.
VarKindType
Framework-wide stuff.
Definition: MooseTypes.h:627
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type It sho...
Definition: MooseEnum.h:31
QuadraturePointMarker(const InputParameters &parameters)
void addMooseVariableDependency(MooseVariableFieldBase *var)
Call this function to add the passed in MooseVariableFieldBase as a variable that this object depends...
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
An interface for accessing Materials.
const std::set< BoundaryID > EMPTY_BOUNDARY_IDS
Definition: MooseTypes.h:597
Interface for objects that need to get values of MooseVariables.
MooseVariableField< Real > & mooseVariableField()
Return the MooseVariableField object that this interface acts on.
virtual MarkerValue computeElementMarker() override
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...
MOOSE now contains C++17 code, so give a reasonable error message stating what the user can do to add...
static InputParameters validParams()
Definition: Marker.C:19
virtual MarkerValue computeQpMarker()=0
Override this to compute a marker value at each quadrature point.