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-38.h
Go to the documentation of this file.
1) const
494 * {
495 * return (-8. * p(0) * p(1));
496 * }
497 *
498 *
499 * template <>
500 * double RightHandSide<3>::value(const Point<3> &p,
501 * const unsigned int /*component*/) const
502 * {
503 * using numbers::PI;
504 *
506 *
507 * hessian[0][0] = -PI * PI * sin(PI * p(0)) * cos(PI * p(1)) * exp(p(2));
508 * hessian[1][1] = -PI * PI * sin(PI * p(0)) * cos(PI * p(1)) * exp(p(2));
509 * hessian[2][2] = sin(PI * p(0)) * cos(PI * p(1)) * exp(p(2));
510 *
511 * hessian[0][1] = -PI * PI * cos(PI * p(0)) * sin(PI * p(1)) * exp(p(2));
512 * hessian[1][0] = -PI * PI * cos(PI * p(0)) * sin(PI * p(1)) * exp(p(2));
513 *
514 * hessian[0][2] = PI * cos(PI * p(0)) * cos(PI * p(1)) * exp(p(2));
515 * hessian[2][0] = PI * cos(PI * p(0)) * cos(PI * p(1)) * exp(p(2));
516 *
517 * hessian[1][2] = -PI * sin(PI * p(0)) * sin(PI * p(1)) * exp(p(2));
518 * hessian[2][1] = -PI * sin(PI * p(0)) * sin(PI * p(1)) * exp(p(2));
519 *
521 * gradient[0] = PI * cos(PI * p(0)) * cos(PI * p(1)) * exp(p(2));
522 * gradient[1] = -PI * sin(PI * p(0)) * sin(PI * p(1)) * exp(p(2));
523 * gradient[2] = sin(PI * p(0)) * cos(PI * p(1)) * exp(p(2));
524 *
525 * Point<3> normal = p;
526 * normal /= p.norm();
527 *
528 * return (-trace(hessian) + 2 * (gradient * normal) +
529 * (hessian * normal) * normal);
530 * }
531 *
532 *
533 * @endcode
534 *
535 *
536 * <a name="ImplementationofthecodeLaplaceBeltramiProblemcodeclass"></a>
537 * <h3>Implementation of the <code>LaplaceBeltramiProblem</code> class</h3>
538 *
539
540 *
541 * The rest of the program is actually quite unspectacular if you know
542 * @ref step_4 "step-4". Our first step is to define the constructor, setting the
543 * polynomial degree of the finite element and mapping, and associating the
544 * DoF handler to the triangulation:
545 *
546 * @code
547 * template <int spacedim>
548 * LaplaceBeltramiProblem<spacedim>::LaplaceBeltramiProblem(
549 * const unsigned degree)
550 * : fe(degree)
551 * , dof_handler(triangulation)
552 * , mapping(degree)
553 * {}
554 *
555 *
556 * @endcode
557 *
558 *
559 * <a name="LaplaceBeltramiProblemmake_grid_and_dofs"></a>
560 * <h4>LaplaceBeltramiProblem::make_grid_and_dofs</h4>
561 *
562
563 *
564 * The next step is to create the mesh, distribute degrees of freedom, and
565 * set up the various variables that describe the linear system. All of
566 * these steps are standard with the exception of how to create a mesh that
567 * describes a surface. We could generate a mesh for the domain we are
568 * interested in, generate a triangulation using a mesh generator, and read
569 * it in using the GridIn class. Or, as we do here, we generate the mesh
570 * using the facilities in the GridGenerator namespace.
571 *
572
573 *
574 * In particular, what we're going to do is this (enclosed between the set
575 * of braces below): we generate a <code>spacedim</code> dimensional mesh
576 * for the half disk (in 2d) or half ball (in 3d), using the
577 * GridGenerator::half_hyper_ball function. This function sets the boundary
578 * indicators of all faces on the outside of the boundary to zero for the
579 * ones located on the perimeter of the disk/ball, and one on the straight
580 * part that splits the full disk/ball into two halves. The next step is the
581 * main point: The GridGenerator::extract_boundary_mesh function creates a
582 * mesh that consists of those cells that are the faces of the previous mesh,
583 * i.e. it describes the <i>surface</i> cells of the original (volume)
584 * mesh. However, we do not want all faces: only those on the perimeter of
585 * the disk or ball which carry boundary indicator zero; we can select these
586 * cells using a set of boundary indicators that we pass to
587 * GridGenerator::extract_boundary_mesh.
588 *
589
590 *
591 * There is one point that needs to be mentioned. In order to refine a
592 * surface mesh appropriately if the manifold is curved (similarly to
593 * refining the faces of cells that are adjacent to a curved boundary), the
594 * triangulation has to have an object attached to it that describes where
595 * new vertices should be located. If you don't attach such a boundary
596 * object, they will be located halfway between existing vertices; this is
597 * appropriate if you have a domain with straight boundaries (e.g. a
598 * polygon) but not when, as here, the manifold has curvature. So for things
599 * to work properly, we need to attach a manifold object to our (surface)
600 * triangulation, in much the same way as we've already done in 1d for the
601 * boundary. We create such an object and attach it to the triangulation.
602 *
603
604 *
605 * The final step in creating the mesh is to refine it a number of
606 * times. The rest of the function is the same as in previous tutorial
607 * programs.
608 *
609 * @code
610 * template <int spacedim>
611 * void LaplaceBeltramiProblem<spacedim>::make_grid_and_dofs()
612 * {
613 * {
614 * Triangulation<spacedim> volume_mesh;
615 * GridGenerator::half_hyper_ball(volume_mesh);
616 *
617 * std::set<types::boundary_id> boundary_ids;
618 * boundary_ids.insert(0);
619 *
620 * GridGenerator::extract_boundary_mesh(volume_mesh,
621 * triangulation,
622 * boundary_ids);
623 * }
624 * triangulation.set_all_manifold_ids(0);
625 * triangulation.set_manifold(0, SphericalManifold<dim, spacedim>());
626 *
627 * triangulation.refine_global(4);
628 *
629 * std::cout << "Surface mesh has " << triangulation.n_active_cells()
630 * << " cells." << std::endl;
631 *
632 * dof_handler.distribute_dofs(fe);
633 *
634 * std::cout << "Surface mesh has " << dof_handler.n_dofs()
635 * << " degrees of freedom." << std::endl;
636 *
637 * DynamicSparsityPattern dsp(dof_handler.n_dofs(), dof_handler.n_dofs());
638 * DoFTools::make_sparsity_pattern(dof_handler, dsp);
639 * sparsity_pattern.copy_from(dsp);
640 *
641 * system_matrix.reinit(sparsity_pattern);
642 *
643 * solution.reinit(dof_handler.n_dofs());
644 * system_rhs.reinit(dof_handler.n_dofs());
645 * }
646 *
647 *
648 * @endcode
649 *
650 *
651 * <a name="LaplaceBeltramiProblemassemble_system"></a>
652 * <h4>LaplaceBeltramiProblem::assemble_system</h4>
653 *
654
655 *
656 * The following is the central function of this program, assembling the
657 * matrix that corresponds to the surface Laplacian (Laplace-Beltrami
658 * operator). Maybe surprisingly, it actually looks exactly the same as for
659 * the regular Laplace operator discussed in, for example, @ref step_4 "step-4". The key
660 * is that the FEValues::shape_grad() function does the magic: It returns
661 * the surface gradient @f$\nabla_K \phi_i(x_q)@f$ of the @f$i@f$th shape function
662 * at the @f$q@f$th quadrature point. The rest then does not need any changes
663 * either:
664 *
665 * @code
666 * template <int spacedim>
667 * void LaplaceBeltramiProblem<spacedim>::assemble_system()
668 * {
669 * system_matrix = 0;
670 * system_rhs = 0;
671 *
672 * const QGauss<dim> quadrature_formula(2 * fe.degree);
673 * FEValues<dim, spacedim> fe_values(mapping,
674 * fe,
675 * quadrature_formula,
676 * update_values | update_gradients |
677 * update_quadrature_points |
678 * update_JxW_values);
679 *
680 * const unsigned int dofs_per_cell = fe.n_dofs_per_cell();
681 * const unsigned int n_q_points = quadrature_formula.size();
682 *
683 * FullMatrix<double> cell_matrix(dofs_per_cell, dofs_per_cell);
684 * Vector<double> cell_rhs(dofs_per_cell);
685 *
686 * std::vector<double> rhs_values(n_q_points);
687 * std::vector<types::global_dof_index> local_dof_indices(dofs_per_cell);
688 *
689 * RightHandSide<spacedim> rhs;
690 *
691 * for (const auto &cell : dof_handler.active_cell_iterators())
692 * {
693 * cell_matrix = 0;
694 * cell_rhs = 0;
695 *
696 * fe_values.reinit(cell);
697 *
698 * rhs.value_list(fe_values.get_quadrature_points(), rhs_values);
699 *
700 * for (unsigned int i = 0; i < dofs_per_cell; ++i)
701 * for (unsigned int j = 0; j < dofs_per_cell; ++j)
702 * for (unsigned int q_point = 0; q_point < n_q_points; ++q_point)
703 * cell_matrix(i, j) += fe_values.shape_grad(i, q_point) *
704 * fe_values.shape_grad(j, q_point) *
705 * fe_values.JxW(q_point);
706 *
707 * for (unsigned int i = 0; i < dofs_per_cell; ++i)
708 * for (unsigned int q_point = 0; q_point < n_q_points; ++q_point)
709 * cell_rhs(i) += fe_values.shape_value(i, q_point) *
710 * rhs_values[q_point] * fe_values.JxW(q_point);
711 *
712 * cell->get_dof_indices(local_dof_indices);
713 * for (unsigned int i = 0; i < dofs_per_cell; ++i)
714 * {
715 * for (unsigned int j = 0; j < dofs_per_cell; ++j)
716 * system_matrix.add(local_dof_indices[i],
717 * local_dof_indices[j],
718 * cell_matrix(i, j));
719 *
720 * system_rhs(local_dof_indices[i]) += cell_rhs(i);
721 * }
722 * }
723 *
724 * std::map<types::global_dof_index, double> boundary_values;
725 * VectorTools::interpolate_boundary_values(
726 * mapping, dof_handler, 0, Solution<spacedim>(), boundary_values);
727 *
728 * MatrixTools::apply_boundary_values(
729 * boundary_values, system_matrix, solution, system_rhs, false);
730 * }
731 *
732 *
733 *
734 * @endcode
735 *
736 *
737 * <a name="LaplaceBeltramiProblemsolve"></a>
738 * <h4>LaplaceBeltramiProblem::solve</h4>
739 *
740
741 *
742 * The next function is the one that solves the linear system. Here, too, no
743 * changes are necessary:
744 *
745 * @code
746 * template <int spacedim>
747 * void LaplaceBeltramiProblem<spacedim>::solve()
748 * {
749 * SolverControl solver_control(solution.size(), 1e-7 * system_rhs.l2_norm());
750 * SolverCG<Vector<double>> cg(solver_control);
751 *
752 * PreconditionSSOR<SparseMatrix<double>> preconditioner;
753 * preconditioner.initialize(system_matrix, 1.2);
754 *
755 * cg.solve(system_matrix, solution, system_rhs, preconditioner);
756 * }
757 *
758 *
759 *
760 * @endcode
761 *
762 *
763 * <a name="LaplaceBeltramiProblemoutput_result"></a>
764 * <h4>LaplaceBeltramiProblem::output_result</h4>
765 *
766
767 *
768 * This is the function that generates graphical output from the
769 * solution. Most of it is boilerplate code, but there are two points worth
770 * pointing out:
771 *
772
773 *
774 * - The DataOut::add_data_vector() function can take two kinds of vectors:
775 * Either vectors that have one value per degree of freedom defined by the
776 * DoFHandler object previously attached via DataOut::attach_dof_handler();
777 * and vectors that have one value for each cell of the triangulation, for
778 * example to output estimated errors for each cell. Typically, the
779 * DataOut class knows to tell these two kinds of vectors apart: there are
780 * almost always more degrees of freedom than cells, so we can
781 * differentiate by the two kinds looking at the length of a vector. We
782 * could do the same here, but only because we got lucky: we use a half
783 * sphere. If we had used the whole sphere as domain and @f$Q_1@f$ elements,
784 * we would have the same number of cells as vertices and consequently the
785 * two kinds of vectors would have the same number of elements. To avoid
786 * the resulting confusion, we have to tell the DataOut::add_data_vector()
787 * function which kind of vector we have: DoF data. This is what the third
788 * argument to the function does.
789 * - The DataOut::build_patches() function can generate output that subdivides
790 * each cell so that visualization programs can resolve curved manifolds
791 * or higher polynomial degree shape functions better. We here subdivide
792 * each element in each coordinate direction as many times as the
793 * polynomial degree of the finite element in use.
794 *
795 * @code
796 * template <int spacedim>
797 * void LaplaceBeltramiProblem<spacedim>::output_results() const
798 * {
799 * DataOut<dim, DoFHandler<dim, spacedim>> data_out;
800 * data_out.attach_dof_handler(dof_handler);
801 * data_out.add_data_vector(
802 * solution,
803 * "solution",
804 * DataOut<dim, DoFHandler<dim, spacedim>>::type_dof_data);
805 * data_out.build_patches(mapping, mapping.get_degree());
806 *
807 * const std::string filename =
808 * "solution-" + std::to_string(spacedim) + "d.vtk";
809 * std::ofstream output(filename);
810 * data_out.write_vtk(output);
811 * }
812 *
813 *
814 *
815 * @endcode
816 *
817 *
818 * <a name="LaplaceBeltramiProblemcompute_error"></a>
819 * <h4>LaplaceBeltramiProblem::compute_error</h4>
820 *
821
822 *
823 * This is the last piece of functionality: we want to compute the error in
824 * the numerical solution. It is a verbatim copy of the code previously
825 * shown and discussed in @ref step_7 "step-7". As mentioned in the introduction, the
826 * <code>Solution</code> class provides the (tangential) gradient of the
827 * solution. To avoid evaluating the error only a superconvergence points,
828 * we choose a quadrature rule of sufficiently high order.
829 *
830 * @code
831 * template <int spacedim>
832 * void LaplaceBeltramiProblem<spacedim>::compute_error() const
833 * {
834 * Vector<float> difference_per_cell(triangulation.n_active_cells());
835 * VectorTools::integrate_difference(mapping,
836 * dof_handler,
837 * solution,
838 * Solution<spacedim>(),
839 * difference_per_cell,
840 * QGauss<dim>(2 * fe.degree + 1),
841 * VectorTools::H1_norm);
842 *
843 * double h1_error = VectorTools::compute_global_error(triangulation,
844 * difference_per_cell,
845 * VectorTools::H1_norm);
846 * std::cout << "H1 error = " << h1_error << std::endl;
847 * }
848 *
849 *
850 *
851 * @endcode
852 *
853 *
854 * <a name="LaplaceBeltramiProblemrun"></a>
855 * <h4>LaplaceBeltramiProblem::run</h4>
856 *
857
858 *
859 * The last function provides the top-level logic. Its contents are
860 * self-explanatory:
861 *
862 * @code
863 * template <int spacedim>
864 * void LaplaceBeltramiProblem<spacedim>::run()
865 * {
866 * make_grid_and_dofs();
867 * assemble_system();
868 * solve();
869 * output_results();
870 * compute_error();
871 * }
872 * } // namespace Step38
873 *
874 *
875 * @endcode
876 *
877 *
878 * <a name="Themainfunction"></a>
879 * <h3>The main() function</h3>
880 *
881
882 *
883 * The remainder of the program is taken up by the <code>main()</code>
884 * function. It follows exactly the general layout first introduced in @ref step_6 "step-6"
885 * and used in all following tutorial programs:
886 *
887 * @code
888 * int main()
889 * {
890 * try
891 * {
892 * using namespace Step38;
893 *
894 * LaplaceBeltramiProblem<3> laplace_beltrami;
895 * laplace_beltrami.run();
896 * }
897 * catch (std::exception &exc)
898 * {
899 * std::cerr << std::endl
900 * << std::endl
901 * << "----------------------------------------------------"
902 * << std::endl;
903 * std::cerr << "Exception on processing: " << std::endl
904 * << exc.what() << std::endl
905 * << "Aborting!" << std::endl
906 * << "----------------------------------------------------"
907 * << std::endl;
908 * return 1;
909 * }
910 * catch (...)
911 * {
912 * std::cerr << std::endl
913 * << std::endl
914 * << "----------------------------------------------------"
915 * << std::endl;
916 * std::cerr << "Unknown exception!" << std::endl
917 * << "Aborting!" << std::endl
918 * << "----------------------------------------------------"
919 * << std::endl;
920 * return 1;
921 * }
922 *
923 * return 0;
924 * }
925 * @endcode
926<a name="Results"></a><h1>Results</h1>
927
928
929When you run the program, the following output should be printed on screen:
930
931@verbatim
932Surface mesh has 1280 cells.
933Surface mesh has 5185 degrees of freedom.
934H1 error = 0.0217136
935@endverbatim
936
937
938By playing around with the number of global refinements in the
939<code>LaplaceBeltrami::make_grid_and_dofs</code> function you increase or decrease mesh
940refinement. For example, doing one more refinement and only running the 3d surface
941problem yields the following
942output:
943
944@verbatim
945Surface mesh has 5120 cells.
946Surface mesh has 20609 degrees of freedom.
947H1 error = 0.00543481
948@endverbatim
949
950This is what we expect: make the mesh size smaller by a factor of two and the
951error goes down by a factor of four (remember that we use bi-quadratic
952elements). The full sequence of errors from one to five refinements looks like
953this, neatly following the theoretically predicted pattern:
954@verbatim
9550.339438
9560.0864385
9570.0217136
9580.00543481
9590.00135913
960@endverbatim
961
962Finally, the program produces graphical output that we can visualize. Here is
963a plot of the results:
964
965<img src="https://www.dealii.org/images/steps/developer/step-38.solution-3d.png" alt="">
966
967The program also works for 1d curves in 2d, not just 2d surfaces in 3d. You
968can test this by changing the template argument in <code>main()</code> like
969so:
970@code
971 LaplaceBeltramiProblem<2> laplace_beltrami;
972@endcode
973The domain is a curve in 2d, and we can visualize the solution by using the
974third dimension (and color) to denote the value of the function @f$u(x)@f$. This
975then looks like so (the white curve is the domain, the colored curve is the
976solution extruded into the third dimension, clearly showing the change in sign
977as the curve moves from one quadrant of the domain into the adjacent one):
978
979<img src="https://www.dealii.org/images/steps/developer/step-38.solution-2d.png" alt="">
980
981
982<a name="extensions"></a>
983<a name="Possibilitiesforextensions"></a><h3>Possibilities for extensions</h3>
984
985
986Computing on surfaces only becomes interesting if the surface is more
987interesting than just a half sphere. To achieve this, deal.II can read
988meshes that describe surfaces through the usual GridIn class. Or, in case you
989have an analytic description, a simple mesh can sometimes be stretched and
990bent into a shape we are interested in.
991
992Let us consider a relatively simple example: we take the half sphere we used
993before, we stretch it by a factor of 10 in the z-direction, and then we jumble
994the x- and y-coordinates a bit. Let's show the computational domain and the
995solution first before we go into details of the implementation below:
996
997<img src="https://www.dealii.org/images/steps/developer/step-38.warp-1.png" alt="">
998
999<img src="https://www.dealii.org/images/steps/developer/step-38.warp-2.png" alt="">
1000
1001The way to produce such a mesh is by using the GridTools::transform()
1002function. It needs a way to transform each individual mesh point to a
1003different position. Let us here use the following, rather simple function
1004(remember: stretch in one direction, jumble in the other two):
1005
1006@code
1007template <int spacedim>
1008Point<spacedim> warp(const Point<spacedim> &p)
1009{
1010 Point<spacedim> q = p;
1011 q[spacedim-1] *= 10;
1012
1013 if (spacedim >= 2)
1014 q[0] += 2*std::sin(q[spacedim-1]);
1015 if (spacedim >= 3)
1016 q[1] += 2*std::cos(q[spacedim-1]);
1017
1018 return q;
1019}
1020@endcode
1021
1022If we followed the <code>LaplaceBeltrami::make_grid_and_dofs</code> function, we would
1023extract the half spherical surface mesh as before, warp it into the shape we
1024want, and refine as often as necessary. This is not quite as simple as we'd
1025like here, though: refining requires that we have an appropriate manifold
1026object attached to the triangulation that describes where new vertices of the
1027mesh should be located upon refinement. I'm sure it's possible to describe
1028this manifold in a not-too-complicated way by simply undoing the
1029transformation above (yielding the spherical surface again), finding the
1030location of a new point on the sphere, and then re-warping the result. But I'm
1031a lazy person, and since doing this is not really the point here, let's just
1032make our lives a bit easier: we'll extract the half sphere, refine it as
1033often as necessary, get rid of the object that describes the manifold since we
1034now no longer need it, and then finally warp the mesh. With the function
1035above, this would look as follows:
1036
1037@code
1038template <int spacedim>
1039void LaplaceBeltrami<spacedim>::make_grid_and_dofs()
1040{
1041 {
1042 Triangulation<spacedim> volume_mesh;
1043 GridGenerator::half_hyper_ball(volume_mesh);
1044
1045 volume_mesh.refine_global(4);
1046
1047 std::set<types::boundary_id> boundary_ids;
1048 boundary_ids.insert(0);
1049
1051 boundary_ids);
1052 GridTools::transform(&warp<spacedim>, triangulation); /* ** */
1053 std::ofstream x("x"), y("y");
1054 GridOut().write_gnuplot(volume_mesh, x);
1056 }
1057
1058 std::cout << "Surface mesh has " << triangulation.n_active_cells()
1059 << " cells."
1060 << std::endl;
1061 ...
1062}
1063@endcode
1064
1065Note that the only essential addition is the line marked with
1066asterisks. It is worth pointing out one other thing here, though: because we
1067detach the manifold description from the surface mesh, whenever we use a
1068mapping object in the rest of the program, it has no curves boundary
1069description to go on any more. Rather, it will have to use the implicit,
1070FlatManifold class that is used on all parts of the domain not
1071explicitly assigned a different manifold object. Consequently, whether we use
1072MappingQ(2), MappingQ(15) or MappingQ1, each cell of our mesh will be mapped
1073using a bilinear approximation.
1074
1075All these drawbacks aside, the resulting pictures are still pretty. The only
1076other differences to what's in @ref step_38 "step-38" is that we changed the right hand side
1077to @f$f(\mathbf x)=\sin x_3@f$ and the boundary values (through the
1078<code>Solution</code> class) to @f$u(\mathbf x)|_{\partial\Omega}=\cos x_3@f$. Of
1079course, we now no longer know the exact solution, so the computation of the
1080error at the end of <code>LaplaceBeltrami::run</code> will yield a meaningless
1081number.
1082 *
1083 *
1084<a name="PlainProg"></a>
1085<h1> The plain program</h1>
1086@include "step-38.cc"
1087*/
void write_gnuplot(const Triangulation< dim, spacedim > &tria, std::ostream &out, const Mapping< dim, spacedim > *mapping=nullptr) const
Definition: grid_out.cc:4572
constexpr Number trace(const SymmetricTensor< 2, dim, Number > &d)
Definition: tensor.h:472
numbers::NumberTraits< Number >::real_type norm() const
void refine_global(const unsigned int times=1)
VectorizedArray< Number, width > exp(const ::VectorizedArray< Number, width > &x)
VectorizedArray< Number, width > cos(const ::VectorizedArray< Number, width > &x)
VectorizedArray< Number, width > sin(const ::VectorizedArray< Number, width > &x)
Point< 3 > vertices[4]
Point< 2 > first
Definition: grid_out.cc:4587
__global__ void set(Number *val, const Number s, const size_type N)
std::map< typename MeshType< dim - 1, spacedim >::cell_iterator, typename MeshType< dim, spacedim >::face_iterator > extract_boundary_mesh(const MeshType< dim, spacedim > &volume_mesh, MeshType< dim - 1, spacedim > &surface_mesh, const std::set< types::boundary_id > &boundary_ids=std::set< types::boundary_id >())
void half_hyper_ball(Triangulation< dim > &tria, const Point< dim > &center=Point< dim >(), const double radius=1.)
void refine(Triangulation< dim, spacedim > &tria, const Vector< Number > &criteria, const double threshold, const unsigned int max_to_mark=numbers::invalid_unsigned_int)
void transform(const Transformation &transformation, Triangulation< dim, spacedim > &triangulation)
static const types::blas_int one
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)
constexpr ReturnType< rank, T >::value_type & extract(T &t, const ArrayType &indices)
VectorType::value_type * end(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
static constexpr double PI
Definition: numbers.h:231
void transform(const InputIterator &begin_in, const InputIterator &end_in, OutputIterator out, const Predicate &predicate, const unsigned int grainsize)
Definition: parallel.h:213
::VectorizedArray< Number, width > cos(const ::VectorizedArray< Number, width > &)
::VectorizedArray< Number, width > sin(const ::VectorizedArray< Number, width > &)
const ::parallel::distributed::Triangulation< dim, spacedim > * triangulation