]> https://gitweb.dealii.org/ - code-gallery.git/commitdiff
Update to deal.II 9.1: Quasi_static_Finite_strain_Compressible_Elasticity
authorJean-Paul Pelteret <jppelteret@gmail.com>
Sat, 23 May 2020 16:24:04 +0000 (18:24 +0200)
committerJean-Paul Pelteret <jppelteret@gmail.com>
Sat, 23 May 2020 17:32:21 +0000 (19:32 +0200)
Quasi_static_Finite_strain_Compressible_Elasticity/CMakeLists.txt
Quasi_static_Finite_strain_Compressible_Elasticity/cook_membrane.cc

index b6b64346cd43e7c6e94b75d5397016153d603866..6f40451d903fb8f9f9f623317462b8a12477cd73 100644 (file)
@@ -20,7 +20,7 @@ SET(CLEAN_UP_FILES
 
 CMAKE_MINIMUM_REQUIRED(VERSION 2.8.8)
 
-FIND_PACKAGE(deal.II 8.5 QUIET
+FIND_PACKAGE(deal.II 9.1 QUIET
   HINTS ${deal.II_DIR} ${DEAL_II_DIR} ../ ../../ $ENV{DEAL_II_DIR}
   )
 IF(NOT ${deal.II_FOUND})
index 7a40efe22ab2584b45de877eb329f78fa131bdad..9b19f0c68f25b8ee9a2db9726206568f5e65963c 100644 (file)
@@ -32,7 +32,6 @@
 #include <deal.II/base/timer.h>
 #include <deal.II/base/work_stream.h>
 #include <deal.II/base/quadrature_point_data.h>
-#include <deal.II/base/std_cxx11/shared_ptr.h>
 
 #include <deal.II/dofs/dof_renumbering.h>
 #include <deal.II/dofs/dof_tools.h>
@@ -50,6 +49,7 @@
 #include <deal.II/fe/mapping_q_eulerian.h>
 #include <deal.II/fe/mapping_q.h>
 
+#include <deal.II/lac/affine_constraints.h>
 #include <deal.II/lac/block_sparse_matrix.h>
 #include <deal.II/lac/block_vector.h>
 #include <deal.II/lac/dynamic_sparsity_pattern.h>
@@ -57,7 +57,6 @@
 #include <deal.II/lac/precondition_selector.h>
 #include <deal.II/lac/solver_cg.h>
 #include <deal.II/lac/sparse_direct.h>
-#include <deal.II/lac/constraint_matrix.h>
 
 #include <deal.II/numerics/data_out.h>
 #include <deal.II/numerics/vector_tools.h>
@@ -76,6 +75,7 @@
 
 #include <iostream>
 #include <fstream>
+#include <memory>
 
 
 // We then stick everything that relates to this tutorial program into a
@@ -892,10 +892,10 @@ namespace Cook_Membrane
     const unsigned int               n_q_points_f;
 
     // Objects that store the converged solution and right-hand side vectors,
-    // as well as the tangent matrix. There is a ConstraintMatrix object used
+    // as well as the tangent matrix. There is a AffineConstraints object used
     // to keep track of constraints.  We make use of a sparsity pattern
     // designed for a block system.
-    ConstraintMatrix                 constraints;
+    AffineConstraints<double>        constraints;
     BlockSparsityPattern             sparsity_pattern;
     BlockSparseMatrix<double>        tangent_matrix;
     BlockVector<double>              system_rhs;
@@ -1591,7 +1591,7 @@ namespace Cook_Membrane
     void
     copy_local_to_global_ASM(const PerTaskData_ASM &data)
     {
-      const ConstraintMatrix &constraints = data.solid->constraints;
+      const AffineConstraints<double> &constraints = data.solid->constraints;
       BlockSparseMatrix<double> &tangent_matrix = const_cast<Solid<dim,NumberType> *>(data.solid)->tangent_matrix;
       BlockVector<double> &system_rhs =  const_cast<Solid<dim,NumberType> *>(data.solid)->system_rhs;
 
@@ -2139,8 +2139,9 @@ namespace Cook_Membrane
       std::cout << " SLV " << std::flush;
       if (parameters.type_lin == "CG")
         {
-          const int solver_its = tangent_matrix.block(u_dof, u_dof).m()
-                                 * parameters.max_iterations_lin;
+          const int solver_its = static_cast<unsigned int>(
+                                    tangent_matrix.block(u_dof, u_dof).m()
+                                    * parameters.max_iterations_lin);
           const double tol_sol = parameters.tol_lin
                                  * system_rhs.block(u_dof).l2_norm();
 

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.