From: Jean-Paul Pelteret Date: Mon, 28 Dec 2015 08:22:11 +0000 (+0100) Subject: Added parameter file; code cleanup; fixed small bug in post-processing X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=11ef1804fb9e3d41a5a50f44d482bdc651e00c95;p=code-gallery.git Added parameter file; code cleanup; fixed small bug in post-processing --- diff --git a/Quasi_static_Finite_strain_Compressible_Elasticity/cook_membrane.cc b/Quasi_static_Finite_strain_Compressible_Elasticity/cook_membrane.cc index 6dffb8b..fbfc824 100644 --- a/Quasi_static_Finite_strain_Compressible_Elasticity/cook_membrane.cc +++ b/Quasi_static_Finite_strain_Compressible_Elasticity/cook_membrane.cc @@ -40,7 +40,6 @@ #include #include #include -#include // TEMP #include #include @@ -48,6 +47,7 @@ #include #include #include +#include #include #include @@ -1747,9 +1747,18 @@ Point grid_y_transform (const Point &pt_in) for (unsigned int v=0; v::vertices_per_cell; ++v) if (cell->vertex(v).distance(soln_pt) < 1e-6) { - // Create some object to help us extract the solution value at - // the desired point - const Quadrature soln_qrule (soln_pt); + // Extract y-component of solution at the given point + // This point is coindicent with a vertex, so we can + // extract it directly as we're using FE_Q finite elements + // that have support at the vertices + vertical_tip_displacement = solution_n(cell->vertex_dof_index(v,u_dof+1)); + + // Sanity check using alternate method to extract the solution + // at the given point. To do this, we must create an FEValues instance + // to help us extract the solution value at the desired point + const MappingQ mapping (parameters.poly_degree); + const Point qp_unit = mapping.transform_real_to_unit_cell(cell,soln_pt); + const Quadrature soln_qrule (qp_unit); AssertThrow(soln_qrule.size() == 1, ExcInternalError()); FEValues fe_values_soln (fe, soln_qrule, update_values); fe_values_soln.reinit(cell); @@ -1758,12 +1767,7 @@ Point grid_y_transform (const Point &pt_in) std::vector< Tensor<1,dim> > soln_values (soln_qrule.size()); fe_values_soln[u_fe].get_function_values(solution_n, soln_values); - vertical_tip_displacement = soln_values[0][u_dof+1]; - - // Sanity Check - for (unsigned int v=0; v::vertices_per_cell; ++v) - if (cell->vertex(v).distance(soln_pt) < 1e-6) - vertical_tip_displacement_check = solution_n(cell->vertex_dof_index(v,u_dof+1)); + vertical_tip_displacement_check = soln_values[0][u_dof+1]; break; } diff --git a/Quasi_static_Finite_strain_Compressible_Elasticity/parameters.prm b/Quasi_static_Finite_strain_Compressible_Elasticity/parameters.prm new file mode 100644 index 0000000..595f06a --- /dev/null +++ b/Quasi_static_Finite_strain_Compressible_Elasticity/parameters.prm @@ -0,0 +1,68 @@ +# Listing of Parameters +# --------------------- +subsection Finite element system + # Displacement system polynomial order + set Polynomial degree = 1 + + # Gauss quadrature order + set Quadrature order = 2 +end + + +subsection Geometry + # Number of elements per long edge of the beam + set Elements per edge = 32 + + # Global grid scaling factor + set Grid scale = 1e-3 +end + + +subsection Linear solver + # Linear solver iterations (multiples of the system matrix size) + set Max iteration multiplier = 1 + + # Linear solver residual (scaled by residual norm) + set Residual = 1e-6 + + # Preconditioner type + set Preconditioner type = ssor + + # Preconditioner relaxation value + set Preconditioner relaxation = 0.65 + + # Type of solver used to solve the linear system + set Solver type = Direct +end + + +subsection Material properties + # Poisson's ratio + set Poisson's ratio = 0.3 + + # Shear modulus + set Shear modulus = 0.4225e6 +end + + +subsection Nonlinear solver + # Number of Newton-Raphson iterations allowed + set Max iterations Newton-Raphson = 10 + + # Displacement error tolerance + set Tolerance displacement = 1.0e-6 + + # Force residual tolerance + set Tolerance force = 1.0e-9 +end + + +subsection Time + # End time + set End time = 1 + + # Time step size + set Time step size = 0.1 +end + +