]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
updated results
authormcbride <mcbride@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 17 Feb 2012 07:12:41 +0000 (07:12 +0000)
committermcbride <mcbride@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 17 Feb 2012 07:12:41 +0000 (07:12 +0000)
git-svn-id: https://svn.dealii.org/trunk@25105 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/examples/step-44/doc/results.dox

index dd5dfe7dd7bb229d66b964d30e1a6bee41f56183..e1540a9d3e5260d78d768ca3df6414cae049dae9 100644 (file)
 <h1>Results</h1>
 
+Firstly, a comparison of a batch of results with that in the literature demonstrates that the program works as expected.
+<table align="center" class="tutorial" cellspacing="3" cellpadding="3">
+  <tr>
+     <td align="center">
+        @image html "Q1-P0_convergence.png"
+       <p align="center">
+        Convergence for the Q1-P0-P0 formulation.
+       </p>
+    </td>
+    <td align="center">
+        @image html "Q2-P1_convergence.png"
+       <p align="center">
+        Convergence for the Q2-P1-P1 formulation.
+       </p>
+    </td>
+  </tr>
+</table>
+Using the appropriate material and loading parameters, a set of results matches that presented in <em>Reese (2000) </em> is produced. Both schemes demonstrate good convergence properties upon grid refinement, with an accurate measure of the centre-point vertical displacement attained within a few grid refinements. The lower order formulation typically overestimates the displacement solution at low levels of refinement, while the higher order interpolation scheme underestimates it, but be a lesser degree. This result gives confidence that the program produces the correct output.
+
+A typical output generated by running the problem looks is shown below. The particular case demonstrated is that of the Q2-P1-P1 element.
+
+@code
+Grid:
+        Reference volume: 1e-09
+Triangulation:
+        Number of active cells: 512
+        Number of degrees of freedom: 18835
+    Setting up quadrature point data...
+
+Timestep 1 @ 0.1s
+___________________________________________________________________________________________________________________________________________________________
+                 SOLVER STEP                   |  LIN_IT   LIN_RES    RES_NORM     RES_U     RES_P      RES_J     NU_NORM      NU_U       NU_P       NU_J 
+___________________________________________________________________________________________________________________________________________________________
+  0  ASM_R  ASM_K  CST  ASM_SC  SLV  PP  UQPH  |    1163  1.319e-06  1.000e+00  1.000e+00  0.000e+00  0.000e+00  1.000e+00  1.000e+00  1.000e+00  1.000e+00  
+  1  ASM_R  ASM_K  CST  ASM_SC  SLV  PP  UQPH  |     551  3.047e-03  3.752e+00  8.173e-02  9.198e-11  5.080e+00  5.181e+00  2.521e+00  5.181e+00  1.812e+04  
+  2  ASM_R  ASM_K  CST  ASM_SC  SLV  PP  UQPH  |     690  5.008e-03  2.918e+00  2.562e+00  2.624e-12  1.892e+00  4.254e+00  2.094e+00  4.254e+00  1.638e+03  
+  3  ASM_R  ASM_K  CST  ASM_SC  SLV  PP  UQPH  |     661  2.053e-03  1.828e+00  1.804e+00  1.654e-12  3.941e-01  1.248e+00  4.299e-01  1.248e+00  1.555e+03  
+  4  ASM_R  ASM_K  CST  ASM_SC  SLV  PP  UQPH  |     805  1.476e-04  1.209e-01  1.197e-01  8.815e-14  2.373e-02  6.740e-02  2.693e-02  6.740e-02  1.539e+02  
+  5  ASM_R  ASM_K  CST  ASM_SC  SLV  PP  UQPH  |    1076  6.240e-07  5.569e-04  5.515e-04  3.003e-16  1.053e-04  3.640e-04  1.310e-04  3.640e-04  6.834e-01  
+  6  ASM_R  ASM_K  CST  ASM_SC  SLV  PP  UQPH  |    1346  1.279e-11  4.610e-07  4.610e-07  8.115e-21  2.062e-09  4.382e-07  4.221e-08  4.382e-07  1.330e-05  
+  7  ASM_R  CONVERGED! 
+___________________________________________________________________________________________________________________________________________________________
+Relative errors:
+Displacement:  4.221e-08
+Force:                 9.444e-12
+Dilatation:    6.731e-08
+v / V_0:       1.000e-09 / 1.000e-09 = 1.000e+00
+
+[...]
+
+Timestep 10 @ 1.000e+00s
+___________________________________________________________________________________________________________________________________________________________
+                 SOLVER STEP                   |  LIN_IT   LIN_RES    RES_NORM     RES_U     RES_P      RES_J     NU_NORM      NU_U       NU_P       NU_J 
+___________________________________________________________________________________________________________________________________________________________
+  0  ASM_R  ASM_K  CST  ASM_SC  SLV  PP  UQPH  |    1141  1.313e-06  1.000e+00  1.000e+00  1.000e+00  1.000e+00  1.000e+00  1.000e+00  1.000e+00  1.000e+00  
+  1  ASM_R  ASM_K  CST  ASM_SC  SLV  PP  UQPH  |     969  3.959e-05  1.200e-01  1.200e-01  2.955e+11  4.575e+07  7.143e-02  7.434e-02  7.143e-02  7.143e-02  
+  2  ASM_R  ASM_K  CST  ASM_SC  SLV  PP  UQPH  |    1153  1.980e-07  2.523e-03  2.523e-03  1.444e+09  5.757e+04  3.790e-03  1.340e-03  3.790e-03  3.794e-03  
+  3  ASM_R  ASM_K  CST  ASM_SC  SLV  PP  UQPH  |    1340  2.440e-10  1.480e-06  1.480e-06  1.234e+06  2.078e+02  1.387e-06  5.401e-07  1.387e-06  1.387e-06  
+  4  ASM_R  CONVERGED! 
+___________________________________________________________________________________________________________________________________________________________
+Relative errors:
+Displacement:  5.401e-07
+Force:                 1.804e-10
+Dilatation:    7.492e-07
+v / V_0:       9.999e-10 / 1.000e-09 = 9.999e-01
+@endcode
+
+It is clear that, using the Newton-Raphson method, quadratic convergence of the solution is obtained. Other than the for the first timestep, solution convergence is achieved within 5 Newton increments. The end displacement L2-norm is several orders of magnitude less than the geometry scale.
+
+Using the Timer class, we can discern which parts of the code require the highest computational expense. For a case with a large number of degrees-of-freedom, a typical output of the Timer may look like
+
+@code
++---------------------------------------------+------------+------------+
+| Total wallclock time elapsed since start    | 9.874e+02s |            |
+|                                             |            |            |
+| Section                         | no. calls |  wall time | % of total |
++---------------------------------+-----------+------------+------------+
+| Assemble system right-hand side |        53 | 1.727e+00s |  1.75e-01% |
+| Assemble tangent matrix         |        43 | 2.707e+01s |  2.74e+00% |
+| Linear solver                   |        43 | 9.248e+02s |  9.37e+01% |
+| Linear solver postprocessing    |        43 | 2.743e-02s |  2.78e-03% |
+| Perform static condensation     |        43 | 1.437e+01s |  1.46e+00% |
+| Setup system                    |         1 | 3.897e-01s |  3.95e-02% |
+| Update QPH data                 |        43 | 5.770e-01s |  5.84e-02% |
++---------------------------------+-----------+------------+------------+
+@endcode
+
+Much of the code in the tutorial has been developed based on the optimisations described, discussed and demonstrated in Step-18 and others. With over 93% of the time being spent in the linear solver, it is obvious that it may be necessary to invest in a better solver for large three-dimensional problems. The SSOR preconditioner is not multi-threaded but is effective for this class of solid problems. It may be beneficial to investigate the use of another solver such as those available through Trilinos.
+
+Using ParaView, we visualise the results for the vertical displacement of the centre-point of the traction surface and pressure field from two of the cases. 
+
+The first case shown is that using the coarsest grid using the lowest order interpolation method, namely the Q1-P0-P0 element. It is clear that the coarse spatial discretisation coupled with large displacements leads to a low quality solution. Additionally, the pressure difference between elements is very large. The constant pressure field on the element ensures that the large pressure gradient is not captured. However, it should be noted that locking that would be present in a standard Q1 displacement formulation still does not arise in this poorly discretised case. The images below show the results for the tri-linear displacement formulation with p/p0=80.
 
 <table align="center" class="tutorial" cellspacing="3" cellpadding="3">
   <tr>
     <td align="center">
-        @image html "step-44.Q1P0_Ref_2.png"
+        @image html "Q1-P0_gr_1_p_ratio_80-displacement.png"
        <p align="center">
-        Two refinement levels for a tri-linear displacement formulation.
+        Z-displacement solution.
        </p>
     </td>
     <td align="center">
-        @image html "step-44.Q2P1_Ref_1.png"
+        @image html "Q1-P0_gr_1_p_ratio_80-pressure.png"
+       <p align="center">
+        Discontinuous piece-wise constant pressure field.
+       </p>
+    </td>
+     <td align="center">
+        @image html "Q1-P0_gr_1_p_ratio_80-dilatation.png"
        <p align="center">
-        One refinement level for a tri-quadratic displacement formulation.
+        Discontinuous piece-wise constant dilatation field.
        </p>
     </td>
   </tr>
+</table>
+Although the displacement solution is coarse, the final vertical displacement of the node at the centre of the top surface of the block is within 12.5% of the converged solution. The pressure solution is very course and has large jumps between adjacent cells. It is clear that the volume nearest the applied traction undergoes compression while the extents of the domain are in a state of dilatation. The dilatation solution field and pressure field are clearly linked, with positive dilatation indicating regions of positive hyprodstatic pressure and negative showing regions placed in compression. It should be noted that in these, a compressive pressure generated in each element has a negative sign while expansive pressure takes a positive sign. This stems from the definition of the volumetric strain energy function and is opposite to the physically realistic interpretation of pressure.
+
+Both spatial refinement and the use of a higher-order interpolation scheme results in a high quality displacement solution. Three grid refinements coupled with a Q2-P1-P1 formulation produce a result that clearly captures the mechanics of the problem. The deformation of the traction surface is well resolved. The images below show the result for tri-quadratic displacement solution with p/p0=80.
+
+<table align="center" class="tutorial" cellspacing="3" cellpadding="3">
 
   <tr>
     <td align="center">
-        @image html "step-44.Q1P0_Ref_4.png"
+        @image html "Q2-P1_gr_3_p_ratio_80-displacement.png"
+       <p align="center">
+        Z-displacement solution. 
+       </p>
+    </td>
+    <td align="center">
+        @image html "Q2-P1_gr_3_p_ratio_80-pressure.png"
        <p align="center">
-        Four refinement levels for a tri-linear displacement formulation.
+        Discontinuous linear pressure field.
        </p>
     </td>
     <td align="center">
-        @image html "step-44.Q2P1_Ref_3.png"
+        @image html "Q2-P1_gr_3_p_ratio_80-dilatation.png"
+       <p align="center">
+        Discontinuous linear dilatation field.
+       </p>
+    </td>
+  </tr>
+</table>
+Here we observe the true extent of the applied traction, with the maximum force being applied at the central point of the surface causing the largest compression. Even though very high strains are experienced in the domain, especially at the junction between the regions of the upper surface when the traction is and isn't applied, the solution remains accurate. The pressure field is captured in far greater detail than before. There is a clear distinction and transition between regions of compression and tension, and the linear approximation of the pressure field allows a refined visualisation of the pressure at a sub-element scale. It should however be noted that the pressure field remains discontinuous and could be smoothed on a continuous grid for the post-processing purposes.
+
+A brief analysis of the results for an selection of input parameters for the described problem demonstrates that the three-field formulation is effective in circumventing volumetric-locking of highly-incompressible media. The mixed formulation is able to accurately simulate the displacement of a near-incompressible block under compression. The command-line output indicates that the volumetric change under extreme compression resulted in less than 0.01% volume change for a Poisson's ratio of 0.4999.
+
+In terms of run-time, the Q2-P1-P1 formulation tends to be a bit more computationally expensive than the Q1-P0-P0 for a similar number of degree-of-freedom (produced by adding an extra grid refinement level for the lower order interpolation method). This is demonstrated in the graph below for a batch of tests run serially on a single 4-core (8 thread) machine. 
+<table align="center" class="tutorial" cellspacing="3" cellpadding="3">
+  <tr>
+     <td align="center">
+        @image html "Normalised_runtime.png"
        <p align="center">
-        Three refinement levels for a tri-quadratic displacement formulation.
+        Runtime on a 4-core machine, normalised against the lowest grid resolution Q1-P0 solution that utilised a SSOR preconditioner.
        </p>
     </td>
   </tr>
 </table>
+The increase in computational time for the higher order method is likely due to the increased band-width required for the higher order elements. As previously mentioned, the use of a better solver and precondtioner may mitigate the expense of using a higher-order formulation. It was observed that using the multithreaded Jacobi preconditioner does reduced computational runtime by up to 25% in comparison to the single-thread SSOR preconditioner. However, it is the author's experience that the Jacobi method of preconditioning may not be suitable for some finite-strain problems involving different constitutive models.
 
 <a name="extensions"></a>
 <h3>Possibilities for extensions</h3>
 
+There are a number of obvious extensions for this work.
+
+Firstly, an additional constraint could be added to the free-energy function in order to enforce a high degree of incompressibility in materials. An additional Lagrange multiplier would be introduced, but this could most easily be dealt with using the principle of augented Lagrange multipliers. This is demonstrated in <em>Simo and Taylor (1991) </em>.
+
+The constitutive relationship used in this model is very basic. It may be beneficial to split the the material class into two separate classes, one dealing with the volumetric response and the other the isochoric response, and produce a generic materials class that would allow for the addition of more complex material models. Such models could include other hyperelastic materials, plasticity and viscoelastic materials and others. 
+
+The program has been developed for solving problems on single-node multicore machines. With a little effort, the program could be extended to a large-scale computing environment through the use of Petsc or Trilinos, using a similar technique to that demonstrated in step-40. This would mostly involve changes to the setup, assembly, Point History and linear solver routines.
+
+As this program solves assumes quasi-static equilibrium, extension to include dynamic effects would be necessary to study problems where inertial effects are important, e.g. problems involving impact.
+
+Load and solution limiting procedures may be necessary for highly nonlinear problems. It is possible to add a linesearch algorithm to limit the step size within a Newton increment to ensure optimum convergence. It may also be necessary to use a load limiting method, such as the Riks method, to solve unstable problems involving geometric nonlinearity such as buckling and snap-through.
+
+Many physical problems involve contact. It is possible to include the effect of frictional or frictionless contact between objects into this program. This would involve that addition of an extra term in the free-energy functional and therefore addition to the assembly routine would be required in addition to managing the contact problem (detection and stress calculations) itself.
+
+Finally, adaptive mesh refinement, as demonstrated in step-18, could provide additional solution accuracy.

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.