forked from idaholab/moose
/
NonlinearSystem.h
122 lines (97 loc) · 3.75 KB
/
NonlinearSystem.h
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
//* This file is part of the MOOSE framework
//* https://www.mooseframework.org
//*
//* All rights reserved, see COPYRIGHT for full restrictions
//* https://github.com/idaholab/moose/blob/master/COPYRIGHT
//*
//* Licensed under LGPL 2.1, please see LICENSE for details
//* https://www.gnu.org/licenses/lgpl-2.1.html
#pragma once
#include "NonlinearSystemBase.h"
#include "ComputeResidualFunctor.h"
#include "ComputeFDResidualFunctor.h"
#include "SubProblem.h"
#include "MooseError.h"
#include "libmesh/transient_system.h"
#include "libmesh/nonlinear_implicit_system.h"
/**
* Nonlinear system to be solved
*
* It is a part of FEProblemBase ;-)
*/
class NonlinearSystem : public NonlinearSystemBase
{
public:
NonlinearSystem(FEProblemBase & problem, const std::string & name);
virtual ~NonlinearSystem();
virtual SparseMatrix<Number> & addMatrix(TagID tag) override;
virtual void solve() override;
void init() override;
/**
* Quit the current solve as soon as possible.
*/
virtual void stopSolve() override;
/**
* Returns the current nonlinear iteration number. In libmesh, this is
* updated during the nonlinear solve, so it should be up-to-date.
*/
virtual unsigned int getCurrentNonlinearIterationNumber() override
{
return _transient_sys.get_current_nonlinear_iteration_number();
}
virtual void setupFiniteDifferencedPreconditioner() override;
/**
* Returns the convergence state
* @return true if converged, otherwise false
*/
virtual bool converged() override;
virtual NumericVector<Number> & RHS() override { return *_transient_sys.rhs; }
virtual NonlinearSolver<Number> * nonlinearSolver() override
{
return _transient_sys.nonlinear_solver.get();
}
virtual TransientNonlinearImplicitSystem & sys() { return _transient_sys; }
void computeScaling() override;
virtual void attachPreconditioner(Preconditioner<Number> * preconditioner) override;
protected:
NumericVector<Number> & solutionOldInternal() const override
{
return *_transient_sys.old_local_solution;
}
NumericVector<Number> & solutionOlderInternal() const override
{
return *_transient_sys.older_local_solution;
}
TransientNonlinearImplicitSystem & _transient_sys;
ComputeResidualFunctor _nl_residual_functor;
ComputeFDResidualFunctor _fd_residual_functor;
private:
/**
* Setup group scaling containers
*/
void setupScalingGrouping();
/**
* Form preconditioning matrix via a standard finite difference method
* column-by-column. This method computes both diagonal and off-diagonal
* entrices regardless of the structure pattern of the Jacobian matrix.
*/
void setupStandardFiniteDifferencedPreconditioner();
/**
* According to the nonzero pattern provided in the matrix, a graph is constructed.
* A coloring algorithm is applied to the graph. The graph is partitioned into several
* independent subgraphs (colors), and a finte difference method is applied color-by-color
* to form a preconditioning matrix. If the number of colors is small, this method is much
* faster than the standard FD. But there is an issue. If the matrix provided by users does not
* represent the actual structure of the true Jacobian, the matrix computed via coloring could
* be wrong or inaccurate. In this case, users should switch to the standard finite difference
* method.
*/
void setupColoringFiniteDifferencedPreconditioner();
bool _use_coloring_finite_difference;
/// Whether we've initialized the automatic scaling data structures
bool _auto_scaling_initd;
/// A map from variable index to group variable index and it's associated (inverse) scaling factor
std::unordered_map<unsigned int, unsigned int> _var_to_group_var;
/// The number of scaling groups
std::size_t _num_scaling_groups;
};