Merge pull request idaholab#44 from jessecarterMOOSE/sinkmapuserobjec…

…t_42

Adds ability to define sink regions on a uniform grid, closes idaholab#42

examples/2d_sink_map_with_refinement/2d_sink_map_with_refinement.i

@@ -0,0 +1,233 @@
+[Mesh]
+  type = GeneratedMesh
+  dim = 2
+  xmin = -1
+  xmax = 1
+  ymin = -1
+  ymax = 1
+  nx = 20
+  ny = 20
+[]
+
+[Variables]
+  [./u]
+    order = FIRST
+    family = LAGRANGE
+
+    [./InitialCondition]
+      type = ConstantIC
+      value = 0
+    [../]
+  [../]
+[]
+
+[AuxVariables]
+  [./sink_map_aux]
+    order = CONSTANT
+    family = MONOMIAL
+  [../]
+[]
+
+[Kernels]
+  [./ie]
+    type = TimeDerivative
+    variable = u
+  [../]
+
+  [./diff]
+    type = MaterialDiffusion
+    variable = u
+    diffusivity_name = diffusivity
+  [../]
+
+  [./sink_map]
+    type = SinkMapKernel
+    variable = u
+    sink_map_user_object = sink_map_uo
+  [../]
+
+  [./event_inserter_source]
+    type = EventInserterSource
+    variable = u
+    inserter = inserter
+    gaussian_user_object = gaussian_uo
+  [../]
+[]
+
+[AuxKernels]
+  [./sink_map]
+    type = SinkMapAux
+    variable = sink_map_aux
+    sink_map_user_object = sink_map_uo
+    execute_on = timestep_end
+  [../]
+[]
+
+[BCs]
+  [./Periodic]
+    [./all]
+      variable = u
+      auto_direction = 'x y'
+    [../]
+  [../]
+[]
+
+[Materials]
+  [./simple]
+    type = GenericConstantMaterial
+    block = 0
+    prop_names = 'diffusivity'
+    prop_values = '2.0'
+  [../]
+[]
+
+[UserObjects]
+  [./inserter_circle_average]
+    type = CircleAverageMaterialProperty
+    mat_prop = diffusivity
+    periodic_variable = u
+    inserter = inserter
+    radius = 0.2
+  [../]
+  [./circle_average]
+    type = CircleAverageMaterialProperty
+    mat_prop = diffusivity
+    periodic_variable = u
+  [../]
+  [./random_point_uo]
+    type = RandomPointUserObject
+    seed = 1
+  [../]
+  [./gaussian_uo]
+    type = GaussianUserObject
+    sigma = 0.05
+    use_random_points = true
+    random_point_user_object = random_point_uo
+    periodic_variable = u
+    scale = 3
+  [../]
+  [./sink_gaussian_uo]
+    type = GaussianUserObject
+    sigma = 0.05
+  [../]
+  [./circle_max_original_element_size_uo]
+    type = CircleMaxOriginalElementSize
+    periodic_variable = u
+  [../]
+  [./inserter]
+    type = EventInserter
+    insert_initial = true
+    random_timing = false
+    mean = 1.2
+    random_point_user_object = random_point_uo
+    seed = 3
+    verbose = true
+    track_old_events = true
+    removal_method = sigma_element_size_ratio
+    removal_sigma_element_size_ratio = 2.0
+    radius = 3.0
+    gaussian_user_object = gaussian_uo
+    circle_average_material_property_user_object = circle_average
+    inserter_circle_average_material_property_user_object = inserter_circle_average
+    circle_max_original_element_size_user_object = circle_max_original_element_size_uo
+  [../]
+  [./sink_map_uo]
+    type = SinkMapUserObject
+    spacing = 1.0
+    strength = 0.05
+    sink_placement = inside
+    gaussian_user_object = sink_gaussian_uo
+  [../]
+[]
+
+[Adaptivity]
+  initial_marker = event_marker
+  initial_steps = 10
+  marker = event_marker
+  cycles_per_step = 10
+  max_h_level = 0
+  recompute_markers_during_cycles = true
+  [./Markers]
+    [./event_marker]
+      type = EventMarker
+      inserter = inserter
+      gaussian_user_object = gaussian_uo
+      marker_radius = 6.0
+      coarsen_events = true
+      verbose = true
+      periodic_variable = u
+      event_sigma_mesh_ratio = 2.0
+      refine_sinks = true
+      sink_marker_radius = 6.0
+      sink_map_user_object = sink_map_uo
+      sink_gaussian_user_object = sink_gaussian_uo
+      sink_sigma_mesh_ratio = 2.0
+    [../]
+  [../]
+[]
+
+[Executioner]
+  type = Transient
+  scheme = 'implicit-euler'
+
+  solve_type = 'NEWTON'
+
+  start_time = 0.0
+  end_time = 50.0
+
+  verbose = true
+
+  petsc_options_iname = '-pc_type -pc_hypre_type'
+  petsc_options_value = 'hypre boomeramg'
+
+  nl_abs_tol = 1.0e-12
+
+  [./TimeStepper]
+    type = EventTimeStepper
+    dt = 0.01
+    event_inserter = inserter
+    growth_factor = 2.0
+  [../]
+[]
+
+[Preconditioning]
+  [./smp]
+    type = SMP
+    full = true
+  [../]
+[]
+
+[Postprocessors]
+  [./dt]
+    type = TimestepSize
+  [../]
+  [./solution_integral]
+    type = ElementIntegralVariablePostprocessor
+    variable = u
+  [../]
+  [./solution_average]
+    type = ElementAverageValue
+    variable = u
+  [../]
+  [./sink_integral]
+    type = ElementIntegralVariablePostprocessor
+    variable = sink_map_aux
+  [../]
+  [./sink_average]
+    type = ElementAverageValue
+    variable = sink_map_aux
+  [../]
+  [./num_elems]
+    type = NumElems
+  [../]
+[]
+
+[Outputs]
+  exodus = true
+  csv = true
+  execute_on = 'initial timestep_end'
+  [./console]
+    type = Console
+    print_mesh_changed_info = true
+  [../]
+[]

examples/2d_sink_map_with_refinement/plot_concentration.py

@@ -0,0 +1,29 @@
+import numpy as np
+import matplotlib.pyplot as plt
+import os
+
+# enter values from simulation here
+diffusivity = 2.0
+sink_strength = 0.05
+mean = 1.2  # mean of EventInserter distribution
+Lx = 2.0    # length of problem in x
+Ly = 2.0    # length of problem in y
+scale = 3.0 # scale of Gaussian source term
+
+# postprocessor output should be the only csv file in directory
+for file in os.listdir("."):
+  if file.endswith(".csv"):
+    csvfile=file
+
+data=np.loadtxt(csvfile, delimiter=',', skiprows=1)
+
+plt.plot(data[:,0], data[:,5], 'r-', label='solution')
+plt.xlabel('time')
+plt.ylabel('average concentration')
+
+analytic = scale/diffusivity/sink_strength/mean/Lx/Ly*(1.0 - np.exp(-sink_strength*diffusivity*data[:,0]))
+plt.plot(data[:,0], analytic, 'k-', label='analytic')
+
+plt.legend(loc="lower right")
+
+plt.show()

examples/2d_sink_map_with_refinement/plot_num_elements.py

@@ -0,0 +1,16 @@
+import numpy as np
+import matplotlib.pyplot as plt
+import os
+
+# postprocessor output should be the only csv file in directory
+for file in os.listdir("."):
+  if file.endswith(".csv"):
+    csvfile=file
+
+data=np.loadtxt(csvfile, delimiter=',', skiprows=1)
+
+plt.plot(data[:,0], data[:,2], 'r-o')
+plt.xlabel('time')
+plt.ylabel('number of elements')
+
+plt.show()

include/auxkernels/SinkMapAux.h

@@ -0,0 +1,42 @@
+/****************************************************************/
+/*               DO NOT MODIFY THIS HEADER                      */
+/* MOOSE - Multiphysics Object Oriented Simulation Environment  */
+/*                                                              */
+/*           (c) 2010 Battelle Energy Alliance, LLC             */
+/*                   ALL RIGHTS RESERVED                        */
+/*                                                              */
+/*          Prepared by Battelle Energy Alliance, LLC           */
+/*            Under Contract No. DE-AC07-05ID14517              */
+/*            With the U. S. Department of Energy               */
+/*                                                              */
+/*            See COPYRIGHT for full restrictions               */
+/****************************************************************/
+
+#ifndef SINKMAPAUX_H
+#define SINKMAPAUX_H
+
+#include "AuxKernel.h"
+
+//Forward Declarations
+class SinkMapAux;
+class SinkMapUserObject;
+
+template<>
+InputParameters validParams<SinkMapAux>();
+
+class SinkMapAux : public AuxKernel
+{
+public:
+  SinkMapAux(const InputParameters & parameters);
+
+  virtual ~SinkMapAux() {}
+
+protected:
+  virtual void precalculateValue();
+  virtual Real computeValue();
+
+  const SinkMapUserObject & _sink_map_uo;
+  std::vector<Real> _element_sink_map;
+};
+
+#endif // SINKMAPAUX_H

include/kernels/SinkMapKernel.h

@@ -0,0 +1,42 @@
+/****************************************************************/
+/*               DO NOT MODIFY THIS HEADER                      */
+/* MOOSE - Multiphysics Object Oriented Simulation Environment  */
+/*                                                              */
+/*           (c) 2010 Battelle Energy Alliance, LLC             */
+/*                   ALL RIGHTS RESERVED                        */
+/*                                                              */
+/*          Prepared by Battelle Energy Alliance, LLC           */
+/*            Under Contract No. DE-AC07-05ID14517              */
+/*            With the U. S. Department of Energy               */
+/*                                                              */
+/*            See COPYRIGHT for full restrictions               */
+/****************************************************************/
+#ifndef SINKMAPKERNEL_H
+#define SINKMAPKERNEL_H
+
+#include "Reaction.h"
+
+//Forward Declarations
+class SinkMapKernel;
+class SinkMapUserObject;
+
+template<>
+InputParameters validParams<SinkMapKernel>();
+
+class SinkMapKernel : public Reaction
+{
+public:
+  SinkMapKernel(const InputParameters & parameters);
+
+protected:
+  virtual void computeResidual() override;
+  virtual Real computeQpResidual() override;
+  virtual void computeJacobian() override;
+  virtual Real computeQpJacobian() override;
+
+  const SinkMapUserObject & _sink_map_uo;
+  const MaterialProperty<Real> & _diffusivity;
+  std::vector<Real> _element_sink_map;
+};
+
+#endif //SINKMAPKERNEL_H