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