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-4.h
Go to the documentation of this file.
1 ) const
406  * {
407  * double return_value = 0.0;
408  * for (unsigned int i = 0; i < dim; ++i)
409  * return_value += 4.0 * std::pow(p(i), 4.0);
410  *
411  * return return_value;
412  * }
413  *
414  *
415  * @endcode
416  *
417  * As boundary values, we choose @f$x^2+y^2@f$ in 2D, and @f$x^2+y^2+z^2@f$ in 3D. This
418  * happens to be equal to the square of the vector from the origin to the
419  * point at which we would like to evaluate the function, irrespective of the
420  * dimension. So that is what we return:
421  *
422  * @code
423  * template <int dim>
424  * double BoundaryValues<dim>::value(const Point<dim> &p,
425  * const unsigned int /*component*/) const
426  * {
427  * return p.square();
428  * }
429  *
430  *
431  *
432  * @endcode
433  *
434  *
435  * <a name="ImplementationofthecodeStep4codeclass"></a>
436  * <h3>Implementation of the <code>Step4</code> class</h3>
437  *
438 
439  *
440  * Next for the implementation of the class template that makes use of the
441  * functions above. As before, we will write everything as templates that have
442  * a formal parameter <code>dim</code> that we assume unknown at the time we
443  * define the template functions. Only later, the compiler will find a
444  * declaration of <code>Step4@<2@></code> (in the <code>main</code> function,
445  * actually) and compile the entire class with <code>dim</code> replaced by 2,
446  * a process referred to as `instantiation of a template'. When doing so, it
447  * will also replace instances of <code>RightHandSide@<dim@></code> by
448  * <code>RightHandSide@<2@></code> and instantiate the latter class from the
449  * class template.
450  *
451 
452  *
453  * In fact, the compiler will also find a declaration <code>Step4@<3@></code>
454  * in <code>main()</code>. This will cause it to again go back to the general
455  * <code>Step4@<dim@></code> template, replace all occurrences of
456  * <code>dim</code>, this time by 3, and compile the class a second time. Note
457  * that the two instantiations <code>Step4@<2@></code> and
458  * <code>Step4@<3@></code> are completely independent classes; their only
459  * common feature is that they are both instantiated from the same general
460  * template, but they are not convertible into each other, for example, and
461  * share no code (both instantiations are compiled completely independently).
462  *
463 
464  *
465  *
466 
467  *
468  *
469  * <a name="Step4Step4"></a>
470  * <h4>Step4::Step4</h4>
471  *
472 
473  *
474  * After this introduction, here is the constructor of the <code>Step4</code>
475  * class. It specifies the desired polynomial degree of the finite elements
476  * and associates the DoFHandler to the triangulation just as in the previous
477  * example program, @ref step_3 "step-3":
478  *
479  * @code
480  * template <int dim>
481  * Step4<dim>::Step4()
482  * : fe(1)
483  * , dof_handler(triangulation)
484  * {}
485  *
486  *
487  * @endcode
488  *
489  *
490  * <a name="Step4make_grid"></a>
491  * <h4>Step4::make_grid</h4>
492  *
493 
494  *
495  * Grid creation is something inherently dimension dependent. However, as long
496  * as the domains are sufficiently similar in 2D or 3D, the library can
497  * abstract for you. In our case, we would like to again solve on the square
498  * @f$[-1,1]\times [-1,1]@f$ in 2D, or on the cube @f$[-1,1] \times [-1,1] \times
499  * [-1,1]@f$ in 3D; both can be termed GridGenerator::hyper_cube(), so we may
500  * use the same function in whatever dimension we are. Of course, the
501  * functions that create a hypercube in two and three dimensions are very much
502  * different, but that is something you need not care about. Let the library
503  * handle the difficult things.
504  *
505  * @code
506  * template <int dim>
507  * void Step4<dim>::make_grid()
508  * {
509  * GridGenerator::hyper_cube(triangulation, -1, 1);
510  * triangulation.refine_global(4);
511  *
512  * std::cout << " Number of active cells: " << triangulation.n_active_cells()
513  * << std::endl
514  * << " Total number of cells: " << triangulation.n_cells()
515  * << std::endl;
516  * }
517  *
518  * @endcode
519  *
520  *
521  * <a name="Step4setup_system"></a>
522  * <h4>Step4::setup_system</h4>
523  *
524 
525  *
526  * This function looks exactly like in the previous example, although it
527  * performs actions that in their details are quite different if
528  * <code>dim</code> happens to be 3. The only significant difference from a
529  * user's perspective is the number of cells resulting, which is much higher
530  * in three than in two space dimensions!
531  *
532  * @code
533  * template <int dim>
534  * void Step4<dim>::setup_system()
535  * {
536  * dof_handler.distribute_dofs(fe);
537  *
538  * std::cout << " Number of degrees of freedom: " << dof_handler.n_dofs()
539  * << std::endl;
540  *
541  * DynamicSparsityPattern dsp(dof_handler.n_dofs());
542  * DoFTools::make_sparsity_pattern(dof_handler, dsp);
543  * sparsity_pattern.copy_from(dsp);
544  *
545  * system_matrix.reinit(sparsity_pattern);
546  *
547  * solution.reinit(dof_handler.n_dofs());
548  * system_rhs.reinit(dof_handler.n_dofs());
549  * }
550  *
551  *
552  * @endcode
553  *
554  *
555  * <a name="Step4assemble_system"></a>
556  * <h4>Step4::assemble_system</h4>
557  *
558 
559  *
560  * Unlike in the previous example, we would now like to use a non-constant
561  * right hand side function and non-zero boundary values. Both are tasks that
562  * are readily achieved with only a few new lines of code in the assemblage of
563  * the matrix and right hand side.
564  *
565 
566  *
567  * More interesting, though, is the way we assemble matrix and right hand side
568  * vector dimension independently: there is simply no difference to the
569  * two-dimensional case. Since the important objects used in this function
570  * (quadrature formula, FEValues) depend on the dimension by way of a template
571  * parameter as well, they can take care of setting up properly everything for
572  * the dimension for which this function is compiled. By declaring all classes
573  * which might depend on the dimension using a template parameter, the library
574  * can make nearly all work for you and you don't have to care about most
575  * things.
576  *
577  * @code
578  * template <int dim>
579  * void Step4<dim>::assemble_system()
580  * {
581  * QGauss<dim> quadrature_formula(fe.degree + 1);
582  *
583  * @endcode
584  *
585  * We wanted to have a non-constant right hand side, so we use an object of
586  * the class declared above to generate the necessary data. Since this right
587  * hand side object is only used locally in the present function, we declare
588  * it here as a local variable:
589  *
590  * @code
591  * RightHandSide<dim> right_hand_side;
592  *
593  * @endcode
594  *
595  * Compared to the previous example, in order to evaluate the non-constant
596  * right hand side function we now also need the quadrature points on the
597  * cell we are presently on (previously, we only required values and
598  * gradients of the shape function from the FEValues object, as well as the
599  * quadrature weights, FEValues::JxW() ). We can tell the FEValues object to
600  * do for us by also giving it the #update_quadrature_points flag:
601  *
602  * @code
603  * FEValues<dim> fe_values(fe,
604  * quadrature_formula,
605  * update_values | update_gradients |
606  * update_quadrature_points | update_JxW_values);
607  *
608  * @endcode
609  *
610  * We then again define the same abbreviation as in the previous program.
611  * The value of this variable of course depends on the dimension which we
612  * are presently using, but the FiniteElement class does all the necessary
613  * work for you and you don't have to care about the dimension dependent
614  * parts:
615  *
616  * @code
617  * const unsigned int dofs_per_cell = fe.dofs_per_cell;
618  *
619  * FullMatrix<double> cell_matrix(dofs_per_cell, dofs_per_cell);
620  * Vector<double> cell_rhs(dofs_per_cell);
621  *
622  * std::vector<types::global_dof_index> local_dof_indices(dofs_per_cell);
623  *
624  * @endcode
625  *
626  * Next, we again have to loop over all cells and assemble local
627  * contributions. Note, that a cell is a quadrilateral in two space
628  * dimensions, but a hexahedron in 3D. In fact, the
629  * <code>active_cell_iterator</code> data type is something different,
630  * depending on the dimension we are in, but to the outside world they look
631  * alike and you will probably never see a difference. In any case, the real
632  * type is hidden by using `auto`:
633  *
634  * @code
635  * for (const auto &cell : dof_handler.active_cell_iterators())
636  * {
637  * fe_values.reinit(cell);
638  * cell_matrix = 0;
639  * cell_rhs = 0;
640  *
641  * @endcode
642  *
643  * Now we have to assemble the local matrix and right hand side. This is
644  * done exactly like in the previous example, but now we revert the
645  * order of the loops (which we can safely do since they are independent
646  * of each other) and merge the loops for the local matrix and the local
647  * vector as far as possible to make things a bit faster.
648  *
649 
650  *
651  * Assembling the right hand side presents the only significant
652  * difference to how we did things in @ref step_3 "step-3": Instead of using a
653  * constant right hand side with value 1, we use the object representing
654  * the right hand side and evaluate it at the quadrature points:
655  *
656  * @code
657  * for (const unsigned int q_index : fe_values.quadrature_point_indices())
658  * for (const unsigned int i : fe_values.dof_indices())
659  * {
660  * for (const unsigned int j : fe_values.dof_indices())
661  * cell_matrix(i, j) +=
662  * (fe_values.shape_grad(i, q_index) * // grad phi_i(x_q)
663  * fe_values.shape_grad(j, q_index) * // grad phi_j(x_q)
664  * fe_values.JxW(q_index)); // dx
665  *
666  * const auto x_q = fe_values.quadrature_point(q_index);
667  * cell_rhs(i) += (fe_values.shape_value(i, q_index) * // phi_i(x_q)
668  * right_hand_side.value(x_q) * // f(x_q)
669  * fe_values.JxW(q_index)); // dx
670  * }
671  * @endcode
672  *
673  * As a final remark to these loops: when we assemble the local
674  * contributions into <code>cell_matrix(i,j)</code>, we have to multiply
675  * the gradients of shape functions @f$i@f$ and @f$j@f$ at point number
676  * q_index and
677  * multiply it with the scalar weights JxW. This is what actually
678  * happens: <code>fe_values.shape_grad(i,q_index)</code> returns a
679  * <code>dim</code> dimensional vector, represented by a
680  * <code>Tensor@<1,dim@></code> object, and the operator* that
681  * multiplies it with the result of
682  * <code>fe_values.shape_grad(j,q_index)</code> makes sure that the
683  * <code>dim</code> components of the two vectors are properly
684  * contracted, and the result is a scalar floating point number that
685  * then is multiplied with the weights. Internally, this operator* makes
686  * sure that this happens correctly for all <code>dim</code> components
687  * of the vectors, whether <code>dim</code> be 2, 3, or any other space
688  * dimension; from a user's perspective, this is not something worth
689  * bothering with, however, making things a lot simpler if one wants to
690  * write code dimension independently.
691  *
692 
693  *
694  * With the local systems assembled, the transfer into the global matrix
695  * and right hand side is done exactly as before, but here we have again
696  * merged some loops for efficiency:
697  *
698  * @code
699  * cell->get_dof_indices(local_dof_indices);
700  * for (const unsigned int i : fe_values.dof_indices())
701  * {
702  * for (const unsigned int j : fe_values.dof_indices())
703  * system_matrix.add(local_dof_indices[i],
704  * local_dof_indices[j],
705  * cell_matrix(i, j));
706  *
707  * system_rhs(local_dof_indices[i]) += cell_rhs(i);
708  * }
709  * }
710  *
711  * @endcode
712  *
713  * As the final step in this function, we wanted to have non-homogeneous
714  * boundary values in this example, unlike the one before. This is a simple
715  * task, we only have to replace the Functions::ZeroFunction used there by an
716  * object of the class which describes the boundary values we would like to
717  * use (i.e. the <code>BoundaryValues</code> class declared above):
718  *
719 
720  *
721  * The function VectorTools::interpolate_boundary_values() will only work
722  * on faces that have been marked with boundary indicator 0 (because that's
723  * what we say the function should work on with the second argument below).
724  * If there are faces with boundary id other than 0, then the function
725  * interpolate_boundary_values will do nothing on these faces. For
726  * the Laplace equation doing nothing is equivalent to assuming that
727  * on those parts of the boundary a zero Neumann boundary condition holds.
728  *
729  * @code
730  * std::map<types::global_dof_index, double> boundary_values;
732  * 0,
733  * BoundaryValues<dim>(),
734  * boundary_values);
735  * MatrixTools::apply_boundary_values(boundary_values,
736  * system_matrix,
737  * solution,
738  * system_rhs);
739  * }
740  *
741  *
742  * @endcode
743  *
744  *
745  * <a name="Step4solve"></a>
746  * <h4>Step4::solve</h4>
747  *
748 
749  *
750  * Solving the linear system of equations is something that looks almost
751  * identical in most programs. In particular, it is dimension independent, so
752  * this function is copied verbatim from the previous example.
753  *
754  * @code
755  * template <int dim>
756  * void Step4<dim>::solve()
757  * {
758  * SolverControl solver_control(1000, 1e-12);
759  * SolverCG<Vector<double>> solver(solver_control);
760  * solver.solve(system_matrix, solution, system_rhs, PreconditionIdentity());
761  *
762  * @endcode
763  *
764  * We have made one addition, though: since we suppress output from the
765  * linear solvers, we have to print the number of iterations by hand.
766  *
767  * @code
768  * std::cout << " " << solver_control.last_step()
769  * << " CG iterations needed to obtain convergence." << std::endl;
770  * }
771  *
772  *
773  * @endcode
774  *
775  *
776  * <a name="Step4output_results"></a>
777  * <h4>Step4::output_results</h4>
778  *
779 
780  *
781  * This function also does what the respective one did in @ref step_3 "step-3". No changes
782  * here for dimension independence either.
783  *
784 
785  *
786  * Since the program will run both 2d and 3d versions of the Laplace solver,
787  * we use the dimension in the filename to generate distinct filenames for
788  * each run (in a better program, one would check whether <code>dim</code> can
789  * have other values than 2 or 3, but we neglect this here for the sake of
790  * brevity).
791  *
792  * @code
793  * template <int dim>
794  * void Step4<dim>::output_results() const
795  * {
796  * DataOut<dim> data_out;
797  *
798  * data_out.attach_dof_handler(dof_handler);
799  * data_out.add_data_vector(solution, "solution");
800  *
801  * data_out.build_patches();
802  *
803  * std::ofstream output(dim == 2 ? "solution-2d.vtk" : "solution-3d.vtk");
804  * data_out.write_vtk(output);
805  * }
806  *
807  *
808  *
809  * @endcode
810  *
811  *
812  * <a name="Step4run"></a>
813  * <h4>Step4::run</h4>
814  *
815 
816  *
817  * This is the function which has the top-level control over everything. Apart
818  * from one line of additional output, it is the same as for the previous
819  * example.
820  *
821  * @code
822  * template <int dim>
823  * void Step4<dim>::run()
824  * {
825  * std::cout << "Solving problem in " << dim << " space dimensions."
826  * << std::endl;
827  *
828  * make_grid();
829  * setup_system();
830  * assemble_system();
831  * solve();
832  * output_results();
833  * }
834  *
835  *
836  * @endcode
837  *
838  *
839  * <a name="Thecodemaincodefunction"></a>
840  * <h3>The <code>main</code> function</h3>
841  *
842 
843  *
844  * And this is the main function. It also looks mostly like in @ref step_3 "step-3", but if
845  * you look at the code below, note how we first create a variable of type
846  * <code>Step4@<2@></code> (forcing the compiler to compile the class template
847  * with <code>dim</code> replaced by <code>2</code>) and run a 2d simulation,
848  * and then we do the whole thing over in 3d.
849  *
850 
851  *
852  * In practice, this is probably not what you would do very frequently (you
853  * probably either want to solve a 2d problem, or one in 3d, but not both at
854  * the same time). However, it demonstrates the mechanism by which we can
855  * simply change which dimension we want in a single place, and thereby force
856  * the compiler to recompile the dimension independent class templates for the
857  * dimension we request. The emphasis here lies on the fact that we only need
858  * to change a single place. This makes it rather trivial to debug the program
859  * in 2d where computations are fast, and then switch a single place to a 3 to
860  * run the much more computing intensive program in 3d for `real'
861  * computations.
862  *
863 
864  *
865  * Each of the two blocks is enclosed in braces to make sure that the
866  * <code>laplace_problem_2d</code> variable goes out of scope (and releases
867  * the memory it holds) before we move on to allocate memory for the 3d
868  * case. Without the additional braces, the <code>laplace_problem_2d</code>
869  * variable would only be destroyed at the end of the function, i.e. after
870  * running the 3d problem, and would needlessly hog memory while the 3d run
871  * could actually use it.
872  *
873  * @code
874  * int main()
875  * {
876  * deallog.depth_console(0);
877  * {
878  * Step4<2> laplace_problem_2d;
879  * laplace_problem_2d.run();
880  * }
881  *
882  * {
883  * Step4<3> laplace_problem_3d;
884  * laplace_problem_3d.run();
885  * }
886  *
887  * return 0;
888  * }
889  * @endcode
890 <a name="Results"></a><h1>Results</h1>
891 
892 
893 
894 The output of the program looks as follows (the number of iterations
895 may vary by one or two, depending on your computer, since this is
896 often dependent on the round-off accuracy of floating point
897 operations, which differs between processors):
898 @code
899 Solving problem in 2 space dimensions.
900  Number of active cells: 256
901  Total number of cells: 341
902  Number of degrees of freedom: 289
903  26 CG iterations needed to obtain convergence.
904 Solving problem in 3 space dimensions.
905  Number of active cells: 4096
906  Total number of cells: 4681
907  Number of degrees of freedom: 4913
908  30 CG iterations needed to obtain convergence.
909 @endcode
910 It is obvious that in three spatial dimensions the number of cells and
911 therefore also the number of degrees of freedom is
912 much higher. What cannot be seen here, is that besides this higher
913 number of rows and columns in the matrix, there are also significantly
914 more entries per row of the matrix in three space
915 dimensions. Together, this leads to a much higher numerical effort for
916 solving the system of equation, which you can feel in the run time of the two
917 solution steps when you actually run the program.
918 
919 
920 
921 The program produces two files: <code>solution-2d.vtk</code> and
922 <code>solution-3d.vtk</code>, which can be viewed using the programs
923 VisIt or Paraview (in case you do not have these programs, you can easily
924 change the
925 output format in the program to something which you can view more
926 easily). Visualizing solutions is a bit of an art, but it can also be fun, so
927 you should play around with your favorite visualization tool to get familiar
928 with its functionality. Here's what I have come up with for the 2d solution:
929 
930 <p align="center">
931  <img src="https://www.dealii.org/images/steps/developer/step-4.solution-2d.png" alt="">
932 </p>
933 
934 (See also <a href="http://www.math.colostate.edu/~bangerth/videos.676.11.html">video lecture 11</a>, <a href="http://www.math.colostate.edu/~bangerth/videos.676.32.html">video lecture 32</a>.)
935 The picture shows the solution of the problem under consideration as
936 a 3D plot. As can be seen, the solution is almost flat in the interior
937 of the domain and has a higher curvature near the boundary. This, of
938 course, is due to the fact that for Laplace's equation the curvature
939 of the solution is equal to the right hand side and that was chosen as
940 a quartic polynomial which is nearly zero in the interior and is only
941 rising sharply when approaching the boundaries of the domain; the
942 maximal values of the right hand side function are at the corners of
943 the domain, where also the solution is moving most rapidly.
944 It is also nice to see that the solution follows the desired quadratic
945 boundary values along the boundaries of the domain.
946 It can also be useful to verify a computed solution against an analytical
947 solution. For an explanation of this technique, see @ref step_7 "step-7".
948 
949 On the other hand, even though the picture does not show the mesh lines
950 explicitly, you can see them as little kinks in the solution. This clearly
951 indicates that the solution hasn't been computed to very high accuracy and
952 that to get a better solution, we may have to compute on a finer mesh.
953 
954 In three spatial dimensions, visualization is a bit more difficult. The left
955 picture shows the solution and the mesh it was computed on on the surface of
956 the domain. This is nice, but it has the drawback that it completely hides
957 what is happening on the inside. The picture on the right is an attempt at
958 visualizing the interior as well, by showing surfaces where the solution has
959 constant values (as indicated by the legend at the top left). Isosurface
960 pictures look best if one makes the individual surfaces slightly transparent
961 so that it is possible to see through them and see what's behind.
962 
963 <table width="60%" align="center">
964  <tr>
965  <td align="center">
966  <img src="https://www.dealii.org/images/steps/developer/step-4.solution-3d.png" alt="">
967  </td>
968  <td align="center">
969  <img src="https://www.dealii.org/images/steps/developer/step-4.contours-3d.png" alt="">
970  </td>
971  </tr>
972 </table>
973 
974 @note
975 A final remark on visualization: the idea of visualization is to give insight,
976 which is not the same as displaying information. In particular, it is easy to
977 overload a picture with information, but while it shows more information it
978 makes it also more difficult to glean insight. As an example, the program I
979 used to generate these pictures, VisIt, by default puts tick marks on every
980 axis, puts a big fat label "X Axis" on the @f$x@f$ axis and similar for the other
981 axes, shows the file name from which the data was taken in the top left and
982 the name of the user doing so and the time and date on the bottom right. None
983 of this is important
984 here: the axes are equally easy to make out because the tripod at the bottom
985 left is still visible, and we know from the program that the domain is
986 @f$[-1,1]^3@f$, so there is no need for tick marks. As a consequence, I have
987 switched off all the extraneous stuff in the picture: the art of visualization
988 is to reduce the picture to those parts that are important to see what one
989 wants to see, but no more.
990 
991 
992 
993 <a name="extensions"></a>
994 <a name="Possibilitiesforextensions"></a><h3>Possibilities for extensions</h3>
995 
996 
997 
998 Essentially the possibilities for playing around with the program are the same
999 as for the previous one, except that they will now also apply to the 3d
1000 case. For inspiration read up on <a href="step_3.html#extensions"
1001 target="body">possible extensions in the documentation of step 3</a>.
1002 
1003  *
1004  *
1005 <a name="PlainProg"></a>
1006 <h1> The plain program</h1>
1007 @include "step-4.cc"
1008 */
SolverCG
Definition: solver_cg.h:98
internal::SymmetricTensorAccessors::merge
constexpr TableIndices< 2 > merge(const TableIndices< 2 > &previous_indices, const unsigned int new_index, const unsigned int position)
Definition: symmetric_tensor.h:146
MatrixTools::apply_boundary_values
void apply_boundary_values(const std::map< types::global_dof_index, number > &boundary_values, SparseMatrix< number > &matrix, Vector< number > &solution, Vector< number > &right_hand_side, const bool eliminate_columns=true)
Definition: matrix_tools.cc:81
Physics::Elasticity::Kinematics::e
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
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
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
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
level
unsigned int level
Definition: grid_out.cc:4355
PreconditionIdentity
Definition: precondition.h:80
LAPACKSupport::one
static const types::blas_int one
Definition: lapack_support.h:183
DoFTools::make_sparsity_pattern
void make_sparsity_pattern(const DoFHandlerType &dof_handler, SparsityPatternType &sparsity_pattern, const AffineConstraints< number > &constraints=AffineConstraints< number >(), const bool keep_constrained_dofs=true, const types::subdomain_id subdomain_id=numbers::invalid_subdomain_id)
Definition: dof_tools_sparsity.cc:63
VectorTools::interpolate_boundary_values
void interpolate_boundary_values(const Mapping< dim, spacedim > &mapping, const DoFHandlerType< dim, spacedim > &dof, const std::map< types::boundary_id, const Function< spacedim, number > * > &function_map, std::map< types::global_dof_index, number > &boundary_values, const ComponentMask &component_mask=ComponentMask())
Tensor
Definition: tensor.h:450
LAPACKSupport::matrix
@ matrix
Contents is actually a matrix.
Definition: lapack_support.h:60
std_cxx17::apply
auto apply(F &&fn, Tuple &&t) -> decltype(apply_impl(std::forward< F >(fn), std::forward< Tuple >(t), std_cxx14::make_index_sequence< std::tuple_size< typename std::remove_reference< Tuple >::type >::value >()))
Definition: tuple.h:40
DynamicSparsityPattern
Definition: dynamic_sparsity_pattern.h:323
SIMDComparison::equal
@ equal
LAPACKSupport::A
static const char A
Definition: lapack_support.h:155
DataOut_DoFData::attach_dof_handler
void attach_dof_handler(const DoFHandlerType &)
value
static const bool value
Definition: dof_tools_constraints.cc:433
internal::assemble
void assemble(const MeshWorker::DoFInfoBox< dim, DOFINFO > &dinfo, A *assembler)
Definition: loop.h:71
LAPACKSupport::zero
static const types::blas_int zero
Definition: lapack_support.h:179
Point< dim >
FullMatrix< double >
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< dim >
Vector< double >
center
Point< 3 > center
Definition: data_out_base.cc:169
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