Skip to content

Commit

Permalink
Create infrastructure in FEProblem and Compute*JacobianThread as well…
Browse files Browse the repository at this point in the history
… as NonlocalIntegratedBC (issue idaholab#7980).
  • Loading branch information
Alex Lindsay committed Nov 3, 2016
1 parent fe584dd commit 15bc94a
Show file tree
Hide file tree
Showing 15 changed files with 851 additions and 1 deletion.
4 changes: 4 additions & 0 deletions framework/include/base/FEProblem.h
Original file line number Diff line number Diff line change
Expand Up @@ -67,6 +67,7 @@ class InternalSideUserObject;
class GeneralUserObject;
class Function;
class KernelBase;
class IntegratedBC;

// libMesh forward declarations
namespace libMesh
Expand Down Expand Up @@ -1092,6 +1093,9 @@ class FEProblem :
/// nonlocal kernels
MooseObjectWarehouse<KernelBase> _nonlocal_kernels;

/// nonlocal integrated_bcs
MooseObjectWarehouse<IntegratedBC> _nonlocal_integrated_bcs;

///@{
/// Initial condition storage
InitialConditionWarehouse _ics;
Expand Down
12 changes: 12 additions & 0 deletions framework/include/bcs/IntegratedBC.h
Original file line number Diff line number Diff line change
Expand Up @@ -53,6 +53,18 @@ class IntegratedBC :
*/
void computeJacobianBlockScalar(unsigned int jvar);

/**
* Compute this IntegratedBC's contribution to the diagonal Jacobian entries
* corresponding to nonlocal dofs of the variable
*/
virtual void computeNonlocalJacobian() {}

/**
* Computes d-residual / d-jvar... corresponding to nonlocal dofs of the jvar
* and stores the result in nonlocal ke
*/
virtual void computeNonlocalOffDiagJacobian(unsigned int /* jvar */) {}

protected:
/// current element
const Elem * & _current_elem;
Expand Down
71 changes: 71 additions & 0 deletions framework/include/bcs/NonlocalIntegratedBC.h
Original file line number Diff line number Diff line change
@@ -0,0 +1,71 @@
/****************************************************************/
/* 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 NONLOCALINTEGRATEDBC_H
#define NONLOCALINTEGRATEDBC_H

#include "IntegratedBC.h"

class NonlocalIntegratedBC;

template<>
InputParameters validParams<NonlocalIntegratedBC>();

/**
* NonlocalIntegratedBC is used for solving integral terms in integro-differential equations.
* Integro-differential equations includes spatial integral terms over variables in the domain.
* In this case the jacobian calculation is not restricted to local dofs of an element, it requires
* additional contributions from all the dofs in the domain. NonlocalIntegratedBC adds capability to consider
* nonlocal jacobians in addition to the local jacobians.
*/
class NonlocalIntegratedBC : public IntegratedBC
{
public:
NonlocalIntegratedBC(const InputParameters & parameters);

/**
* computeJacobian and computeQpOffDiagJacobian methods are almost same
* as IntegratedBC except for few additional optimization options regarding the integral terms.
* Looping order of _i & _j are reversed for providing opimization options and
* additional getUserObjectJacobian method is provided as an option to obtain
* jocobians of the integral term from userobject once per dof
*/
virtual void computeJacobian();
virtual void computeOffDiagJacobian(unsigned int jvar);

/**
* computeNonlocalJacobian and computeNonlocalOffDiagJacobian methods are
* introduced for providing the jacobian contribution corresponding to nolocal dofs.
* Nonlocal jacobian block is sized with respect to the all the dofs of jacobian variable.
* Looping order of _i & _k are reversed for opimization purpose and
* additional globalDoFEnabled method is provided as an option to execute nonlocal
* jocobian calculations only for nonlocal dofs that has nonzero jacobian contribution.
*/
virtual void computeNonlocalJacobian();
virtual void computeNonlocalOffDiagJacobian(unsigned int jvar);

protected:
/// Compute this IntegratedBC's contribution to the Jacobian corresponding to nolocal dof at the current quadrature point
virtual Real computeQpNonlocalJacobian(dof_id_type /*dof_index*/) { return 0; }
virtual Real computeQpNonlocalOffDiagJacobian(unsigned int /*jvar*/, dof_id_type /*dof_index*/) { return 0; }

/// Optimization option for getting jocobinas from userobject once per dof
virtual void getUserObjectJacobian(unsigned int /*jvar*/, dof_id_type /*dof_index*/) {}
/// optimization option for executing nonlocal jacobian calculation only for nonzero elements
virtual bool globalDoFEnabled(MooseVariable & /*var*/, dof_id_type /*dof_index*/) { return true; }

unsigned int _k;
};

#endif /* NONLOCALINTEGRATEDBC_H */
30 changes: 30 additions & 0 deletions framework/src/base/ComputeFullJacobianThread.C
Original file line number Diff line number Diff line change
Expand Up @@ -20,6 +20,7 @@
#include "DGKernel.h"
#include "InterfaceKernel.h"
#include "NonlocalKernel.h"
#include "NonlocalIntegratedBC.h"
// libmesh includes
#include "libmesh/threads.h"

Expand Down Expand Up @@ -148,6 +149,35 @@ ComputeFullJacobianThread::computeFaceJacobian(BoundaryID bnd_id)
}
}

/// done only when nonlocal integrated_bcs exist in the system
if (_fe_problem.checkNonlocalCouplingRequirement())
{
std::vector<std::pair<MooseVariable *, MooseVariable *> > & cne = _fe_problem.nonlocalCouplingEntries(_tid);
for (const auto & it : cne)
{
MooseVariable & ivariable = *(it.first);
MooseVariable & jvariable = *(it.second);

unsigned int ivar = ivariable.number();
unsigned int jvar = jvariable.number();

if (ivariable.activeOnSubdomain(_subdomain) && jvariable.activeOnSubdomain(_subdomain) && _integrated_bcs.hasActiveBoundaryObjects(bnd_id, _tid))
{
const std::vector<MooseSharedPointer<IntegratedBC> > & integrated_bcs = _integrated_bcs.getBoundaryObjects(bnd_id, _tid);
for (const auto & integrated_bc : integrated_bcs)
{
MooseSharedPointer<NonlocalIntegratedBC> nonlocal_integrated_bc = MooseSharedNamespace::dynamic_pointer_cast<NonlocalIntegratedBC>(integrated_bc);
if (nonlocal_integrated_bc)
if ((integrated_bc->variable().number() == ivar) && integrated_bc->isImplicit())
{
integrated_bc->subProblem().prepareShapes(jvar, _tid);
integrated_bc->computeNonlocalOffDiagJacobian(jvar);
}
}
}
}
}

const std::vector<MooseVariableScalar *> & scalar_vars = _nl.getScalarVariables(_tid);
if (scalar_vars.size() > 0)
{
Expand Down
8 changes: 8 additions & 0 deletions framework/src/base/ComputeJacobianThread.C
Original file line number Diff line number Diff line change
Expand Up @@ -21,6 +21,7 @@
#include "InterfaceKernel.h"
#include "KernelWarehouse.h"
#include "NonlocalKernel.h"
#include "NonlocalIntegratedBC.h"

// libmesh includes
#include "libmesh/threads.h"
Expand Down Expand Up @@ -85,6 +86,13 @@ ComputeJacobianThread::computeFaceJacobian(BoundaryID bnd_id)
{
bc->subProblem().prepareFaceShapes(bc->variable().number(), _tid);
bc->computeJacobian();
/// done only when nonlocal integrated_bcs exist in the system
if (_fe_problem.checkNonlocalCouplingRequirement())
{
MooseSharedPointer<NonlocalIntegratedBC> nonlocal_integrated_bc = MooseSharedNamespace::dynamic_pointer_cast<NonlocalIntegratedBC>(bc);
if (nonlocal_integrated_bc)
bc->computeNonlocalJacobian();
}
}
}

Expand Down
32 changes: 31 additions & 1 deletion framework/src/base/FEProblem.C
Original file line number Diff line number Diff line change
Expand Up @@ -74,6 +74,7 @@
#include "XFEMInterface.h"
#include "ConsoleUtils.h"
#include "NonlocalKernel.h"
#include "NonlocalIntegratedBC.h"
#include "ShapeElementUserObject.h"
#include "ShapeSideUserObject.h"

Expand Down Expand Up @@ -701,7 +702,7 @@ void FEProblem::timestepSetup()
_app.getOutputWarehouse().timestepSetup();

if (_requires_nonlocal_coupling)
if (_nonlocal_kernels.hasActiveObjects())
if (_nonlocal_kernels.hasActiveObjects() || _nonlocal_integrated_bcs.hasActiveObjects())
_has_nonlocal_coupling = true;
}

Expand Down Expand Up @@ -744,6 +745,18 @@ FEProblem::checkNonlocalCoupling()
_nonlocal_kernels.addObject(kernel, tid);
}
}
const MooseObjectWarehouse<IntegratedBC> & all_integrated_bcs = _nl.getIntegratedBCWarehouse();
const std::vector<MooseSharedPointer<IntegratedBC> > & integrated_bcs = all_integrated_bcs.getObjects(tid);
for (const auto & integrated_bc : integrated_bcs)
{
MooseSharedPointer<NonlocalIntegratedBC> nonlocal_integrated_bc = MooseSharedNamespace::dynamic_pointer_cast<NonlocalIntegratedBC>(integrated_bc);
if (nonlocal_integrated_bc)
{
if (_calculate_jacobian_in_uo)
_requires_nonlocal_coupling = true;
_nonlocal_integrated_bcs.addObject(integrated_bc, tid);
}
}
}
}

Expand Down Expand Up @@ -3122,7 +3135,9 @@ FEProblem::setNonlocalCouplingMatrix()
_nonlocal_cm.resize(n_vars);
const std::vector<MooseVariable *> & vars = _nl.getVariables(0);
const std::vector<MooseSharedPointer<KernelBase> > & nonlocal_kernel = _nonlocal_kernels.getObjects();
const std::vector<MooseSharedPointer<IntegratedBC> > & nonlocal_integrated_bc = _nonlocal_integrated_bcs.getObjects();
for (const auto & ivar : vars)
{
for (const auto & kernel : nonlocal_kernel)
{
unsigned int i = ivar->number();
Expand All @@ -3137,6 +3152,21 @@ FEProblem::setNonlocalCouplingMatrix()
}
}
}
for (const auto & integrated_bc : nonlocal_integrated_bc)
{
unsigned int i = ivar->number();
if (i == integrated_bc->variable().number())
for (const auto & jvar : vars)
{
const auto it = _var_dof_map.find(jvar->name());
if (it != _var_dof_map.end())
{
unsigned int j = jvar->number();
_nonlocal_cm(i,j) = 1;
}
}
}
}
}

bool
Expand Down
151 changes: 151 additions & 0 deletions framework/src/bcs/NonlocalIntegratedBC.C
Original file line number Diff line number Diff line change
@@ -0,0 +1,151 @@
/****************************************************************/
/* 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 */
/****************************************************************/

#include "NonlocalIntegratedBC.h"
#include "Assembly.h"
#include "MooseVariable.h"
#include "Problem.h"
#include "SubProblem.h"
#include "SystemBase.h"
#include "MooseMesh.h"

// libmesh includes
#include "libmesh/threads.h"
#include "libmesh/quadrature.h"

template<>
InputParameters validParams<NonlocalIntegratedBC>()
{
InputParameters params = validParams<IntegratedBC>();
return params;
}

NonlocalIntegratedBC::NonlocalIntegratedBC(const InputParameters & parameters) :
IntegratedBC(parameters)
{
_mesh.errorIfDistributedMesh("NonlocalIntegratedBC");
mooseWarning("NonlocalIntegratedBC is a computationally expensive experimental capability used only for integral terms.");
}

void
NonlocalIntegratedBC::computeJacobian()
{
DenseMatrix<Number> & ke = _assembly.jacobianBlock(_var.number(), _var.number());
_local_ke.resize(ke.m(), ke.n());
_local_ke.zero();

for (_j = 0; _j < _phi.size(); _j++) // looping order for _i & _j are reversed for performance improvement
{
getUserObjectJacobian(_var.number(), _var.dofIndices()[_j]);
for (_i = 0; _i < _test.size(); _i++)
for (_qp = 0; _qp < _qrule->n_points(); _qp++)
_local_ke(_i, _j) += _JxW[_qp] * _coord[_qp] * computeQpJacobian();
}

ke += _local_ke;

if (_has_diag_save_in)
{
unsigned int rows = ke.m();
DenseVector<Number> diag(rows);
for (unsigned int i=0; i<rows; i++)
diag(i) = _local_ke(i,i);

Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
for (const auto & var : _diag_save_in)
var->sys().solution().add_vector(diag, var->dofIndices());
}
}

void
NonlocalIntegratedBC::computeOffDiagJacobian(unsigned int jvar)
{
if (jvar == _var.number())
computeJacobian();
else
{
MooseVariable & jv = _sys.getVariable(_tid, jvar);
DenseMatrix<Number> & ke = _assembly.jacobianBlock(_var.number(), jvar);

for (_j = 0; _j < _phi.size(); _j++) // looping order for _i & _j are reversed for performance improvement
{
getUserObjectJacobian(jvar, jv.dofIndices()[_j]);
for (_i = 0; _i < _test.size(); _i++)
for (_qp = 0; _qp < _qrule->n_points(); _qp++)
ke(_i, _j) += _JxW[_qp] * _coord[_qp] * computeQpOffDiagJacobian(jvar);
}
}
}

void
NonlocalIntegratedBC::computeNonlocalJacobian()
{
DenseMatrix<Number> & keg = _assembly.jacobianBlockNonlocal(_var.number(), _var.number());
// compiling set of global IDs for the local DOFs on the element
std::set<dof_id_type> local_dofindices(_var.dofIndices().begin(), _var.dofIndices().end());
// storing the global IDs for all the DOFs of the variable
const std::vector<dof_id_type> & var_alldofindices = _var.allDofIndices();
unsigned int n_total_dofs = var_alldofindices.size();

for (_k = 0; _k < n_total_dofs; _k++) // looping order for _i & _k are reversed for performance improvement
{
// eliminating the local components
auto it = local_dofindices.find(var_alldofindices[_k]);
if (it == local_dofindices.end())
{
getUserObjectJacobian(_var.number(), var_alldofindices[_k]);
// skip global DOFs that do not contribute to the jacobian
if (!globalDoFEnabled(_var, var_alldofindices[_k]))
continue;

for (_i = 0; _i < _test.size(); _i++)
for (_qp = 0; _qp < _qrule->n_points(); _qp++)
keg(_i, _k) += _JxW[_qp] * _coord[_qp] * computeQpNonlocalJacobian(var_alldofindices[_k]);
}
}
}

void
NonlocalIntegratedBC::computeNonlocalOffDiagJacobian(unsigned int jvar)
{
if (jvar == _var.number())
computeNonlocalJacobian();
else
{
MooseVariable & jv = _sys.getVariable(_tid, jvar);
DenseMatrix<Number> & keg = _assembly.jacobianBlockNonlocal(_var.number(), jvar);
// compiling set of global IDs for the local DOFs on the element
std::set<dof_id_type> local_dofindices(jv.dofIndices().begin(), jv.dofIndices().end());
// storing the global IDs for all the DOFs of the variable
const std::vector<dof_id_type> & jv_alldofindices = jv.allDofIndices();
unsigned int n_total_dofs = jv_alldofindices.size();

for (_k = 0; _k < n_total_dofs; _k++) // looping order for _i & _k are reversed for performance improvement
{
// eliminating the local components
auto it = local_dofindices.find(jv_alldofindices[_k]);
if (it == local_dofindices.end())
{
getUserObjectJacobian(jvar, jv_alldofindices[_k]);
// skip global DOFs that do not contribute to the jacobian
if (!globalDoFEnabled(jv, jv_alldofindices[_k]))
continue;

for (_i = 0; _i < _test.size(); _i++)
for (_qp = 0; _qp < _qrule->n_points(); _qp++)
keg(_i, _k) += _JxW[_qp] * _coord[_qp] * computeQpNonlocalOffDiagJacobian(jvar, jv_alldofindices[_k]);
}
}
}
}

0 comments on commit 15bc94a

Please sign in to comment.