Skip to content

Commit

Permalink
Add variable shim for FFT (idaholab#401)
Browse files Browse the repository at this point in the history
  • Loading branch information
dschwen committed Jan 20, 2020
1 parent e16e74f commit 077ca17
Show file tree
Hide file tree
Showing 2 changed files with 196 additions and 0 deletions.
175 changes: 175 additions & 0 deletions include/variables/MooseFFTVariable.h
Original file line number Diff line number Diff line change
@@ -0,0 +1,175 @@
/**********************************************************************/
/* DO NOT MODIFY THIS HEADER */
/* MAGPIE - Mesoscale Atomistic Glue Program for Integrated Execution */
/* */
/* Copyright 2017 Battelle Energy Alliance, LLC */
/* ALL RIGHTS RESERVED */
/**********************************************************************/

#pragma once

#include "MooseVariableFEBase.h"

class MooseFFTVariable : public MooseVariableFEBase
{
public:
static InputParameters validParams();

MooseFFTVariable(const InputParameters & parameters);

virtual bool isNodal() const { return false; }

/**
* Clear out the dof indices. We do this in case this variable is not going to be prepared at
* all...
*/
virtual void clearDofIndices() {}

/**
* Prepare the elemental degrees of freedom
*/
virtual void prepare() {}

/**
* Prepare the neighbor element degrees of freedom
*/
virtual void prepareNeighbor() { mooseError("Neighbor FFT variables are not supported yet."); }

/**
* Prepare a lower dimensional element's degrees of freedom
*/
virtual void prepareLowerD()
{
mooseError("Lower dimensional FFT variables are not supported yet.");
}

virtual void prepareAux() {}

virtual void reinitNode() {}
virtual void reinitAux() {}
virtual void reinitAuxNeighbor() { mooseError("Neighbor FFT variables are not supported yet."); }

virtual void reinitNodes(const std::vector<dof_id_type> &)
{
mooseError("Nodal FFT variables are not supported.");
}
virtual void reinitNodesNeighbor(const std::vector<dof_id_type> & nodes)
{
mooseError("Nodal FFT variables are not supported.");
}

/**
* Field type of this variable
*/
virtual Moose::VarFieldType fieldType() const { return Moose::VarFieldType::VAR_FIELD_STANDARD; }

/**
* @returns true if this is a vector-valued element, false otherwise.
*/
virtual bool isVector() const { return false; };

/**
* Is this variable defined at nodes
* @return true if it the variable is defined at nodes, otherwise false
*/
virtual bool isNodalDefined() const { return false; }

virtual const dof_id_type & nodalDofIndex() const
{
mooseError("Nodal FFT variables are not supported.");
}

virtual const dof_id_type & nodalDofIndexNeighbor() const
{
mooseError("Nodal FFT variables are not supported.");
}

/**
* Current element this variable is evaluated at
*/
virtual const Elem * const & currentElem() const
{
mooseError("Current element access not supported.");
}

/**
* The subdomains the variable is active on
*/
virtual const std::set<SubdomainID> & activeSubdomains() const
{
mooseError("Active subdomain access not supported.");
}
/**
* Is the variable active on the subdomain?
* @param subdomain The subdomain id in question
* @return true if active on subdomain, false otherwise
*/
virtual bool activeOnSubdomain(SubdomainID subdomain) const { return true; }

/**
* Prepare the initial condition
*/
virtual void prepareIC() {}

/**
* Compute values at interior quadrature points
*/
virtual void computeElemValues()
{ // meat!
}

/**
* Compute values at facial quadrature points
*/
virtual void computeElemValuesFace() {}
/**
* Compute values at facial quadrature points for the neighbor
*/
virtual void computeNeighborValuesFace() {}
/**
* Compute values at quadrature points for the neighbor
*/
virtual void computeNeighborValues() {}
/**
* compute values at quadrature points on the lower dimensional element
*/
virtual void computeLowerDValues() {}
/**
* Compute nodal values of this variable in the neighbor
*/
virtual void computeNodalNeighborValues() {}
/**
* Compute nodal values of this variable
*/
virtual void computeNodalValues() {}

virtual void getDofIndices(const Elem * elem, std::vector<dof_id_type> & dof_indices) const {}
/**
* Get neighbor DOF indices for currently selected element
* @return the neighbor degree of freedom indices
*/
virtual const std::vector<dof_id_type> & dofIndicesNeighbor() const { return _no_dofs; }

/**
* Get dof indices for the current lower dimensional element (this is meaningful when performing
* mortar FEM)
* @return the lower dimensional element's dofs
*/
virtual const std::vector<dof_id_type> & dofIndicesLower() const { return _no_dofs; }

virtual unsigned int numberOfDofsNeighbor() { return 0; }

virtual void insert(NumericVector<Number> & residual) {}
virtual void add(NumericVector<Number> & residual) {}

///@{ FFT variables have no shape functions
virtual size_t phiSize() const { return 0; }
virtual size_t phiFaceSize() const { return 0; }
virtual size_t phiNeighborSize() const { return 0; }
virtual size_t phiFaceNeighborSize() const { return 0; }
virtual size_t phiLowerSize() const { return 0; }
///@}

protected:
std::vector<dof_id_type> _no_dofs;
};
21 changes: 21 additions & 0 deletions src/variables/MooseFFTVariable.C
Original file line number Diff line number Diff line change
@@ -0,0 +1,21 @@
/**********************************************************************/
/* DO NOT MODIFY THIS HEADER */
/* MAGPIE - Mesoscale Atomistic Glue Program for Integrated Execution */
/* */
/* Copyright 2017 Battelle Energy Alliance, LLC */
/* ALL RIGHTS RESERVED */
/**********************************************************************/

#include <MooseFFTVariable.h>

InputParameters
MooseFFTVariable::validParams()
{
InputParameters params = MooseVariableFEBase::validParams();
return params;
}

MooseFFTVariable::MooseFFTVariable(const InputParameters & parameters)
: MooseVariableFEBase(parameters)
{
}

0 comments on commit 077ca17

Please sign in to comment.