SpectralExecutionerDiffusion

Executioner for spectral solve of diffusion equation

Diffusion spectral solve

SpectralExecutionerDiffusion inherits from SpectralExecutionerBase and implements problem-specific methods. For example, it computes the diffusion equation's Green's function for a time step _dt through the method getGreensFunction. This function is used in the execute method to explicitly obtain the diffusion solution step by step. The Green's function is transformed to the reciprocal space where it is point-wise multiplied by the diffused variable in the Fourier. The transformation of its result into the real space carries the proper grid-related scale factor through the FFTBufferBase method backward.

As shown below, this diffusion execution does not need to iterate due to its linear nature. Spectral executioners are, however, not restricted to linear solves.

void
SpectralExecutionerDiffusion::executeExplicit()
{
  unsigned int thisStep = 0;

  _time_step = thisStep;
  _time = _time_step;
  _fe_problem.outputStep(EXEC_INITIAL);
  _fe_problem.advanceState();

  auto & u_buffer = getFFTBuffer<Real>("u");
  auto u = u_buffer.realSpace();

  for (unsigned int step_no = 0; step_no < _nsteps; step_no++)
  {
    // Forward transform to get \hat{u}_{k}
    u_buffer.forwardRaw();

    advanceDiffusionStep(u_buffer, _diff_coeff);

    u_buffer.backward();

    thisStep++;
    _t_current += _dt;
    _time_step = thisStep;

    _fe_problem.execute(EXEC_FINAL);
    _time = _t_current;

    Moose::out << "_t_current: " << _t_current << ". \n";

    _fe_problem.outputStep(EXEC_FINAL);

    if (step_no != _nsteps - 1)
      _fe_problem.advanceState();
  }
}

void
SpectralExecutionerDiffusion::getGreensFunction(FFTBufferBase<Real> & greens,
                                                Real time,
                                                const Real D)
{
  Real accGreens = 0.0;

  const Point & box_size = greens.getBoxSize();

  const auto & grid = greens.grid();
  switch (greens.dim())
  {
    case 1:
    {
      mooseError("Error: Dimension 1 not implemented yet.");
      break;
    }

    case 2:
    {
      std::size_t index = 0;

      const int ni = grid[0];
      const int nj = grid[1]; // Real space.

      for (int i = 0; i < ni; ++i)
        for (int j = 0; j < nj; ++j)
        {

          Real x1 = Real(i) / Real(ni) * box_size(0);
          Real y1 = Real(j) / Real(nj) * box_size(1);

          Real x2 = box_size(0) - x1;
          Real y2 = box_size(1) - y1;

          Moose::out << "box size: " << box_size(0) << ", " << box_size(1) << "\n";
          Moose::out << "number of cells: " << ni << ", " << nj << "\n";

          if (time == 0)
            mooseError("Greens function undefined at t = 0s in the FFT solver.");
          else
          {
            Real sum = std::exp(-(x1 * x1 + y1 * y1) / 4 / D / time) +
                       std::exp(-(x1 * x1 + y2 * y2) / 4 / D / time) +
                       std::exp(-(x2 * x2 + y1 * y1) / 4 / D / time) +
                       std::exp(-(x2 * x2 + y2 * y2) / 4 / D / time);

            if (sum == 0 && i != 0 && j != 0)
              mooseError("Precision lost in the analytical integration of heat equation");

            Real correction_factor = (Real(ni) / box_size(0)) * (Real(nj) / box_size(1));

            greens.realSpace()[index] =
                (1.0 / (4 * libMesh::pi * D * time)) * sum / correction_factor;
          }
          accGreens += greens.realSpace()[index];
          index++;
        }
      mooseInfo("Accumulated value of Greens function: ", accGreens, "\n");
      break;
    }

    case 3:
    {
      mooseError("Error: Dimension 3 not implemented yet.");
      break;
    }

    default:
      mooseError("Invalid buffer dimension");
  }
}
void
SpectralExecutionerDiffusion::execute()
{
  unsigned int thisStep = 0;

  _time_step = thisStep;
  _time = _time_step;
  _fe_problem.outputStep(EXEC_INITIAL);
  _fe_problem.advanceState();

  auto & u_buffer = getFFTBuffer<Real>("u");
  auto & greens = getFFTBuffer<Real>("greens_buffer");

  // Get Greens function of
  getGreensFunction(greens, _dt, _diff_coeff);
  greens.forwardRaw();

  for (unsigned int step_no = 0; step_no < _nsteps; step_no++)
  {
    u_buffer.forwardRaw();

    u_buffer.reciprocalSpace() *= greens.reciprocalSpace();
    // scalarMultiplyBuffer(u_buffer, greens);
    u_buffer.backward();

    // End of diffusion computations
    thisStep++;
    _t_current += _dt;
    _time_step = thisStep;

    _fe_problem.execute(EXEC_FINAL);
    _time = _t_current;

    Moose::out << "_t_current: " << _t_current << ". \n";

    _fe_problem.outputStep(EXEC_FINAL);

    if (step_no != _nsteps - 1)
      _fe_problem.advanceState();
  }
}

void
SpectralExecutionerDiffusion::executeTotalDiffusion()
{
  unsigned int thisStep = 0;

  _time_step = thisStep;
  _time = _time_step;
  _fe_problem.outputStep(EXEC_INITIAL);
  _fe_problem.advanceState();

  FFTBufferBase<Real> & u_buffer = getFFTBuffer<Real>("u");
  FFTBufferBase<Real> & u_init_buffer = getFFTBuffer<Real>("u_initial");

  u_init_buffer.forwardRaw();
  // auto &greens = getFFTBuffer<Real>("greens_buffer");

  for (unsigned int step_no = 0; step_no < _nsteps; step_no++)
  {
    // Forward transform to get \hat{u}_{k}
    u_buffer.forwardRaw();

    advanceDiffusionStepTotal(u_init_buffer, u_buffer, _diff_coeff, _time);

    u_buffer.backward();

    thisStep++;
    _t_current += _dt;
    _time_step = thisStep;

    _fe_problem.execute(EXEC_FINAL);
    _time = _t_current;

    Moose::out << "_t_current: " << _t_current << ". \n";

    _fe_problem.outputStep(EXEC_FINAL);

    if (step_no != _nsteps - 1)
      _fe_problem.advanceState();
  }
}

void
SpectralExecutionerDiffusion::scalarMultiplyBuffer(FFTBufferBase<Real> & u,
                                                   const FFTBufferBase<Real> & greens)
{
  const auto & grid = u.grid();

  switch (u.dim())
  {
    case 1:
    {
      mooseError("Error: Dimension 1 not implemented yet.");
      return;
    }

    case 2:
    {
      std::size_t index = 0;
      const int ni = grid[0];
      const int nj = (grid[1] >> 1) + 1;

      for (int i = 0; i < ni; ++i)
        for (int j = 0; j < nj; ++j)
        {
          u.reciprocalSpace()[index] = greens.reciprocalSpace()[index] * u.reciprocalSpace()[index];
          index++;
        }
      return;
    }

    case 3:
    {
      mooseError("Error: Dimension 3 not implemented yet.");
      return;
    }

    default:
      mooseError("Invalid buffer dimension");
  }
}

void
SpectralExecutionerDiffusion::advanceDiffusionStepTotal(const FFTBufferBase<Real> & u_initial,
                                                        FFTBufferBase<Real> & u,
                                                        const Real D,
                                                        const Real time)
{
  const auto & grid = u.grid();
  switch (u.dim())
  {
    case 1:
    {
      mooseError("Error: Dimension 1 not implemented yet.");
      return;
    }

    case 2:
    {
      std::size_t index = 0;
      const int ni = grid[0];
      const int nj = (grid[1] >> 1) + 1;

      const auto & ivec = u.kTable(0);
      const auto & jvec = u.kTable(1);

      for (int i = 0; i < ni; ++i)
        for (int j = 0; j < nj; ++j)
        {
          u.reciprocalSpace()[index] =
              u_initial.reciprocalSpace()[index] *
              std::exp(-D * (ivec[i] * ivec[i] + jvec[j] * jvec[j]) * time);
          index++;
        }
      return;
    }

    case 3:
    {
      mooseError("Error: Dimension 3 not implemented yet.");
      return;
    }

    default:
      mooseError("Invalid buffer dimension");
  }
}

void
SpectralExecutionerDiffusion::advanceDiffusionStep(FFTBufferBase<Real> & inOut, const Real D)
{
  const Real sizePixel = 1.0; // Assumed same size in all directions

  const auto & grid = inOut.grid();
  switch (inOut.dim())
  {
    case 1:
    {
      mooseError("Error: Dimension 1 not implemented yet.");
      return;
    }

    case 2:
    {
      std::size_t index = 0;
      const int ni = grid[0];
      const int nj = (grid[1] >> 1) + 1;

      for (int i = 0; i < ni; ++i)
        for (int j = 0; j < nj; ++j)
        {
          Real r1 = D * _dt / (sizePixel * sizePixel);
          Real r2 = D * _dt / (sizePixel * sizePixel);

          Real hk1 = (1 - 2 * r1 * (1 - std::cos(libMesh::pi * i / ((ni) / 2))));
          Real hk2 = (1 - 2 * r2 * (1 - std::cos(libMesh::pi * j / ((nj) / 2))));
          inOut.reciprocalSpace()[index] *= hk1 * hk2;
          index++;
        }
      return;
    }

    case 3:
    {
      mooseError("Error: Dimension 3 not implemented yet.");
      return;
    }

    default:
      mooseError("Invalid buffer dimension");
  }
}
(src/executioners/SpectralExecutionerDiffusion.C)

Input Parameters

  • diffusion_coefficient1Diffusion coefficient

    Default:1

    C++ Type:double

    Unit:(no unit assumed)

    Controllable:No

    Description:Diffusion coefficient

  • number_steps0Time step for ODE integration

    Default:0

    C++ Type:unsigned int

    Controllable:No

    Description:Time step for ODE integration

  • time0System time

    Default:0

    C++ Type:double

    Unit:(no unit assumed)

    Controllable:No

    Description:System time

  • time_step1Time step for ODE integration

    Default:1

    C++ Type:double

    Unit:(no unit assumed)

    Controllable:No

    Description:Time step for ODE integration

  • verboseFalseSet to true to print additional information

    Default:False

    C++ Type:bool

    Controllable:No

    Description:Set to true to print additional information

Optional Parameters

  • accept_on_max_fixed_point_iterationFalseTrue to treat reaching the maximum number of fixed point iterations as converged.

    Default:False

    C++ Type:bool

    Controllable:No

    Description:True to treat reaching the maximum number of fixed point iterations as converged.

  • auto_advanceFalseWhether to automatically advance sub-applications regardless of whether their solve converges, for transient executioners only.

    Default:False

    C++ Type:bool

    Controllable:No

    Description:Whether to automatically advance sub-applications regardless of whether their solve converges, for transient executioners only.

  • custom_abs_tol1e-50The absolute nonlinear residual to shoot for during fixed point iterations. This check is performed based on postprocessor defined by the custom_pp residual.

    Default:1e-50

    C++ Type:double

    Unit:(no unit assumed)

    Controllable:No

    Description:The absolute nonlinear residual to shoot for during fixed point iterations. This check is performed based on postprocessor defined by the custom_pp residual.

  • custom_ppPostprocessor for custom fixed point convergence check.

    C++ Type:PostprocessorName

    Unit:(no unit assumed)

    Controllable:No

    Description:Postprocessor for custom fixed point convergence check.

  • custom_rel_tol1e-08The relative nonlinear residual drop to shoot for during fixed point iterations. This check is performed based on the postprocessor defined by custom_pp residual.

    Default:1e-08

    C++ Type:double

    Unit:(no unit assumed)

    Controllable:No

    Description:The relative nonlinear residual drop to shoot for during fixed point iterations. This check is performed based on the postprocessor defined by custom_pp residual.

  • direct_pp_valueFalseTrue to use direct postprocessor value (scaled by value on first iteration). False (default) to use difference in postprocessor value between fixed point iterations.

    Default:False

    C++ Type:bool

    Controllable:No

    Description:True to use direct postprocessor value (scaled by value on first iteration). False (default) to use difference in postprocessor value between fixed point iterations.

  • disable_fixed_point_residual_norm_checkFalseDisable the residual norm evaluation thus the three parameters fixed_point_rel_tol, fixed_point_abs_tol and fixed_point_force_norms.

    Default:False

    C++ Type:bool

    Controllable:No

    Description:Disable the residual norm evaluation thus the three parameters fixed_point_rel_tol, fixed_point_abs_tol and fixed_point_force_norms.

  • fixed_point_abs_tol1e-50The absolute nonlinear residual to shoot for during fixed point iterations. This check is performed based on the main app's nonlinear residual.

    Default:1e-50

    C++ Type:double

    Unit:(no unit assumed)

    Controllable:No

    Description:The absolute nonlinear residual to shoot for during fixed point iterations. This check is performed based on the main app's nonlinear residual.

  • fixed_point_algorithmpicardThe fixed point algorithm to converge the sequence of problems.

    Default:picard

    C++ Type:MooseEnum

    Options:picard, secant, steffensen

    Controllable:No

    Description:The fixed point algorithm to converge the sequence of problems.

  • fixed_point_force_normsFalseForce the evaluation of both the TIMESTEP_BEGIN and TIMESTEP_END norms regardless of the existence of active MultiApps with those execute_on flags, default: false.

    Default:False

    C++ Type:bool

    Controllable:No

    Description:Force the evaluation of both the TIMESTEP_BEGIN and TIMESTEP_END norms regardless of the existence of active MultiApps with those execute_on flags, default: false.

  • fixed_point_max_its1Specifies the maximum number of fixed point iterations.

    Default:1

    C++ Type:unsigned int

    Controllable:No

    Description:Specifies the maximum number of fixed point iterations.

  • fixed_point_min_its1Specifies the minimum number of fixed point iterations.

    Default:1

    C++ Type:unsigned int

    Controllable:No

    Description:Specifies the minimum number of fixed point iterations.

  • fixed_point_rel_tol1e-08The relative nonlinear residual drop to shoot for during fixed point iterations. This check is performed based on the main app's nonlinear residual.

    Default:1e-08

    C++ Type:double

    Unit:(no unit assumed)

    Controllable:No

    Description:The relative nonlinear residual drop to shoot for during fixed point iterations. This check is performed based on the main app's nonlinear residual.

  • multiapp_fixed_point_convergenceName of the Convergence object to use to assess convergence of the MultiApp fixed point solve. If not provided, a default Convergence will be constructed internally from the executioner parameters.

    C++ Type:ConvergenceName

    Controllable:No

    Description:Name of the Convergence object to use to assess convergence of the MultiApp fixed point solve. If not provided, a default Convergence will be constructed internally from the executioner parameters.

  • relaxation_factor1Fraction of newly computed value to keep.Set between 0 and 2.

    Default:1

    C++ Type:double

    Unit:(no unit assumed)

    Controllable:No

    Description:Fraction of newly computed value to keep.Set between 0 and 2.

  • transformed_postprocessorsList of main app postprocessors to transform during fixed point iterations

    C++ Type:std::vector<PostprocessorName>

    Unit:(no unit assumed)

    Controllable:No

    Description:List of main app postprocessors to transform during fixed point iterations

  • transformed_variablesList of main app variables to transform during fixed point iterations

    C++ Type:std::vector<std::string>

    Controllable:No

    Description:List of main app variables to transform during fixed point iterations

Fixed Point Iterations Parameters

  • control_tagsAdds user-defined labels for accessing object parameters via control logic.

    C++ Type:std::vector<std::string>

    Controllable:No

    Description:Adds user-defined labels for accessing object parameters via control logic.

  • enableTrueSet the enabled status of the MooseObject.

    Default:True

    C++ Type:bool

    Controllable:No

    Description:Set the enabled status of the MooseObject.

  • outputsVector of output names where you would like to restrict the output of variables(s) associated with this object

    C++ Type:std::vector<OutputName>

    Controllable:No

    Description:Vector of output names where you would like to restrict the output of variables(s) associated with this object

Advanced Parameters

  • max_xfem_update4294967295Maximum number of times to update XFEM crack topology in a step due to evolving cracks

    Default:4294967295

    C++ Type:unsigned int

    Controllable:No

    Description:Maximum number of times to update XFEM crack topology in a step due to evolving cracks

  • update_xfem_at_timestep_beginFalseShould XFEM update the mesh at the beginning of the timestep

    Default:False

    C++ Type:bool

    Controllable:No

    Description:Should XFEM update the mesh at the beginning of the timestep

Xfem Fixed Point Iterations Parameters

    Restart Parameters

    Input Files

    References

    No citations exist within this document.