#include <matrixFreePDE.h>
Public Member Functions | |
MatrixFreePDE (userInputParameters< dim >) | |
~MatrixFreePDE () | |
virtual void | init () |
virtual void | makeTriangulation (parallel::distributed::Triangulation< dim > &) const |
void | initForTests () |
void | solve () |
void | vmult (vectorType &dst, const vectorType &src) const |
void | buildFields () |
virtual void | setInitialCondition (const dealii::Point< dim > &p, const unsigned int index, double &scalar_IC, dealii::Vector< double > &vector_IC)=0 |
virtual void | setNonUniformDirichletBCs (const dealii::Point< dim > &p, const unsigned int index, const unsigned int direction, const double time, double &scalar_BC, dealii::Vector< double > &vector_BC)=0 |
Public Attributes | |
std::vector< Field< dim > > | fields |
ConditionalOStream | pcout |
Protected Member Functions | |
void | reinit () |
void | reassignGrains () |
virtual void | solveIncrement (bool skip_time_dependent) |
void | outputResults () |
void | computeInvM () |
void | refineGrid () |
void | adaptiveRefine (unsigned int _currentIncrement) |
virtual void | adaptiveRefineCriterion () |
void | computeExplicitRHS () |
void | computeNonexplicitRHS () |
void | getLHS (const MatrixFree< dim, double > &data, vectorType &dst, const vectorType &src, const std::pair< unsigned int, unsigned int > &cell_range) const |
void | getLaplaceLHS (const MatrixFree< dim, double > &data, vectorType &dst, const vectorType &src, const std::pair< unsigned int, unsigned int > &cell_range) const |
void | setNonlinearEqInitialGuess () |
void | computeLaplaceRHS (unsigned int fieldIndex) |
void | getLaplaceRHS (const MatrixFree< dim, double > &data, vectorType &dst, const vectorType &src, const std::pair< unsigned int, unsigned int > &cell_range) const |
void | getExplicitRHS (const MatrixFree< dim, double > &data, std::vector< vectorType *> &dst, const std::vector< vectorType *> &src, const std::pair< unsigned int, unsigned int > &cell_range) const |
void | getNonexplicitRHS (const MatrixFree< dim, double > &data, std::vector< vectorType *> &dst, const std::vector< vectorType *> &src, const std::pair< unsigned int, unsigned int > &cell_range) const |
virtual void | explicitEquationRHS (variableContainer< dim, degree, dealii::VectorizedArray< double > > &variable_list, dealii::Point< dim, dealii::VectorizedArray< double > > q_point_loc) const =0 |
virtual void | nonExplicitEquationRHS (variableContainer< dim, degree, dealii::VectorizedArray< double > > &variable_list, dealii::Point< dim, dealii::VectorizedArray< double > > q_point_loc) const =0 |
virtual void | equationLHS (variableContainer< dim, degree, dealii::VectorizedArray< double > > &variable_list, dealii::Point< dim, dealii::VectorizedArray< double > > q_point_loc) const =0 |
virtual void | postProcessedFields (const variableContainer< dim, degree, dealii::VectorizedArray< double > > &variable_list, variableContainer< dim, degree, dealii::VectorizedArray< double > > &pp_variable_list, const dealii::Point< dim, dealii::VectorizedArray< double > > q_point_loc) const |
void | computePostProcessedFields (std::vector< vectorType *> &postProcessedSet) |
void | getPostProcessedFields (const dealii::MatrixFree< dim, double > &data, std::vector< vectorType *> &dst, const std::vector< vectorType *> &src, const std::pair< unsigned int, unsigned int > &cell_range) |
void | markBoundaries (parallel::distributed::Triangulation< dim > &) const |
void | applyDirichletBCs () |
void | applyNeumannBCs () |
void | setPeriodicity () |
void | setPeriodicityConstraints (ConstraintMatrix *, const DoFHandler< dim > *) const |
void | getComponentsWithRigidBodyModes (std::vector< int > &) const |
void | setRigidBodyModeConstraints (const std::vector< int >, ConstraintMatrix *, const DoFHandler< dim > *) const |
void | applyInitialConditions () |
void | save_checkpoint () |
void | load_checkpoint_triangulation () |
void | load_checkpoint_fields () |
void | load_checkpoint_time_info () |
void | move_file (const std::string &, const std::string &) |
void | verify_checkpoint_file_exists (const std::string filename) |
void | updateNucleiList () |
std::vector< nucleus< dim > > | getNewNuclei () |
void | getLocalNucleiList (std::vector< nucleus< dim > > &newnuclei) const |
void | safetyCheckNewNuclei (std::vector< nucleus< dim > > newnuclei, std::vector< unsigned int > &conflict_ids) |
void | refineMeshNearNuclei (std::vector< nucleus< dim > > newnuclei) |
double | weightedDistanceFromNucleusCenter (const dealii::Point< dim, double > center, const std::vector< double > semiaxes, const dealii::Point< dim, double > q_point_loc, const unsigned int var_index) const |
dealii::VectorizedArray< double > | weightedDistanceFromNucleusCenter (const dealii::Point< dim, double > center, const std::vector< double > semiaxes, const dealii::Point< dim, dealii::VectorizedArray< double > > q_point_loc, const unsigned int var_index) const |
virtual double | getNucleationProbability (variableValueContainer, double, dealii::Point< dim >, unsigned int variable_index) const |
unsigned int | getFieldIndex (std::string _name) |
void | outputFreeEnergy (const std::vector< double > &freeEnergyValues) const |
void | computeIntegral (double &integratedField, int index, std::vector< vectorType *> postProcessedSet) |
void | computeIntegralMF (double &integratedField, int index, const std::vector< vectorType *> postProcessedSet) |
void | getIntegralMF (const MatrixFree< dim, double > &data, std::vector< vectorType *> &dst, const std::vector< vectorType *> &src, const std::pair< unsigned int, unsigned int > &cell_range) |
Protected Attributes | |
userInputParameters< dim > | userInputs |
unsigned int | totalDOFs |
variableAttributeLoader | var_attributes |
std::vector< SimplifiedGrainRepresentation< dim > > | simplified_grain_representations |
parallel::distributed::Triangulation< dim > | triangulation |
std::vector< FESystem< dim > * > | FESet |
std::vector< const ConstraintMatrix * > | constraintsDirichletSet |
std::vector< const ConstraintMatrix * > | constraintsOtherSet |
std::vector< const DoFHandler< dim > * > | dofHandlersSet |
std::vector< const IndexSet * > | locally_relevant_dofsSet |
std::vector< ConstraintMatrix * > | constraintsDirichletSet_nonconst |
std::vector< ConstraintMatrix * > | constraintsOtherSet_nonconst |
std::vector< DoFHandler< dim > * > | dofHandlersSet_nonconst |
std::vector< IndexSet * > | locally_relevant_dofsSet_nonconst |
std::vector< vectorType * > | solutionSet |
std::vector< vectorType * > | residualSet |
std::vector< parallel::distributed::SolutionTransfer< dim, vectorType > * > | soltransSet |
DoFHandler< dim > * | vector_dofHandler |
FESystem< dim > * | vector_fe |
MatrixFree< dim, double > | vector_matrixFreeObject |
MatrixFree< dim, double > | matrixFreeObject |
vectorType | invM |
vectorType | dU_vector |
vectorType | dU_scalar |
unsigned int | currentFieldIndex |
bool | generatingInitialGuess |
std::vector< std::map< dealii::types::global_dof_index, double > * > | valuesDirichletSet |
std::vector< nucleus< dim > > | nuclei |
std::vector< double > | freeEnergyValues |
bool | isTimeDependentBVP |
bool | isEllipticBVP |
bool | hasExplicitEquation |
bool | hasNonExplicitEquation |
bool | has_Dirichlet_BCs |
unsigned int | parabolicFieldIndex |
unsigned int | ellipticFieldIndex |
double | currentTime |
unsigned int | currentIncrement |
unsigned int | currentOutput |
unsigned int | currentCheckpoint |
unsigned int | current_grain_reassignment |
TimerOutput | computing_timer |
std::vector< double > | integrated_postprocessed_fields |
bool | first_integrated_var_output_complete |
double | integrated_var |
unsigned int | integral_index |
dealii::Threads::Mutex | assembler_lock |
Static Protected Attributes | |
static const unsigned int | CIJ_tensor_size = 2*dim-1+dim/3 |
This is the abstract base class for the matrix free implementation of Parabolic and Elliptic BVP's, and supports MPI, Threads and Vectorization (Hybrid Parallel). This class contains the parallel data structures, mesh (referred to as triangulation), parallel degrees of freedom distribution, constraints, and general utility methods.
All the physical models in this package inherit this base class.
Definition at line 67 of file matrixFreePDE.h.
MatrixFreePDE | ( | userInputParameters< dim > | _userInputs | ) |
Class contructor
Definition at line 8 of file matrixFreePDE.cc.
~MatrixFreePDE | ( | ) |
Definition at line 33 of file matrixFreePDE.cc.
|
protected |
|
protected |
Method for applying Dirichlet boundary conditions.
Definition at line 75 of file boundaryConditions.cc.
|
protected |
Definition at line 38 of file initialConditions.cc.
|
protected |
Method for applying Neumann boundary conditions.
Definition at line 12 of file boundaryConditions.cc.
void buildFields | ( | ) |
Definition at line 16 of file buildFields.cc.
|
protected |
Definition at line 8 of file computeRHS.cc.
|
protected |
Definition at line 5 of file computeIntegral.cc.
|
protected |
Definition at line 45 of file computeIntegral.cc.
|
protected |
Definition at line 93 of file setNonlinearEqInitialGuess.cc.
|
protected |
Definition at line 58 of file computeRHS.cc.
|
protected |
Definition at line 5 of file postprocessor.cc.
|
protectedpure virtual |
Implemented in testInvM< dim, degree >, and setRigidBodyModeConstraintsTest< dim, degree >.
|
protectedpure virtual |
Implemented in testInvM< dim, degree >, and setRigidBodyModeConstraintsTest< dim, degree >.
|
protected |
Definition at line 197 of file boundaryConditions.cc.
|
protected |
Definition at line 27 of file computeRHS.cc.
|
protected |
Definition at line 7 of file utilities.cc.
|
protected |
Definition at line 62 of file computeIntegral.cc.
|
protected |
Definition at line 131 of file setNonlinearEqInitialGuess.cc.
|
protected |
Definition at line 109 of file setNonlinearEqInitialGuess.cc.
|
protected |
Definition at line 37 of file computeLHS.cc.
|
protected |
Definition at line 93 of file nucleation.cc.
|
protected |
Definition at line 65 of file nucleation.cc.
|
protected |
Definition at line 77 of file computeRHS.cc.
|
inlineprotectedvirtual |
Definition at line 329 of file matrixFreePDE.h.
|
protected |
Definition at line 27 of file postprocessor.cc.
|
virtual |
void initForTests | ( | ) |
Initializes the data structures for enabling unit tests.
This method initializes the MatrixFreePDE object with a fixed geometry, discretization and other custom selected options specifically to help with unit tests, and should not be called in any of the physical models.
Definition at line 7 of file initForTests.cc.
|
protected |
Definition at line 124 of file checkpoint.cc.
|
protected |
Definition at line 165 of file checkpoint.cc.
|
protected |
Definition at line 103 of file checkpoint.cc.
|
virtual |
|
protected |
Definition at line 9 of file markBoundaries.cc.
|
protected |
Definition at line 187 of file checkpoint.cc.
|
protectedpure virtual |
Implemented in testInvM< dim, degree >, and setRigidBodyModeConstraintsTest< dim, degree >.
|
protected |
|
protected |
Definition at line 8 of file outputResults.cc.
|
inlineprotectedvirtual |
Definition at line 266 of file matrixFreePDE.h.
|
protected |
Method to reassign grains when multiple grains are stored in a single order parameter.
Definition at line 10 of file reassignGrains.cc.
|
protected |
Definition at line 309 of file nucleation.cc.
|
protected |
Definition at line 255 of file nucleation.cc.
|
protected |
Definition at line 11 of file checkpoint.cc.
|
pure virtual |
Implemented in testInvM< dim, degree >, and setRigidBodyModeConstraintsTest< dim, degree >.
|
protected |
Definition at line 8 of file setNonlinearEqInitialGuess.cc.
|
pure virtual |
Implemented in testInvM< dim, degree >, and setRigidBodyModeConstraintsTest< dim, degree >.
|
protected |
Definition at line 151 of file boundaryConditions.cc.
|
protected |
Definition at line 172 of file boundaryConditions.cc.
|
protected |
Definition at line 238 of file boundaryConditions.cc.
void solve | ( | ) |
This method implements the time stepping algorithm and invokes the solveIncrement() method.
|
protectedvirtual |
Method to solve each time increment of a time-dependent problem. For time-independent problems this method is called only once. This method solves for all the fields in a staggered manner (one after another) and also invokes the corresponding solvers: Explicit solver for Parabolic problems, Implicit (matrix-free) solver for Elliptic problems.
Definition at line 8 of file solveIncrement.cc.
|
protected |
Definition at line 12 of file nucleation.cc.
|
protected |
Definition at line 216 of file checkpoint.cc.
void vmult | ( | vectorType & | dst, |
const vectorType & | src | ||
) | const |
This method essentially converts the MatrixFreePDE object into a matrix object which can be used with matrix free iterative solvers. Provides the A*x functionality for solving the system of equations AX=b.
Definition at line 7 of file computeLHS.cc.
|
protected |
Definition at line 375 of file nucleation.cc.
|
protected |
Definition at line 399 of file nucleation.cc.
|
protected |
Definition at line 365 of file matrixFreePDE.h.
|
staticprotected |
Definition at line 131 of file matrixFreePDE.h.
|
mutableprotected |
Definition at line 356 of file matrixFreePDE.h.
|
protected |
Definition at line 168 of file matrixFreePDE.h.
|
protected |
Definition at line 178 of file matrixFreePDE.h.
|
protected |
Definition at line 168 of file matrixFreePDE.h.
|
protected |
Definition at line 178 of file matrixFreePDE.h.
|
protected |
Definition at line 353 of file matrixFreePDE.h.
|
protected |
Definition at line 353 of file matrixFreePDE.h.
|
protected |
Definition at line 207 of file matrixFreePDE.h.
|
protected |
Definition at line 353 of file matrixFreePDE.h.
|
protected |
Definition at line 353 of file matrixFreePDE.h.
|
protected |
Definition at line 352 of file matrixFreePDE.h.
|
protected |
Definition at line 171 of file matrixFreePDE.h.
|
protected |
Definition at line 180 of file matrixFreePDE.h.
|
protected |
Definition at line 203 of file matrixFreePDE.h.
|
protected |
Definition at line 203 of file matrixFreePDE.h.
|
protected |
Definition at line 351 of file matrixFreePDE.h.
|
protected |
Definition at line 163 of file matrixFreePDE.h.
std::vector<Field<dim> > fields |
Vector of all the physical fields in the problem. Fields are identified by dimentionality (SCALAR/VECTOR), the kind of PDE (ELLIPTIC/PARABOLIC) used to compute them and a character identifier (e.g.: "c" for composition) which is used to write the fields to the output files.
Definition at line 107 of file matrixFreePDE.h.
|
protected |
Definition at line 360 of file matrixFreePDE.h.
|
protected |
Definition at line 335 of file matrixFreePDE.h.
|
protected |
Definition at line 231 of file matrixFreePDE.h.
|
protected |
Definition at line 349 of file matrixFreePDE.h.
|
protected |
Definition at line 347 of file matrixFreePDE.h.
|
protected |
Definition at line 348 of file matrixFreePDE.h.
|
protected |
Definition at line 364 of file matrixFreePDE.h.
|
protected |
Definition at line 358 of file matrixFreePDE.h.
|
protected |
Definition at line 363 of file matrixFreePDE.h.
|
protected |
Definition at line 201 of file matrixFreePDE.h.
|
protected |
Definition at line 345 of file matrixFreePDE.h.
|
protected |
Definition at line 343 of file matrixFreePDE.h.
|
protected |
Definition at line 176 of file matrixFreePDE.h.
|
protected |
Definition at line 182 of file matrixFreePDE.h.
|
protected |
Definition at line 199 of file matrixFreePDE.h.
|
protected |
Definition at line 316 of file matrixFreePDE.h.
|
protected |
Definition at line 351 of file matrixFreePDE.h.
ConditionalOStream pcout |
Definition at line 112 of file matrixFreePDE.h.
|
protected |
Definition at line 186 of file matrixFreePDE.h.
|
protected |
Definition at line 141 of file matrixFreePDE.h.
|
protected |
Definition at line 188 of file matrixFreePDE.h.
|
protected |
Definition at line 184 of file matrixFreePDE.h.
|
protected |
Definition at line 123 of file matrixFreePDE.h.
|
protected |
Definition at line 159 of file matrixFreePDE.h.
|
protected |
Definition at line 121 of file matrixFreePDE.h.
|
protected |
Definition at line 278 of file matrixFreePDE.h.
|
protected |
Definition at line 128 of file matrixFreePDE.h.
|
protected |
Definition at line 191 of file matrixFreePDE.h.
|
protected |
Definition at line 192 of file matrixFreePDE.h.
|
protected |
Definition at line 193 of file matrixFreePDE.h.