]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Fixed AS comments. 5167/head
authorLuca Heltai <luca.heltai@sissa.it>
Thu, 5 Oct 2017 10:25:44 +0000 (12:25 +0200)
committerLuca Heltai <luca.heltai@sissa.it>
Thu, 5 Oct 2017 10:36:31 +0000 (12:36 +0200)
include/deal.II/sundials/kinsol.h
source/sundials/kinsol.cc

index b31a0affc3e703227eb1e01826082b2081a16de0..0054d646e3174e1163aa517a58979964aad740ef 100644 (file)
@@ -275,7 +275,7 @@ namespace SUNDIALS
         dq_relative_error(dq_relative_error),
         maximum_beta_failures(maximum_beta_failures),
         anderson_subspace_size(anderson_subspace_size)
-      {};
+      {}
 
       /**
        * Add all AdditionalData() parameters to the given ParameterHandler
@@ -560,15 +560,9 @@ namespace SUNDIALS
      * successive calls to solve_jacobian_systems() lead to better convergence
      * in the Newton process.
      *
-     * If you do not specify a solve_jacobian_system() function, then a fixed
-     * point iteration is used instead of a Newton method. Notice that this may
-     * not converge, or may converge very slowly.
-     *
-     * The Jacobian $J$ should be (an approximation of) the system Jacobian
-     * \f[
-     *   J = M - \gamma \frac{\partial f_I}{\partial y}
-     * \f]
-     * evaluated at `t`, `ycur`. `fcur` is $f_I(t,ycur)$.
+     * If you do not specify a solve_jacobian_system() function, then only a
+     * fixed point iteration strategy can be used. Notice that this may not
+     * converge, or may converge very slowly.
      *
      * A call to this function should store in `dst` the result of $J^{-1}$
      * applied to `src`, i.e., `J*dst = src`. It is the users responsability to
@@ -577,12 +571,11 @@ namespace SUNDIALS
      *
      * Arguments to the function are
      *
-     * @param[in] t  the current time
-     * @param[in] gamma  the current factor to use in the Jacobian computation
      * @param[in] ycur  is the current $y$ vector for the current KINSOL internal step
      * @param[in] fcur  is the current value of the implicit right-hand side at ycur,
      *        $f_I (t_n, ypred)$.
-     *
+     * @param[in] rhs  the system right hand side to solve for
+     * @param[out] dst the solution of $A^{-1} * src$
      *
      * This function should return:
      * - 0: Success
@@ -661,12 +654,10 @@ namespace SUNDIALS
      */
     N_Vector f_scale;
 
-#ifdef DEAL_II_WITH_MPI
     /**
      * MPI communicator. SUNDIALS solver runs happily in parallel.
      */
     MPI_Comm communicator;
-#endif
 
     /**
      * Memory pool of vectors.
index fbffd73efbe15f95d2a6cbd2e7710f06bd72a4b8..178940a478283e0c6a10c482608ec3b18f18037d 100644 (file)
@@ -153,7 +153,9 @@ namespace SUNDIALS
   {
     if (kinsol_mem)
       KINFree(&kinsol_mem);
+#ifdef DEAL_II_WITH_MPI
     MPI_Comm_free(&communicator);
+#endif
   }
 
 
@@ -162,7 +164,6 @@ namespace SUNDIALS
   unsigned int KINSOL<VectorType>::solve(VectorType &initial_guess_and_solution)
   {
     unsigned int system_size = initial_guess_and_solution.size();
-    unsigned int local_system_size = system_size;
 
     // The solution is stored in
     // solution. Here we take only a
@@ -170,8 +171,8 @@ namespace SUNDIALS
 #ifdef DEAL_II_WITH_MPI
     if (is_serial_vector<VectorType>::value == false)
       {
-        IndexSet is = initial_guess_and_solution.locally_owned_elements();
-        local_system_size = is.n_elements();
+        const IndexSet is = initial_guess_and_solution.locally_owned_elements();
+        const unsigned int local_system_size = is.n_elements();
 
         solution        = N_VNew_Parallel(communicator,
                                           local_system_size,
@@ -213,6 +214,7 @@ namespace SUNDIALS
     kinsol_mem = KINCreate();
 
     int status = KINInit(kinsol_mem, t_kinsol_function<VectorType> , solution);
+    (void) status;
     AssertKINSOL(status);
 
     status = KINSetUserData(kinsol_mem, (void *) this);

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.