From: Bruno Turcksin Date: Mon, 30 Aug 2021 18:58:46 +0000 (+0000) Subject: Fix a bug in EmbeddedExplicitRungeKutta default constructor X-Git-Tag: v9.4.0-rc1~1015^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=1c4698521f433107315733c2e85f2c5deef9ed49;p=dealii.git Fix a bug in EmbeddedExplicitRungeKutta default constructor - Fix a bug where some variables were not initialized. - Fix documentation in RungeKutta - Remove unused variable in ImplicitRungeKutta --- diff --git a/doc/news/changes/minor/20210830Turcksin b/doc/news/changes/minor/20210830Turcksin new file mode 100644 index 0000000000..0980c837a7 --- /dev/null +++ b/doc/news/changes/minor/20210830Turcksin @@ -0,0 +1,6 @@ +Fixed: There was a bug in the default constructor of +TimeStepping::EmbeddedExplicitRungeKutta(). Some parameters were left +uninitialized which triggered a segmentation fault when calling +`evolve_one_time_step()`. +
+(Praveen Chandrashekar, Bruno Turcksin, 2021/08/30) diff --git a/include/deal.II/base/time_stepping.h b/include/deal.II/base/time_stepping.h index 0699e1a543..d3a5d9fd23 100644 --- a/include/deal.II/base/time_stepping.h +++ b/include/deal.II/base/time_stepping.h @@ -261,7 +261,7 @@ namespace TimeStepping * the output is value of f at this point. @p id_minus_tau_J_inverse is a * function that computes $ inv(I-\tau J)$ where $ I $ is the identity * matrix, $ \tau $ is given, and $ J $ is the Jacobian $ \frac{\partial - * J}{\partial y} $. The input parameters are the time, $ \tau $, and a + * f}{\partial y} $. The input parameters are the time, $ \tau $, and a * vector. The output is the value of function at this point. * evolve_one_time_step returns the time at the end of the time step. */ @@ -333,10 +333,12 @@ namespace TimeStepping * parameters are the time t and the vector y and the output is value of f * at this point. @p id_minus_tau_J_inverse is a function that computes $ * inv(I-\tau J)$ where $ I $ is the identity matrix, $ \tau $ is given, - * and $ J $ is the Jacobian $ \frac{\partial J}{\partial y} $. The input + * and $ J $ is the Jacobian $ \frac{\partial f}{\partial y} $. The input * parameter are the time, $ \tau $, and a vector. The output is the value * of function at this point. evolve_one_time_step returns the time at the * end of the time step. + * + * @note @p id_minus_tau_J_inverse is ignored since the method is explicit. */ double evolve_one_time_step( @@ -435,7 +437,7 @@ namespace TimeStepping * parameters are the time t and the vector y and the output is value of f * at this point. @p id_minus_tau_J_inverse is a function that computes $ * inv(I-\tau J)$ where $ I $ is the identity matrix, $ \tau $ is given, - * and $ J $ is the Jacobian $ \frac{\partial J}{\partial y} $. The input + * and $ J $ is the Jacobian $ \frac{\partial f}{\partial y} $. The input * parameters are the time, $ \tau $, and a vector. The output is the value * of function at this point. evolve_one_time_step returns the time at the * end of the time step. @@ -558,7 +560,7 @@ namespace TimeStepping * parameters are the time t and the vector y and the output is value of f * at this point. @p id_minus_tau_J_inverse is a function that computes $ * (I-\tau J)^{-1}$ where $ I $ is the identity matrix, $ \tau $ is given, - * and $ J $ is the Jacobian $ \frac{\partial J}{\partial y} $. The input + * and $ J $ is the Jacobian $ \frac{\partial f}{\partial y} $. The input * parameters this function receives are the time, $ \tau $, and a vector. * The output is the value of function at this point. evolve_one_time_step * returns the time at the end of the time step. @@ -642,13 +644,6 @@ namespace TimeStepping VectorType & tendency, VectorType & residual) const; - /** - * When using SDIRK, there is no need to compute the linear combination of - * the stages. Thus, when this flag is true, the linear combination is - * skipped. - */ - bool skip_linear_combi; - /** * Maximum number of iterations of the Newton solver. */ @@ -722,10 +717,12 @@ namespace TimeStepping * parameters are the time t and the vector y and the output is value of f * at this point. @p id_minus_tau_J_inverse is a function that computes $ * inv(I-\tau J)$ where $ I $ is the identity matrix, $ \tau $ is given, - * and $ J $ is the Jacobian $ \frac{\partial J}{\partial y} $. The input + * and $ J $ is the Jacobian $ \frac{\partial f}{\partial y} $. The input * parameters are the time, $ \tau $, and a vector. The output is the * value of function at this point. evolve_one_time_step returns the time * at the end of the time step. + * + * @note @p id_minus_tau_J_inverse is ignored since the method is explicit. */ double evolve_one_time_step( @@ -833,7 +830,7 @@ namespace TimeStepping * If the flag is true, the last stage is the same as the first stage and * one evaluation of f can be saved. */ - bool last_same_as_first; + bool last_same_as_first = false; /** * Butcher tableau coefficients. @@ -849,7 +846,7 @@ namespace TimeStepping * If the last_same_as_first flag is set to true, the last stage is saved * and reused as the first stage of the next time step. */ - VectorType *last_stage; + VectorType *last_stage = nullptr; /** * Status structure of the object. diff --git a/include/deal.II/base/time_stepping.templates.h b/include/deal.II/base/time_stepping.templates.h index 6bb8b682d8..cfc0ff43f1 100644 --- a/include/deal.II/base/time_stepping.templates.h +++ b/include/deal.II/base/time_stepping.templates.h @@ -479,7 +479,6 @@ namespace TimeStepping const unsigned int max_it, const double tolerance) : RungeKutta() - , skip_linear_combi(false) , max_it(max_it) , tolerance(tolerance) { @@ -572,13 +571,9 @@ namespace TimeStepping // Compute the different stages needed. compute_stages(f, id_minus_tau_J_inverse, t, delta_t, y, f_stages); - // If necessary, compute the linear combinations of the stages. - if (skip_linear_combi == false) - { - y = old_y; - for (unsigned int i = 0; i < this->n_stages; ++i) - y.sadd(1., delta_t * this->b[i], f_stages[i]); - } + y = old_y; + for (unsigned int i = 0; i < this->n_stages; ++i) + y.sadd(1., delta_t * this->b[i], f_stages[i]); return (t + delta_t); }