Assemble and solve the equation system.
132{
135
138
139 std::vector<GlobalVector*> x_new{x};
140 x_new[process_id] =
143
144 bool error_norms_met = false;
145
147
148 int iteration = 1;
150 {
152 double time_dirichlet = 0.0;
153
155 time_iteration.
start();
156
157 timer_dirichlet.
start();
158 auto& x_new_process = *x_new[process_id];
160 sys.computeKnownSolutions(x_new_process, process_id);
161 sys.applyKnownSolutions(x_new_process);
162 time_dirichlet += timer_dirichlet.
elapsed();
163
164 sys.preIteration(iteration, x_new_process);
165
167 time_assembly.
start();
168 sys.assemble(x_new, x_prev, process_id);
169 sys.getA(A);
170 sys.getRhs(*x_prev[process_id], rhs);
171
172
173 if (sys.requiresNormalization() &&
175 {
176 sys.getAandRhsNormalized(A, rhs);
178 "The equation system is rectangular, but the current linear "
179 "solver only supports square systems. "
180 "The system will be normalized, which lead to a squared "
181 "condition number and potential numerical issues. "
182 "It is recommended to use a solver that supports rectangular "
183 "equation systems for better numerical stability.");
184 }
185
186 INFO(
"[time] Assembly took {:g} s.", time_assembly.
elapsed());
187
188
190 {
192 }
193
194 timer_dirichlet.
start();
195 sys.applyKnownSolutionsPicard(A, rhs, x_new_process);
196 time_dirichlet += timer_dirichlet.
elapsed();
197 INFO(
"[time] Applying Dirichlet BCs took {:g} s.", time_dirichlet);
198
200 {
205 }
206
207 bool iteration_succeeded =
209 sys.linearSolverNeedsToCompute());
210
211 if (iteration_succeeded)
212 {
213 if (postIterationCallback)
214 {
215 postIterationCallback(iteration, x_new);
216 }
217
218 switch (sys.postIteration(x_new_process))
219 {
221
222
223 break;
225 ERR(
"Picard: The postIteration() hook reported a "
226 "non-recoverable error.");
227 iteration_succeeded = false;
228
229
230
232 break;
235 "Picard: The postIteration() hook decided that this "
236 "iteration has to be repeated.");
238 *x[process_id],
239 x_new_process);
240 continue;
241 }
242 }
243
244 if (!iteration_succeeded)
245 {
246
247 error_norms_met = false;
248 break;
249 }
250
251 if (sys.isLinear())
252 {
253 error_norms_met = true;
254 }
255 else
256 {
258 {
261 x_new_process);
263 x_new_process);
264 }
265
267 }
268
269
271
272 INFO(
"[time] Iteration #{:d} took {:g} s.", iteration,
274
275 if (error_norms_met)
276 {
277 break;
278 }
279
280
281
283 {
284 break;
285 }
286 }
287
289 {
290 ERR(
"Picard: Could not solve the given nonlinear system within {:d} "
291 "iterations",
293 }
294
298
299 return {error_norms_met, iteration};
300}
void ERR(fmt::format_string< Args... > fmt, Args &&... args)
void WARN(fmt::format_string< Args... > fmt, Args &&... args)
double elapsed() const
Get the elapsed time in seconds.
void start()
Start the timer.
bool canSolveRectangular() const
Get, if the solver can handle rectangular equation systems.
Global vector based on Eigen vector.
virtual void preFirstIteration()
virtual void checkResidual(GlobalVector const &residual)=0
Check if the residual satisfies the convergence criterion.
virtual bool hasResidualCheck() const =0
virtual bool isSatisfied() const
Tell if the convergence criterion is satisfied.
virtual void checkDeltaX(GlobalVector const &minus_delta_x, GlobalVector const &x)=0
virtual bool hasDeltaXCheck() const =0
void copy(PETScVector const &x, PETScVector &y)
void setLocalAccessibleVector(PETScVector const &x)
bool solvePicard(GlobalLinearSolver &linear_solver, GlobalMatrix &A, GlobalVector &rhs, GlobalVector &x, MathLib::LinearSolverBehaviour const linear_solver_behaviour)