Skip to content

Commit

Permalink
Must compute time derivatives ahead of Jacobian for coloring
Browse files Browse the repository at this point in the history
With the change to not compute a pre-SMO residual, when we are
computing a finite difference Jacobian via coloring, the initial
Jacobian computation happens before any residual evaluations. In
this case we must ensure that we compute time derivatives on the
`NONLINEAR` exec flag

Refs idaholab#23472
  • Loading branch information
lindsayad committed Apr 7, 2024
1 parent 4781fca commit c5af664
Show file tree
Hide file tree
Showing 4 changed files with 12 additions and 3 deletions.
2 changes: 2 additions & 0 deletions framework/include/systems/NonlinearSystem.h
Original file line number Diff line number Diff line change
Expand Up @@ -100,5 +100,7 @@ class NonlinearSystem : public NonlinearSystemBase
*/
void setupColoringFiniteDifferencedPreconditioner();

virtual bool matrixFromColoring() const override { return _use_coloring_finite_difference; }

bool _use_coloring_finite_difference;
};
7 changes: 7 additions & 0 deletions framework/include/systems/SolverSystem.h
Original file line number Diff line number Diff line change
Expand Up @@ -85,6 +85,13 @@ class SolverSystem : public SystemBase

/// Boolean to see if solution is invalid
bool _solution_is_invalid;

private:
/**
* Whether a system matrix is formed from coloring. This influences things like when to compute
* time derivatives
*/
virtual bool matrixFromColoring() const { return false; }
};

inline const NumericVector<Number> * const &
Expand Down
4 changes: 2 additions & 2 deletions framework/src/systems/NonlinearSystem.C
Original file line number Diff line number Diff line change
Expand Up @@ -232,14 +232,14 @@ NonlinearSystem::setupFiniteDifferencedPreconditioner()

if (fdp->finiteDifferenceType() == "coloring")
{
setupColoringFiniteDifferencedPreconditioner();
_use_coloring_finite_difference = true;
setupColoringFiniteDifferencedPreconditioner();
}

else if (fdp->finiteDifferenceType() == "standard")
{
setupStandardFiniteDifferencedPreconditioner();
_use_coloring_finite_difference = false;
setupStandardFiniteDifferencedPreconditioner();
}
else
mooseError("Unknown finite difference type");
Expand Down
2 changes: 1 addition & 1 deletion framework/src/systems/SolverSystem.C
Original file line number Diff line number Diff line change
Expand Up @@ -135,7 +135,7 @@ SolverSystem::compute(const ExecFlagType type)
compute_tds = true;
else if (type == EXEC_NONLINEAR)
{
if (_fe_problem.computingScalingJacobian())
if (_fe_problem.computingScalingJacobian() || matrixFromColoring())
compute_tds = true;
}
else if ((type == EXEC_TIMESTEP_END) || (type == EXEC_FINAL))
Expand Down

0 comments on commit c5af664

Please sign in to comment.