Reference documentation for deal.II version 9.2.0
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
step-8.h
Go to the documentation of this file.
1  false);
750  * sparsity_pattern.copy_from(dsp);
751  *
752  * system_matrix.reinit(sparsity_pattern);
753  * }
754  *
755  *
756  * @endcode
757  *
758  *
759  * <a name="ElasticProblemassemble_system"></a>
760  * <h4>ElasticProblem::assemble_system</h4>
761  *
762 
763  *
764  * The big changes in this program are in the creation of matrix and right
765  * hand side, since they are problem-dependent. We will go through that
766  * process step-by-step, since it is a bit more complicated than in previous
767  * examples.
768  *
769 
770  *
771  * The first parts of this function are the same as before, however: setting
772  * up a suitable quadrature formula, initializing an FEValues object for the
773  * (vector-valued) finite element we use as well as the quadrature object,
774  * and declaring a number of auxiliary arrays. In addition, we declare the
775  * ever same two abbreviations: <code>n_q_points</code> and
776  * <code>dofs_per_cell</code>. The number of degrees of freedom per cell we
777  * now obviously ask from the composed finite element rather than from the
778  * underlying scalar Q1 element. Here, it is <code>dim</code> times the
779  * number of degrees of freedom per cell of the Q1 element, though this is
780  * not explicit knowledge we need to care about:
781  *
782  * @code
783  * template <int dim>
784  * void ElasticProblem<dim>::assemble_system()
785  * {
786  * QGauss<dim> quadrature_formula(fe.degree + 1);
787  *
788  * FEValues<dim> fe_values(fe,
789  * quadrature_formula,
792  *
793  * const unsigned int dofs_per_cell = fe.dofs_per_cell;
794  * const unsigned int n_q_points = quadrature_formula.size();
795  *
796  * FullMatrix<double> cell_matrix(dofs_per_cell, dofs_per_cell);
797  * Vector<double> cell_rhs(dofs_per_cell);
798  *
799  * std::vector<types::global_dof_index> local_dof_indices(dofs_per_cell);
800  *
801  * @endcode
802  *
803  * As was shown in previous examples as well, we need a place where to
804  * store the values of the coefficients at all the quadrature points on a
805  * cell. In the present situation, we have two coefficients, lambda and
806  * mu.
807  *
808  * @code
809  * std::vector<double> lambda_values(n_q_points);
810  * std::vector<double> mu_values(n_q_points);
811  *
812  * @endcode
813  *
814  * Well, we could as well have omitted the above two arrays since we will
815  * use constant coefficients for both lambda and mu, which can be declared
816  * like this. They both represent functions always returning the constant
817  * value 1.0. Although we could omit the respective factors in the
818  * assemblage of the matrix, we use them here for purpose of
819  * demonstration.
820  *
821  * @code
823  *
824  * @endcode
825  *
826  * Like the two constant functions above, we will call the function
827  * right_hand_side just once per cell to make things simpler.
828  *
829  * @code
830  * std::vector<Tensor<1, dim>> rhs_values(n_q_points);
831  *
832  * @endcode
833  *
834  * Now we can begin with the loop over all cells:
835  *
836  * @code
837  * for (const auto &cell : dof_handler.active_cell_iterators())
838  * {
839  * cell_matrix = 0;
840  * cell_rhs = 0;
841  *
842  * fe_values.reinit(cell);
843  *
844  * @endcode
845  *
846  * Next we get the values of the coefficients at the quadrature
847  * points. Likewise for the right hand side:
848  *
849  * @code
850  * lambda.value_list(fe_values.get_quadrature_points(), lambda_values);
851  * mu.value_list(fe_values.get_quadrature_points(), mu_values);
852  * right_hand_side(fe_values.get_quadrature_points(), rhs_values);
853  *
854  * @endcode
855  *
856  * Then assemble the entries of the local stiffness matrix and right
857  * hand side vector. This follows almost one-to-one the pattern
858  * described in the introduction of this example. One of the few
859  * comments in place is that we can compute the number
860  * <code>comp(i)</code>, i.e. the index of the only nonzero vector
861  * component of shape function <code>i</code> using the
862  * <code>fe.system_to_component_index(i).first</code> function call
863  * below.
864  *
865 
866  *
867  * (By accessing the <code>first</code> variable of the return value
868  * of the <code>system_to_component_index</code> function, you might
869  * already have guessed that there is more in it. In fact, the
870  * function returns a <code>std::pair@<unsigned int, unsigned
871  * int@></code>, of which the first element is <code>comp(i)</code>
872  * and the second is the value <code>base(i)</code> also noted in the
873  * introduction, i.e. the index of this shape function within all the
874  * shape functions that are nonzero in this component,
875  * i.e. <code>base(i)</code> in the diction of the introduction. This
876  * is not a number that we are usually interested in, however.)
877  *
878 
879  *
880  * With this knowledge, we can assemble the local matrix
881  * contributions:
882  *
883  * @code
884  * for (const unsigned int i : fe_values.dof_indices())
885  * {
886  * const unsigned int component_i =
887  * fe.system_to_component_index(i).first;
888  *
889  * for (const unsigned int j : fe_values.dof_indices())
890  * {
891  * const unsigned int component_j =
892  * fe.system_to_component_index(j).first;
893  *
894  * for (const unsigned int q_point :
895  * fe_values.quadrature_point_indices())
896  * {
897  * cell_matrix(i, j) +=
898  * @endcode
899  *
900  * The first term is @f$\lambda \partial_i u_i, \partial_j
901  * v_j) + (\mu \partial_i u_j, \partial_j v_i)@f$. Note
902  * that <code>shape_grad(i,q_point)</code> returns the
903  * gradient of the only nonzero component of the i-th
904  * shape function at quadrature point q_point. The
905  * component <code>comp(i)</code> of the gradient, which
906  * is the derivative of this only nonzero vector
907  * component of the i-th shape function with respect to
908  * the comp(i)th coordinate is accessed by the appended
909  * brackets.
910  *
911  * @code
912  * (
913  * (fe_values.shape_grad(i, q_point)[component_i] *
914  * fe_values.shape_grad(j, q_point)[component_j] *
915  * lambda_values[q_point])
916  * +
917  * (fe_values.shape_grad(i, q_point)[component_j] *
918  * fe_values.shape_grad(j, q_point)[component_i] *
919  * mu_values[q_point])
920  * +
921  * @endcode
922  *
923  * The second term is @f$(\mu \nabla u_i, \nabla
924  * v_j)@f$. We need not access a specific component of
925  * the gradient, since we only have to compute the
926  * scalar product of the two gradients, of which an
927  * overloaded version of <tt>operator*</tt> takes
928  * care, as in previous examples.
929  *
930 
931  *
932  * Note that by using the <tt>?:</tt> operator, we only
933  * do this if <tt>component_i</tt> equals
934  * <tt>component_j</tt>, otherwise a zero is added
935  * (which will be optimized away by the compiler).
936  *
937  * @code
938  * ((component_i == component_j) ?
939  * (fe_values.shape_grad(i, q_point) *
940  * fe_values.shape_grad(j, q_point) *
941  * mu_values[q_point]) :
942  * 0)
943  * ) *
944  * fe_values.JxW(q_point);
945  * }
946  * }
947  * }
948  *
949  * @endcode
950  *
951  * Assembling the right hand side is also just as discussed in the
952  * introduction:
953  *
954  * @code
955  * for (const unsigned int i : fe_values.dof_indices())
956  * {
957  * const unsigned int component_i =
958  * fe.system_to_component_index(i).first;
959  *
960  * for (const unsigned int q_point :
961  * fe_values.quadrature_point_indices())
962  * cell_rhs(i) += fe_values.shape_value(i, q_point) *
963  * rhs_values[q_point][component_i] *
964  * fe_values.JxW(q_point);
965  * }
966  *
967  * @endcode
968  *
969  * The transfer from local degrees of freedom into the global matrix
970  * and right hand side vector does not depend on the equation under
971  * consideration, and is thus the same as in all previous
972  * examples.
973  *
974  * @code
975  * cell->get_dof_indices(local_dof_indices);
976  * constraints.distribute_local_to_global(
977  * cell_matrix, cell_rhs, local_dof_indices, system_matrix, system_rhs);
978  * }
979  * }
980  *
981  *
982  *
983  * @endcode
984  *
985  *
986  * <a name="ElasticProblemsolve"></a>
987  * <h4>ElasticProblem::solve</h4>
988  *
989 
990  *
991  * The solver does not care about where the system of equations comes, as
992  * long as it stays positive definite and symmetric (which are the
993  * requirements for the use of the CG solver), which the system indeed
994  * is. Therefore, we need not change anything.
995  *
996  * @code
997  * template <int dim>
998  * void ElasticProblem<dim>::solve()
999  * {
1000  * SolverControl solver_control(1000, 1e-12);
1001  * SolverCG<Vector<double>> cg(solver_control);
1002  *
1003  * PreconditionSSOR<SparseMatrix<double>> preconditioner;
1004  * preconditioner.initialize(system_matrix, 1.2);
1005  *
1006  * cg.solve(system_matrix, solution, system_rhs, preconditioner);
1007  *
1008  * constraints.distribute(solution);
1009  * }
1010  *
1011  *
1012  * @endcode
1013  *
1014  *
1015  * <a name="ElasticProblemrefine_grid"></a>
1016  * <h4>ElasticProblem::refine_grid</h4>
1017  *
1018 
1019  *
1020  * The function that does the refinement of the grid is the same as in the
1021  * @ref step_6 "step-6" example. The quadrature formula is adapted to the linear elements
1022  * again. Note that the error estimator by default adds up the estimated
1023  * obtained from all components of the finite element solution, i.e., it
1024  * uses the displacement in all directions with the same weight. If we would
1025  * like the grid to be adapted to the x-displacement only, we could pass the
1026  * function an additional parameter which tells it to do so and do not
1027  * consider the displacements in all other directions for the error
1028  * indicators. However, for the current problem, it seems appropriate to
1029  * consider all displacement components with equal weight.
1030  *
1031  * @code
1032  * template <int dim>
1033  * void ElasticProblem<dim>::refine_grid()
1034  * {
1035  * Vector<float> estimated_error_per_cell(triangulation.n_active_cells());
1036  *
1037  * KellyErrorEstimator<dim>::estimate(dof_handler,
1038  * QGauss<dim - 1>(fe.degree + 1),
1039  * {},
1040  * solution,
1041  * estimated_error_per_cell);
1042  *
1044  * estimated_error_per_cell,
1045  * 0.3,
1046  * 0.03);
1047  *
1048  * triangulation.execute_coarsening_and_refinement();
1049  * }
1050  *
1051  *
1052  * @endcode
1053  *
1054  *
1055  * <a name="ElasticProblemoutput_results"></a>
1056  * <h4>ElasticProblem::output_results</h4>
1057  *
1058 
1059  *
1060  * The output happens mostly as has been shown in previous examples
1061  * already. The only difference is that the solution function is vector
1062  * valued. The DataOut class takes care of this automatically, but we have
1063  * to give each component of the solution vector a different name.
1064  *
1065 
1066  *
1067  * To do this, the DataOut::add_vector() function wants a vector of
1068  * strings. Since the number of components is the same as the number
1069  * of dimensions we are working in, we use the <code>switch</code>
1070  * statement below.
1071  *
1072 
1073  *
1074  * We note that some graphics programs have restriction on what
1075  * characters are allowed in the names of variables. deal.II therefore
1076  * supports only the minimal subset of these characters that is supported
1077  * by all programs. Basically, these are letters, numbers, underscores,
1078  * and some other characters, but in particular no whitespace and
1079  * minus/hyphen. The library will throw an exception otherwise, at least
1080  * if in debug mode.
1081  *
1082 
1083  *
1084  * After listing the 1d, 2d, and 3d case, it is good style to let the
1085  * program die if we run upon a case which we did not consider. Remember
1086  * that the Assert macro generates an exception if the condition in the
1087  * first parameter is not satisfied. Of course, the condition
1088  * <code>false</code> can never be satisfied, so the program will always
1089  * abort whenever it gets to the default statement:
1090  *
1091  * @code
1092  * template <int dim>
1093  * void ElasticProblem<dim>::output_results(const unsigned int cycle) const
1094  * {
1095  * DataOut<dim> data_out;
1096  * data_out.attach_dof_handler(dof_handler);
1097  *
1098  * std::vector<std::string> solution_names;
1099  * switch (dim)
1100  * {
1101  * case 1:
1102  * solution_names.emplace_back("displacement");
1103  * break;
1104  * case 2:
1105  * solution_names.emplace_back("x_displacement");
1106  * solution_names.emplace_back("y_displacement");
1107  * break;
1108  * case 3:
1109  * solution_names.emplace_back("x_displacement");
1110  * solution_names.emplace_back("y_displacement");
1111  * solution_names.emplace_back("z_displacement");
1112  * break;
1113  * default:
1114  * Assert(false, ExcNotImplemented());
1115  * }
1116  *
1117  * @endcode
1118  *
1119  * After setting up the names for the different components of the
1120  * solution vector, we can add the solution vector to the list of
1121  * data vectors scheduled for output. Note that the following
1122  * function takes a vector of strings as second argument, whereas
1123  * the one which we have used in all previous examples accepted a
1124  * string there. (In fact, the function we had used before would
1125  * convert the single string into a vector with only one element
1126  * and forwards that to the other function.)
1127  *
1128  * @code
1129  * data_out.add_data_vector(solution, solution_names);
1130  * data_out.build_patches();
1131  *
1132  * std::ofstream output("solution-" + std::to_string(cycle) + ".vtk");
1133  * data_out.write_vtk(output);
1134  * }
1135  *
1136  *
1137  *
1138  * @endcode
1139  *
1140  *
1141  * <a name="ElasticProblemrun"></a>
1142  * <h4>ElasticProblem::run</h4>
1143  *
1144 
1145  *
1146  * The <code>run</code> function does the same things as in @ref step_6 "step-6", for
1147  * example. This time, we use the square [-1,1]^d as domain, and we refine
1148  * it globally four times before starting the first iteration.
1149  *
1150 
1151  *
1152  * The reason for refining is a bit accidental: we use the QGauss
1153  * quadrature formula with two points in each direction for integration of the
1154  * right hand side; that means that there are four quadrature points on each
1155  * cell (in 2D). If we only refine the initial grid once globally, then there
1156  * will be only four quadrature points in each direction on the
1157  * domain. However, the right hand side function was chosen to be rather
1158  * localized and in that case, by pure chance, it happens that all quadrature
1159  * points lie at points where the right hand side function is zero (in
1160  * mathematical terms, the quadrature points happen to be at points outside
1161  * the <i>support</i> of the right hand side function). The right hand side
1162  * vector computed with quadrature will then contain only zeroes (even though
1163  * it would of course be nonzero if we had computed the right hand side vector
1164  * exactly using the integral) and the solution of the system of
1165  * equations is the zero vector, i.e., a finite element function that is zero
1166  * everywhere. In a sense, we
1167  * should not be surprised that this is happening since we have chosen
1168  * an initial grid that is totally unsuitable for the problem at hand.
1169  *
1170 
1171  *
1172  * The unfortunate thing is that if the discrete solution is constant, then
1173  * the error indicators computed by the KellyErrorEstimator class are zero
1174  * for each cell as well, and the call to
1175  * Triangulation::refine_and_coarsen_fixed_number() will not flag any cells
1176  * for refinement (why should it if the indicated error is zero for each
1177  * cell?). The grid in the next iteration will therefore consist of four
1178  * cells only as well, and the same problem occurs again.
1179  *
1180 
1181  *
1182  * The conclusion needs to be: while of course we will not choose the
1183  * initial grid to be well-suited for the accurate solution of the problem,
1184  * we must at least choose it such that it has the chance to capture the
1185  * important features of the solution. In this case, it needs to be able to
1186  * see the right hand side. Thus, we refine globally four times. (Any larger
1187  * number of global refinement steps would of course also work.)
1188  *
1189  * @code
1190  * template <int dim>
1191  * void ElasticProblem<dim>::run()
1192  * {
1193  * for (unsigned int cycle = 0; cycle < 8; ++cycle)
1194  * {
1195  * std::cout << "Cycle " << cycle << ':' << std::endl;
1196  *
1197  * if (cycle == 0)
1198  * {
1200  * triangulation.refine_global(4);
1201  * }
1202  * else
1203  * refine_grid();
1204  *
1205  * std::cout << " Number of active cells: "
1206  * << triangulation.n_active_cells() << std::endl;
1207  *
1208  * setup_system();
1209  *
1210  * std::cout << " Number of degrees of freedom: " << dof_handler.n_dofs()
1211  * << std::endl;
1212  *
1213  * assemble_system();
1214  * solve();
1215  * output_results(cycle);
1216  * }
1217  * }
1218  * } // namespace Step8
1219  *
1220  * @endcode
1221  *
1222  *
1223  * <a name="Thecodemaincodefunction"></a>
1224  * <h3>The <code>main</code> function</h3>
1225  *
1226 
1227  *
1228  * After closing the <code>Step8</code> namespace in the last line above, the
1229  * following is the main function of the program and is again exactly like in
1230  * @ref step_6 "step-6" (apart from the changed class names, of course).
1231  *
1232  * @code
1233  * int main()
1234  * {
1235  * try
1236  * {
1237  * Step8::ElasticProblem<2> elastic_problem_2d;
1238  * elastic_problem_2d.run();
1239  * }
1240  * catch (std::exception &exc)
1241  * {
1242  * std::cerr << std::endl
1243  * << std::endl
1244  * << "----------------------------------------------------"
1245  * << std::endl;
1246  * std::cerr << "Exception on processing: " << std::endl
1247  * << exc.what() << std::endl
1248  * << "Aborting!" << std::endl
1249  * << "----------------------------------------------------"
1250  * << std::endl;
1251  *
1252  * return 1;
1253  * }
1254  * catch (...)
1255  * {
1256  * std::cerr << std::endl
1257  * << std::endl
1258  * << "----------------------------------------------------"
1259  * << std::endl;
1260  * std::cerr << "Unknown exception!" << std::endl
1261  * << "Aborting!" << std::endl
1262  * << "----------------------------------------------------"
1263  * << std::endl;
1264  * return 1;
1265  * }
1266  *
1267  * return 0;
1268  * }
1269  * @endcode
1270 <a name="Results"></a><h1>Results</h1>
1271 
1272 
1273 
1274 There is not much to be said about the results of this program, other than
1275 that they look nice. All images were made using VisIt from the
1276 output files that the program wrote to disk. The first two pictures show
1277 the @f$x@f$- and @f$y@f$-displacements as scalar components:
1278 
1279 <table width="100%">
1280 <tr>
1281 <td>
1282 <img src="https://www.dealii.org/images/steps/developer/step-8.x.png" alt="">
1283 </td>
1284 <td>
1285 <img src="https://www.dealii.org/images/steps/developer/step-8.y.png" alt="">
1286 </td>
1287 </tr>
1288 </table>
1289 
1290 
1291 You can clearly see the sources of @f$x@f$-displacement around @f$x=0.5@f$ and
1292 @f$x=-0.5@f$, and of @f$y@f$-displacement at the origin.
1293 
1294 What one frequently would like to do is to show the displacement as a vector
1295 field, i.e., vectors that for each point illustrate the direction and magnitude
1296 of displacement. Unfortunately, that's a bit more involved. To understand why
1297 this is so, remember that we have just defined our finite element as a
1298 collection of two components (in <code>dim=2</code> dimensions). Nowhere have
1299 we said that this is not just a pressure and a concentration (two scalar
1300 quantities) but that the two components actually are the parts of a
1301 vector-valued quantity, namely the displacement. Absent this knowledge, the
1302 DataOut class assumes that all individual variables we print are separate
1303 scalars, and VisIt and Paraview then faithfully assume that this is indeed what it is. In
1304 other words, once we have written the data as scalars, there is nothing in
1305 these programs that allows us to paste these two scalar fields back together as a
1306 vector field. Where we would have to attack this problem is at the root,
1307 namely in <code>ElasticProblem::output_results()</code>. We won't do so here but
1308 instead refer the reader to the @ref step_22 "step-22" program where we show how to do this
1309 for a more general situation. That said, we couldn't help generating the data
1310 anyway that would show how this would look if implemented as discussed in
1311 @ref step_22 "step-22". The vector field then looks like this (VisIt and Paraview
1312 randomly select a few
1313 hundred vertices from which to draw the vectors; drawing them from each
1314 individual vertex would make the picture unreadable):
1315 
1316 <img src="https://www.dealii.org/images/steps/developer/step-8.vectors.png" alt="">
1317 
1318 
1319 We note that one may have intuitively expected the
1320 solution to be symmetric about the @f$x@f$- and @f$y@f$-axes since the @f$x@f$- and
1321 @f$y@f$-forces are symmetric with respect to these axes. However, the force
1322 considered as a vector is not symmetric and consequently neither is
1323 the solution.
1324 
1325 
1326  *
1327  *
1328 <a name="PlainProg"></a>
1329 <h1> The plain program</h1>
1330 @include "step-8.cc"
1331 */
DoFTools::nonzero
@ nonzero
Definition: dof_tools.h:244
update_quadrature_points
@ update_quadrature_points
Transformed quadrature points.
Definition: fe_update_flags.h:122
SolverCG
Definition: solver_cg.h:98
Differentiation::SD::OptimizerType::lambda
@ lambda
StandardExceptions::ExcNotImplemented
static ::ExceptionBase & ExcNotImplemented()
mu
Definition: function_parser.h:32
parallel::distributed::GridRefinement::refine_and_coarsen_fixed_number
void refine_and_coarsen_fixed_number(parallel::distributed::Triangulation< dim, spacedim > &tria, const ::Vector< Number > &criteria, const double top_fraction_of_cells, const double bottom_fraction_of_cells, const types::global_cell_index max_n_cells=std::numeric_limits< types::global_cell_index >::max())
Definition: grid_refinement.cc:432
GridRefinement::refine_and_coarsen_fixed_number
void refine_and_coarsen_fixed_number(Triangulation< dim, spacedim > &triangulation, const Vector< Number > &criteria, const double top_fraction_of_cells, const double bottom_fraction_of_cells, const unsigned int max_n_cells=std::numeric_limits< unsigned int >::max())
Definition: grid_refinement.cc:189
Functions::ConstantFunction
Definition: function.h:412
Patterns::Tools::to_string
std::string to_string(const T &t)
Definition: patterns.h:2360
Physics::Elasticity::Kinematics::e
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
DoFTools::always
@ always
Definition: dof_tools.h:236
DataOutInterface::write_vtk
void write_vtk(std::ostream &out) const
Definition: data_out_base.cc:6853
LocalIntegrators::Advection::cell_matrix
void cell_matrix(FullMatrix< double > &M, const FEValuesBase< dim > &fe, const FEValuesBase< dim > &fetest, const ArrayView< const std::vector< double >> &velocity, const double factor=1.)
Definition: advection.h:80
update_values
@ update_values
Shape function values.
Definition: fe_update_flags.h:78
second
Point< 2 > second
Definition: grid_out.cc:4353
DataOut::build_patches
virtual void build_patches(const unsigned int n_subdivisions=0)
Definition: data_out.cc:1071
Physics::Elasticity::Kinematics::d
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
internal::p4est::functions
int(&) functions(const void *v1, const void *v2)
Definition: p4est_wrappers.cc:339
OpenCASCADE::point
Point< spacedim > point(const gp_Pnt &p, const double tolerance=1e-10)
Definition: utilities.cc:188
FEValues
Definition: fe.h:38
KellyErrorEstimator
Definition: error_estimator.h:262
WorkStream::run
void run(const std::vector< std::vector< Iterator >> &colored_iterators, Worker worker, Copier copier, const ScratchData &sample_scratch_data, const CopyData &sample_copy_data, const unsigned int queue_length=2 *MultithreadInfo::n_threads(), const unsigned int chunk_size=8)
Definition: work_stream.h:1185
PreconditionSSOR::initialize
void initialize(const MatrixType &A, const typename BaseClass::AdditionalData &parameters=typename BaseClass::AdditionalData())
LAPACKSupport::one
static const types::blas_int one
Definition: lapack_support.h:183
LAPACKSupport::symmetric
@ symmetric
Matrix is symmetric.
Definition: lapack_support.h:115
Algorithms::Events::initial
const Event initial
Definition: event.cc:65
deal_II_exceptions::internals::abort
void abort(const ExceptionBase &exc) noexcept
Definition: exceptions.cc:408
TrilinosWrappers::internal::begin
VectorType::value_type * begin(VectorType &V)
Definition: trilinos_sparse_matrix.cc:51
update_gradients
@ update_gradients
Shape function gradients.
Definition: fe_update_flags.h:84
KellyErrorEstimator::estimate
static void estimate(const Mapping< dim, spacedim > &mapping, const DoFHandlerType &dof, const Quadrature< dim - 1 > &quadrature, const std::map< types::boundary_id, const Function< spacedim, typename InputVector::value_type > * > &neumann_bc, const InputVector &solution, Vector< float > &error, const ComponentMask &component_mask=ComponentMask(), const Function< spacedim > *coefficients=nullptr, const unsigned int n_threads=numbers::invalid_unsigned_int, const types::subdomain_id subdomain_id=numbers::invalid_subdomain_id, const types::material_id material_id=numbers::invalid_material_id, const Strategy strategy=cell_diameter_over_24)
LAPACKSupport::matrix
@ matrix
Contents is actually a matrix.
Definition: lapack_support.h:60
Threads::internal::call
void call(const std::function< RT()> &function, internal::return_value< RT > &ret_val)
Definition: thread_management.h:607
LAPACKSupport::general
@ general
No special properties.
Definition: lapack_support.h:113
numbers
Definition: numbers.h:207
QGauss
Definition: quadrature_lib.h:40
SIMDComparison::equal
@ equal
DataOut_DoFData::attach_dof_handler
void attach_dof_handler(const DoFHandlerType &)
value
static const bool value
Definition: dof_tools_constraints.cc:433
update_JxW_values
@ update_JxW_values
Transformed quadrature weights.
Definition: fe_update_flags.h:129
Assert
#define Assert(cond, exc)
Definition: exceptions.h:1419
internal::assemble
void assemble(const MeshWorker::DoFInfoBox< dim, DOFINFO > &dinfo, A *assembler)
Definition: loop.h:71
GridGenerator::hyper_cube
void hyper_cube(Triangulation< dim, spacedim > &tria, const double left=0., const double right=1., const bool colorize=false)
LAPACKSupport::zero
static const types::blas_int zero
Definition: lapack_support.h:179
PreconditionSSOR
Definition: precondition.h:665
FullMatrix< double >
triangulation
const typename ::parallel::distributed::Triangulation< dim, spacedim > * triangulation
Definition: p4est_wrappers.cc:69
SolverControl
Definition: solver_control.h:67
MeshWorker::loop
void loop(ITERATOR begin, typename identity< ITERATOR >::type end, DOFINFO &dinfo, INFOBOX &info, const std::function< void(DOFINFO &, typename INFOBOX::CellInfo &)> &cell_worker, const std::function< void(DOFINFO &, typename INFOBOX::CellInfo &)> &boundary_worker, const std::function< void(DOFINFO &, DOFINFO &, typename INFOBOX::CellInfo &, typename INFOBOX::CellInfo &)> &face_worker, ASSEMBLER &assembler, const LoopControl &lctrl=LoopControl())
Definition: loop.h:443
first
Point< 2 > first
Definition: grid_out.cc:4352
DataOut
Definition: data_out.h:148
Vector< double >
GridRefinement::refine
void refine(Triangulation< dim, spacedim > &tria, const Vector< Number > &criteria, const double threshold, const unsigned int max_to_mark=numbers::invalid_unsigned_int)
Definition: grid_refinement.cc:41
int
DataOut_DoFData::add_data_vector
void add_data_vector(const VectorType &data, const std::vector< std::string > &names, const DataVectorType type=type_automatic, const std::vector< DataComponentInterpretation::DataComponentInterpretation > &data_component_interpretation=std::vector< DataComponentInterpretation::DataComponentInterpretation >())
Definition: data_out_dof_data.h:1090