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
parallel_in_time.h
Go to the documentation of this file.
1
337 * double tstop; /* evolve to this time*/
338 * int level;
339 * double deltaT;
340 *
341 * int index;
342 * braid_StepStatusGetLevel(status, &level);
343 * braid_StepStatusGetTstartTstop(status, &tstart, &tstop);
344 * braid_StepStatusGetTIndex(status, &index);
345 *
346 * deltaT = tstop - tstart;
347 *
348 * ::Vector<double>& solution = u->data;
349 *
350 * HeatEquation<2>& heateq = app->eq;
351 *
352 * heateq.step(solution, deltaT, tstart, index);
353 *
354 * return 0;
355 * }
356 *
357 *
358 * @endcode
359 *
360 * In this function we initialize a vector at an arbitrary time.
361 * At this point we don't know anything about what the solution
362 * looks like, and we can really initialize to anything, so in
363 * this case use reinit to initialize the memory and set the
364 * values to zero.
365 *
366 * @code
367 * int
368 * my_Init(braid_App app,
369 * double t,
370 * braid_Vector *u_ptr)
371 * {
372 * my_Vector *u = new(my_Vector);
373 * int size = app->eq.size();
374 * u->data.reinit(size);
375 *
376 * app->eq.initialize(t, u->data);
377 *
378 * *u_ptr = u;
379 *
380 * return 0;
381 * }
382 *
383 * @endcode
384 *
385 * Here we need to copy the vector u into the vector v. We do this
386 * by allocating a new vector, then reinitializing the deal.ii
387 * vector to the correct size. The deal.ii reinitialization sets
388 * every value to zero, so next we need to iterate over the vector
389 * u and copy the values to the new vector v.
390 *
391 * @code
392 * int
393 * my_Clone(braid_App app,
394 * braid_Vector u,
395 * braid_Vector *v_ptr)
396 * {
397 * UNUSED(app);
398 * my_Vector *v = new(my_Vector);
399 * int size = u->data.size();
400 * v->data.reinit(size);
401 * for(size_t i=0, end=v->data.size(); i != end; ++i)
402 * {
403 * v->data[i] = u->data[i];
404 * }
405 * *v_ptr = v;
406 *
407 * return 0;
408 * }
409 *
410 * @endcode
411 *
412 * Here we need to free the memory used by vector u. This is
413 * pretty simple since the deal.ii vector is stored inside the
414 * XBraid vector, so we just delete the XBraid vector u and it
415 * puts the deal.ii vector out of scope and releases its memory.
416 *
417 * @code
418 * int
419 * my_Free(braid_App app,
420 * braid_Vector u)
421 * {
422 * UNUSED(app);
423 * delete u;
424 *
425 * return 0;
426 * }
427 *
428 * @endcode
429 *
430 * This is to perform an axpy type operation. That is to say we
431 * do @f$y = \alpha x + \beta y@f$. Fortunately deal.ii already has
432 * this operation built in to its vector class, so we get the
433 * reference to the vector y and call the sadd method.
434 *
435 * @code
436 * int my_Sum(braid_App app,
437 * double alpha,
438 * braid_Vector x,
439 * double beta,
440 * braid_Vector y)
441 * {
442 * UNUSED(app);
443 * Vector<double>& vec = y->data;
444 * vec.sadd(beta, alpha, x->data);
445 *
446 * return 0;
447 * }
448 *
449 * @endcode
450 *
451 * This calculates the spatial norm using the l2 norm. According
452 * to XBraid, this could be just about any spatial norm but we'll
453 * keep it simple and used deal.ii vector's built in l2_norm method.
454 *
455 * @code
456 * int
457 * my_SpatialNorm(braid_App app,
458 * braid_Vector u,
459 * double *norm_ptr)
460 * {
461 * UNUSED(app);
462 * double dot = 0.0;
463 * dot = u->data.l2_norm();
464 * *norm_ptr = dot;
465 *
466 * return 0;
467 * }
468 *
469 * @endcode
470 *
471 * This function is called at various points depending on the access
472 * level specified when configuring the XBraid struct. This function
473 * is used to print out data during the run time, such as plots of the
474 * data. The status struct contains a ton of information about the
475 * simulation run. Here we get the current time and timestep number.
476 * The output_results function is called to plot the solution data.
477 * If the method of manufactured solutions is being used, then the
478 * error of this time step is computed and processed.
479 *
480 * @code
481 * int
482 * my_Access(braid_App app,
483 * braid_Vector u,
484 * braid_AccessStatus astatus)
485 * {
486 * double t;
487 * int index;
488 *
489 * braid_AccessStatusGetT(astatus, &t);
490 * braid_AccessStatusGetTIndex(astatus, &index);
491 *
492 * app->eq.output_results(index, t, u->data);
493 *
494 * #if DO_MFG
495 * if(index == app->final_step)
496 * {
497 * app->eq.process_solution(t, index, u->data);
498 * }
499 * #endif
500 *
501 * return 0;
502 * }
503 *
504 * @endcode
505 *
506 * This calculates the size of buffer needed to pack the solution
507 * data into a linear buffer for transfer to another processor via
508 * MPI. We query the size of the data from the HeatEquation class
509 * and return the buffer size.
510 *
511 * @code
512 * int
513 * my_BufSize(braid_App app,
514 * int *size_ptr,
515 * braid_BufferStatus bstatus)
516 * {
517 * UNUSED(bstatus);
518 * int size = app->eq.size();
519 * *size_ptr = (size+1)*sizeof(double);
520 *
521 * return 0;
522 * }
523 *
524 * @endcode
525 *
526 * This function packs a linear buffer with data so that the buffer
527 * may be sent to another processor via MPI. The buffer is cast to
528 * a type we can work with. The first element of the buffer is the
529 * size of the buffer. Then we iterate over soltuion vector u and
530 * fill the buffer with our solution data. Finally we tell XBraid
531 * how much data we wrote.
532 *
533 * @code
534 * int
535 * my_BufPack(braid_App app,
536 * braid_Vector u,
537 * void *buffer,
538 * braid_BufferStatus bstatus)
539 * {
540 *
541 * UNUSED(app);
542 * double *dbuffer = (double*)buffer;
543 * int size = u->data.size();
544 * dbuffer[0] = size;
545 * for(int i=0; i != size; ++i)
546 * {
547 * dbuffer[i+1] = (u->data)[i];
548 * }
549 * braid_BufferStatusSetSize(bstatus, (size+1)*sizeof(double));
550 *
551 * return 0;
552 * }
553 *
554 * @endcode
555 *
556 * This function unpacks a buffer that was recieved from a different
557 * processor via MPI. The size of the buffer is read from the first
558 * element, then we iterate over the size of the buffer and fill
559 * the values of solution vector u with the data in the buffer.
560 *
561 * @code
562 * int
563 * my_BufUnpack(braid_App app,
564 * void *buffer,
565 * braid_Vector *u_ptr,
566 * braid_BufferStatus bstatus)
567 * {
568 * UNUSED(app);
569 * UNUSED(bstatus);
570 *
571 * my_Vector *u = NULL;
572 * double *dbuffer = (double*)buffer;
573 * int size = static_cast<int>(dbuffer[0]);
574 * u = new(my_Vector);
575 * u->data.reinit(size);
576 *
577 * for(int i = 0; i != size; ++i)
578 * {
579 * (u->data)[i] = dbuffer[i+1];
580 * }
581 * *u_ptr = u;
582 *
583 * return 0;
584 * }
585 * @endcode
586
587
588<a name="ann-src/BraidFuncs.hh"></a>
589<h1>Annotated version of src/BraidFuncs.hh</h1>
590 *
591 *
592 *
593 *
594 * @code
595 * #ifndef _BRAIDFUNCS_H_
596 * #define _BRAIDFUNCS_H_
597 *
598 * /**
599 * * \file BraidFuncs.cc
600 * * \brief Contains the implementation of the mandatory X-Braid functions
601 * *
602 * * X-Braid mandates several functions in order to drive the solution.
603 * * This file contains the implementation of said mandatory functions.
604 * * See the X-Braid documentation for more information.
605 * * There are several functions that are optional in X-Braid that may
606 * * or may not be implemented in here.
607 * *
608 * */
609 *
610 *
611 * /*-------- Third Party --------*/
612 * #include <deal.II/numerics/vector_tools.h>
613 *
614 * #include <braid.h>
615 * #include <braid_test.h>
616 *
617 * /*-------- Project --------*/
618 * #include "HeatEquation.hh"
619 *
620 * @endcode
621 *
622 * This struct contains all data that changes with time. For now
623 * this is just the solution data. When doing AMR this should
624 * probably include the triangulization, the sparsity patter,
625 * constraints, etc.
626 *
627 * @code
628 * /**
629 * * \brief Struct that contains the deal.ii vector.
630 * */
631 * typedef struct _braid_Vector_struct
632 * {
633 * ::Vector<double> data;
634 * } my_Vector;
635 *
636 * @endcode
637 *
638 * This struct contains all the data that is unchanging with time.
639 *
640 * @code
641 * /**
642 * * \brief Struct that contains the HeatEquation and final
643 * * time step number.
644 * */
645 * typedef struct _braid_App_struct
646 * {
647 * HeatEquation<2> eq;
648 * int final_step;
649 * } my_App;
650 *
651 *
652 * /**
653 * * @brief my_Step - Takes a step in time, advancing the u vector
654 * *
655 * * @param app - The braid app struct
656 * * @param ustop - The solution data at the end of this time step
657 * * @param fstop - RHS data (such as forcing function?)
658 * * @param u - The solution data at the beginning of this time step
659 * * @param status - Status structure that contains various info of this time
660 * *
661 * * @return Success (0) or failure (1)
662 * **/
663 * int my_Step(braid_App app,
664 * braid_Vector ustop,
665 * braid_Vector fstop,
666 * braid_Vector u,
667 * braid_StepStatus status);
668 *
669 *
670 * /**
671 * * @brief my_Init - Initializes a solution data at the given time
672 * * For now, initializes the solution to zero no matter what time we are at
673 * *
674 * * @param app - The braid app struct containing user data
675 * * @param t - Time at which the solution is initialized
676 * * @param u_ptr - The solution data that needs to be filled
677 * *
678 * * @return Success (0) or failure (1)
679 * **/
680 * int
681 * my_Init(braid_App app,
682 * double t,
683 * braid_Vector *u_ptr);
684 *
685 *
686 * /**
687 * * @brief my_Clone - Clones a vector into a new vector
688 * *
689 * * @param app - The braid app struct containing user data
690 * * @param u - The existing vector containing data
691 * * @param v_ptr - The empty vector that needs to be filled
692 * *
693 * * @return Success (0) or failure (1)
694 * **/
695 * int
696 * my_Clone(braid_App app,
697 * braid_Vector u,
698 * braid_Vector *v_ptr);
699 *
700 *
701 * /**
702 * * @brief my_Free - Deletes a vector
703 * *
704 * * @param app - The braid app struct containing user data
705 * * @param u - The vector that needs to be deleted
706 * *
707 * * @return Success (0) or failure (1)
708 * **/
709 * int
710 * my_Free(braid_App app,
711 * braid_Vector u);
712 *
713 *
714 * /**
715 * * @brief my_Sum - Sums two vectors in an AXPY operation
716 * * The operation is y = alpha*x + beta*y
717 * *
718 * * @param app - The braid app struct containing user data
719 * * @param alpha - The coefficient in front of x
720 * * @param x - A vector that is multiplied by alpha then added to y
721 * * @param beta - The coefficient of y
722 * * @param y - A vector that is multiplied by beta then summed with x
723 * *
724 * * @return Success (0) or failure (1)
725 * **/
726 * int
727 * my_Sum(braid_App app,
728 * double alpha,
729 * braid_Vector x,
730 * double beta,
731 * braid_Vector y);
732 *
733 * /**
734 * * \brief Returns the spatial norm of the provided vector
735 * *
736 * * Calculates and returns the spatial norm of the provided vector.
737 * * Interestingly enough, X-Braid does not specify a particular norm.
738 * * to keep things simple, we implement the Euclidean norm.
739 * *
740 * * \param app - The braid app struct containing user data
741 * * \param u - The vector we need to take the norm of
742 * * \param norm_ptr - Pointer to the norm that was calculated, need to modify this
743 * * \return Success (0) or failure (1)
744 * */
745 * int
746 * my_SpatialNorm(braid_App app,
747 * braid_Vector u,
748 * double *norm_ptr);
749 *
750 * /**
751 * * \brief Allows the user to output details
752 * *
753 * * The Access function is called at various points to allow the user to output
754 * * information to the screen or to files.
755 * * The astatus parameter provides various information about the simulation,
756 * * see the XBraid documentation for details on what information you can get.
757 * * Example information is what the current timestep number and current time is.
758 * * If the access level (in parallel_in_time.cc) is set to 0, this function is
759 * * never called.
760 * * If the access level is set to 1, the function is called after the last
761 * * XBraid cycle.
762 * * If the access level is set to 2, it is called every XBraid cycle.
763 * *
764 * * \param app - The braid app struct containing user data
765 * * \param u - The vector containing the data at the status provided
766 * * \param astatus - The Braid status structure
767 * * \return Success (0) or failure (1)
768 * */
769 * int
770 * my_Access(braid_App app,
771 * braid_Vector u,
772 * braid_AccessStatus astatus);
773 *
774 * /**
775 * * \brief Calculates the size of a buffer for MPI data transfer
776 * *
777 * * Calculates the size of the buffer that is needed to transfer
778 * * a solution vector to another processor.
779 * * The bstatus parameter provides various information on the
780 * * simulation, see the XBraid documentation for all possible
781 * * fields.
782 * *
783 * * \param app - The braid app struct containing user data
784 * * \param size_ptr A pointer to the calculated size
785 * * \param bstatus The XBraid status structure
786 * * \return Success (0) or failure (1)
787 * */
788 * int
789 * my_BufSize(braid_App app,
790 * int *size_ptr,
791 * braid_BufferStatus bstatus);
792 *
793 * /**
794 * * \brief Linearizes a vector to be sent to another processor
795 * *
796 * * Linearizes (packs) a data buffer with the contents of
797 * * some solution state u.
798 * *
799 * * \param app - The braid app struct containing user data
800 * * \param u The vector that must be packed into buffer
801 * * \param buffer The buffer that must be filled with u
802 * * \param bstatus The XBraid status structure
803 * * \return Success (0) or failure (1)
804 * */
805 * int
806 * my_BufPack(braid_App app,
807 * braid_Vector u,
808 * void *buffer,
809 * braid_BufferStatus bstatus);
810 *
811 * /**
812 * * \brief Unpacks a vector that was sent from another processor
813 * *
814 * * Unpacks a linear data buffer into the vector pointed to by
815 * * u_ptr.
816 * *
817 * * \param app - The braid app struct containing user data
818 * * \param buffer The buffer that must be unpacked
819 * * \param u_ptr The pointer to the vector that is filled
820 * * \param bstatus The XBraid status structure
821 * * \return Success (0) or failure (1)
822 * */
823 * int
824 * my_BufUnpack(braid_App app,
825 * void *buffer,
826 * braid_Vector *u_ptr,
827 * braid_BufferStatus bstatus);
828 *
829 * #endif // _BRAIDFUNCS_H_
830 * @endcode
831
832
833<a name="ann-src/HeatEquation.hh"></a>
834<h1>Annotated version of src/HeatEquation.hh</h1>
835 *
836 *
837 *
838 *
839 * @code
840 * #ifndef _HEATEQUATION_H_
841 * #define _HEATEQUATION_H_
842 *
843 * #include <deal.II/base/utilities.h>
844 * #include <deal.II/base/quadrature_lib.h>
845 * #include <deal.II/base/function.h>
846 * #include <deal.II/base/logstream.h>
847 * #include <deal.II/lac/vector.h>
848 * #include <deal.II/lac/full_matrix.h>
849 * #include <deal.II/lac/dynamic_sparsity_pattern.h>
850 * #include <deal.II/lac/sparse_matrix.h>
851 * #include <deal.II/lac/solver_cg.h>
852 * #include <deal.II/lac/precondition.h>
853 * #include <deal.II/lac/affine_constraints.h>
854 * #include <deal.II/grid/tria.h>
855 * #include <deal.II/grid/grid_generator.h>
856 * #include <deal.II/grid/grid_refinement.h>
857 * #include <deal.II/grid/grid_out.h>
858 * #include <deal.II/grid/tria_accessor.h>
859 * #include <deal.II/grid/tria_iterator.h>
860 * #include <deal.II/dofs/dof_handler.h>
861 * #include <deal.II/dofs/dof_accessor.h>
862 * #include <deal.II/dofs/dof_tools.h>
863 * #include <deal.II/fe/fe_q.h>
864 * #include <deal.II/fe/fe_values.h>
865 * #include <deal.II/numerics/data_out.h>
866 * #include <deal.II/numerics/vector_tools.h>
867 * #include <deal.II/numerics/error_estimator.h>
868 * #include <deal.II/numerics/solution_transfer.h>
869 * #include <deal.II/numerics/matrix_tools.h>
870 * #include <deal.II/base/convergence_table.h>
871 *
872 * #include <fstream>
873 *
874 * using namespace dealii;
875 *
876 * @endcode
877 *
878 * The HeatEquation class is describes the finite element
879 * solver for the heat equation. It contains all the functions
880 * needed to define the problem domain and advance the solution
881 * in time.
882 *
883 * @code
884 * template <int dim>
885 * class HeatEquation
886 * {
887 * public:
888 * HeatEquation();
889 * void define();
890 * void step(Vector<double>& braid_data,
891 * double deltaT,
892 * double a_time,
893 * int a_time_idx);
894 *
895 * int size() const; /// Returns the size of the solution vector
896 *
897 * void output_results(int a_time_idx,
898 * double a_time,
899 * Vector<double>& a_solution) const;
900 *
901 * void initialize(double a_time,
902 * Vector<double>& a_vector) const;
903 *
904 * void process_solution(double a_time,
905 * int a_index,
906 * const Vector<double>& a_vector);
907 *
908 * private:
909 * void setup_system();
910 * void solve_time_step(Vector<double>& a_solution);
911 *
912 * Triangulation<dim> triangulation;
913 * FE_Q<dim> fe;
914 * DoFHandler<dim> dof_handler;
915 *
916 * AffineConstraints<double> constraints;
917 *
918 * SparsityPattern sparsity_pattern;
919 * SparseMatrix<double> mass_matrix;
920 * SparseMatrix<double> laplace_matrix;
921 * SparseMatrix<double> system_matrix;
922 *
923 * Vector<double> system_rhs;
924 *
925 * std::ofstream myfile;
926 *
927 * const double theta;
928 *
929 * @endcode
930 *
931 * These were originally in the run() function but because
932 * I am splitting the run() function up into define and step
933 * they need to become member data
934 *
935 * @code
936 * Vector<double> tmp;
937 * Vector<double> forcing_terms;
938 *
939 * ConvergenceTable convergence_table;
940 * };
941 *
942 * @endcode
943 *
944 * The RightHandSide class describes the RHS of the governing
945 * equations. In this case, it is the forcing function.
946 *
947 * @code
948 * template <int dim>
949 * class RightHandSide : public Function<dim>
950 * {
951 * public:
952 * RightHandSide ()
953 * :
954 * Function<dim>(),
955 * period (0.2)
956 * {}
957 *
958 * virtual double value (const Point<dim> &p,
959 * const unsigned int component = 0) const;
960 *
961 * private:
962 * const double period;
963 * };
964 *
965 * @endcode
966 *
967 * The BoundaryValues class describes the boundary conditions
968 * of the governing equations.
969 *
970 * @code
971 * template <int dim>
972 * class BoundaryValues : public Function<dim>
973 * {
974 * public:
975 * virtual double value (const Point<dim> &p,
976 * const unsigned int component = 0) const;
977 * };
978 *
979 * @endcode
980 *
981 * The RightHandSideMFG class describes the right hand side
982 * function when doing the method of manufactured solutions.
983 *
984 * @code
985 * template <int dim>
986 * class RightHandSideMFG : public Function<dim>
987 * {
988 * public:
989 * virtual double value (const Point<dim> &p,
990 * const unsigned int component = 0) const;
991 * };
992 *
993 * @endcode
994 *
995 * The InitialValuesMFG class describes the initial values
996 * when doing the method of manufactured solutions.
997 *
998 * @code
999 * template <int dim>
1000 * class InitialValuesMFG : public Function<dim>
1001 * {
1002 * public:
1003 * virtual double value (const Point<dim> &p,
1004 * const unsigned int component = 0) const;
1005 * };
1006 *
1007 * @endcode
1008 *
1009 * Provides the exact value for the manufactured solution. This
1010 * is used for the boundary conditions as well.
1011 *
1012 * @code
1013 * template <int dim>
1014 * class ExactValuesMFG : public Function<dim>
1015 * {
1016 * public:
1017 * /**
1018 * * \brief Computes the value at the given point and member data time
1019 * *
1020 * * Computes the exact value of the manufactured solution at point p and
1021 * * the member data time. See the class documentation and the design doc
1022 * * for details on what the exact solution is.
1023 * *
1024 * * \param p The point that the exact solution is computed at
1025 * * \param component The component of the exact solution (always 0 for now)
1026 * * \return double The exact value that was computed
1027 * */
1028 * virtual double value (const Point<dim> &p,
1029 * const unsigned int component = 0) const;
1030 *
1031 * /**
1032 * * \brief Computes the gradient of the exact solution at the given point
1033 * *
1034 * * Computes the gradient of the exact/manufactured solution value at
1035 * * point p and member data time. See the design doc for details on
1036 * * what the gradient of the exact solution is
1037 * *
1038 * * \param p The point that the gradient is calculated at
1039 * * \param component The component of the system of equations this gradient is for
1040 * * \return Tensor<1,dim> A rank 1 tensor that contains the gradient
1041 * * in each spatial dimension
1042 * */
1043 * virtual Tensor<1,dim> gradient (const Point<dim> &p,
1044 * const unsigned int component = 0) const;
1045 * };
1046 *
1047 *
1048 * #include "HeatEquationImplem.hh"
1049 *
1050 * #endif // _HEATEQUATION_H_
1051 * @endcode
1052
1053
1054<a name="ann-src/HeatEquationImplem.hh"></a>
1055<h1>Annotated version of src/HeatEquationImplem.hh</h1>
1056 *
1057 *
1058 *
1059 *
1060 * @code
1061 * #include "Utilities.hh"
1062 *
1063 * #include <iomanip>
1064 * #include <math.h>
1065 *
1066 * @endcode
1067 *
1068 * Calculates the forcing function for the RightHandSide. See the
1069 * documentation for the math.
1070 *
1071 * @code
1072 * template <int dim>
1073 * double RightHandSide<dim>::value (const Point<dim> &p,
1074 * const unsigned int component) const
1075 * {
1076 * (void) component;
1077 * Assert (component == 0, ExcIndexRange(component, 0, 1));
1078 * Assert (dim == 2, ExcNotImplemented());
1079 *
1080 * double time = this->get_time();
1081 *
1082 * if ((p[0] > 0.5) && (p[1] > -0.5))
1083 * {
1084 * return std::exp(-0.5*(time-0.125)*(time-0.125)/(0.005));
1085 * }
1086 * else if ((p[0] > -0.5) && (p[1] > 0.5))
1087 * {
1088 * return std::exp(-0.5*(time-0.375)*(time-0.375)/(0.005));
1089 * }
1090 * else
1091 * {
1092 * return 0;
1093 * }
1094 *
1095 * return 0; // No forcing function
1096 * }
1097 *
1098 * @endcode
1099 *
1100 * Calculates the forcing function for the method of manufactured
1101 * solutions. See the documentation for the math.
1102 *
1103 * @code
1104 * template <int dim>
1105 * double RightHandSideMFG<dim>::value (const Point<dim> &p,
1106 * const unsigned int component) const
1107 * {
1108 * (void) component;
1109 * Assert (component == 0, ExcIndexRange(component, 0, 1));
1110 * Assert (dim == 2, ExcNotImplemented());
1111 *
1112 * double time = this->get_time();
1113 *
1114 * double pi = numbers::PI;
1115 * return 4*pi*pi*std::exp(-4*pi*pi*time)*std::cos(2*pi*p[0])*std::cos(2*pi*p[1]);
1116 * }
1117 *
1118 * @endcode
1119 *
1120 * Calculates the boundary conditions, essentially zero everywhere.
1121 *
1122 * @code
1123 * template <int dim>
1124 * double BoundaryValues<dim>::value (const Point<dim> &p,
1125 * const unsigned int component) const
1126 * {
1127 * UNUSED(p);
1128 * (void) component;
1129 * Assert (component == 0, ExcIndexRange(component, 0, 1));
1130 * return 0;
1131 * }
1132 *
1133 * @endcode
1134 *
1135 * Calculates the exact solution (and thus also boundary conditions)
1136 * for the method of manufactured solutions.
1137 *
1138 * @code
1139 * template <int dim>
1140 * double ExactValuesMFG<dim>::value (const Point<dim> &p,
1141 * const unsigned int component) const
1142 * {
1143 * (void) component;
1144 * Assert (component == 0, ExcIndexRange(component, 0, 1));
1145 *
1146 * double time = this->get_time();
1147 * const double pi = numbers::PI;
1148 *
1149 * return std::exp(-4*pi*pi*time)*std::cos(2*pi*p[0])*std::cos(2*pi*p[1]);
1150 * }
1151 *
1152 * @endcode
1153 *
1154 * Calculates the gradient of the exact solution for the method of manufactured
1155 * solutions. See the documentation for the math.
1156 *
1157 * @code
1158 * template <int dim>
1159 * Tensor<1,dim> ExactValuesMFG<dim>::gradient (const Point<dim> &p,
1160 * const unsigned int) const
1161 * {
1162 * Assert (dim == 2, ExcNotImplemented());
1163 *
1164 * Tensor<1,dim> return_value;
1165 * const double pi = numbers::PI;
1166 * double time = this->get_time();
1167 * return_value[0] = -2*pi*std::exp(-4*pi*pi*time)*std::cos(2*pi*p[1])*std::sin(2*pi*p[0]);
1168 * return_value[1] = -2*pi*std::exp(-4*pi*pi*time)*std::cos(2*pi*p[0])*std::sin(2*pi*p[1]);
1169 * return return_value;
1170 * }
1171 *
1172 * @endcode
1173 *
1174 * Calculates the initial values for the method of manufactured solutions.
1175 * See the documentation for the math.
1176 *
1177 * @code
1178 * template <int dim>
1179 * double InitialValuesMFG<dim>::value (const Point<dim> &p,
1180 * const unsigned int component) const
1181 * {
1182 * (void) component;
1183 * Assert (component == 0, ExcIndexRange(component, 0, 1));
1184 * const double pi = numbers::PI;
1185 *
1186 * return std::cos(2*pi*p[0])*std::cos(2*pi*p[1]);
1187 * }
1188 *
1189 * template <int dim>
1190 * HeatEquation<dim>::HeatEquation ()
1191 * :
1192 * fe(1),
1193 * dof_handler(triangulation),
1194 * theta(0.5)
1195 * {
1196 * }
1197 *
1198 * template <int dim>
1199 * void HeatEquation<dim>::initialize(double a_time,
1200 * Vector<double>& a_vector) const
1201 * {
1202 * #if DO_MFG
1203 * @endcode
1204 *
1205 * We only initialize values in the manufactured solution case
1206 *
1207 * @code
1208 * InitialValuesMFG<dim> iv_function;
1209 * iv_function.set_time(a_time);
1210 * VectorTools::project (dof_handler, constraints,
1211 * QGauss<dim>(fe.degree+1), iv_function,
1212 * a_vector);
1213 * #else
1214 * UNUSED(a_time);
1215 * UNUSED(a_vector);
1216 * #endif // DO_MFG
1217 * @endcode
1218 *
1219 * If not the MFG solution case, a_vector is already zero'd so do nothing
1220 *
1221 * @code
1222 * }
1223 *
1224 * template <int dim>
1225 * void HeatEquation<dim>::setup_system()
1226 * {
1227 * dof_handler.distribute_dofs(fe);
1228 *
1229 * constraints.clear ();
1231 * constraints);
1232 * constraints.close();
1233 *
1234 * DynamicSparsityPattern dsp(dof_handler.n_dofs());
1235 * DoFTools::make_sparsity_pattern(dof_handler,
1236 * dsp,
1237 * constraints,
1238 * /*keep_constrained_dofs = */ true);
1239 * sparsity_pattern.copy_from(dsp);
1240 *
1241 * mass_matrix.reinit(sparsity_pattern);
1242 * laplace_matrix.reinit(sparsity_pattern);
1243 * system_matrix.reinit(sparsity_pattern);
1244 *
1246 * QGauss<dim>(fe.degree+1),
1247 * mass_matrix);
1249 * QGauss<dim>(fe.degree+1),
1250 * laplace_matrix);
1251 *
1252 * system_rhs.reinit(dof_handler.n_dofs());
1253 * }
1254 *
1255 *
1256 * template <int dim>
1257 * void HeatEquation<dim>::solve_time_step(Vector<double>& a_solution)
1258 * {
1259 * SolverControl solver_control(1000, 1e-8 * system_rhs.l2_norm());
1260 * SolverCG<> cg(solver_control);
1261 *
1262 * PreconditionSSOR<> preconditioner;
1263 * preconditioner.initialize(system_matrix, 1.0);
1264 *
1265 * cg.solve(system_matrix, a_solution, system_rhs,
1266 * preconditioner);
1267 *
1268 * constraints.distribute(a_solution);
1269 * }
1270 *
1271 *
1272 *
1273 * template <int dim>
1274 * void HeatEquation<dim>::output_results(int a_time_idx,
1275 * double a_time,
1276 * Vector<double>& a_solution) const
1277 * {
1278 *
1279 * DataOutBase::VtkFlags vtk_flags;
1280 * vtk_flags.time = a_time;
1281 * vtk_flags.cycle = a_time_idx;
1282 *
1283 * DataOut<dim> data_out;
1284 * data_out.set_flags(vtk_flags);
1285 *
1286 * data_out.attach_dof_handler(dof_handler);
1287 * data_out.add_data_vector(a_solution, "U");
1288 *
1289 * data_out.build_patches();
1290 *
1291 * const std::string filename = "solution-"
1292 * + Utilities::int_to_string(a_time_idx, 3) +
1293 * ".vtk";
1294 * std::ofstream output(filename.c_str());
1295 * data_out.write_vtk(output);
1296 * }
1297 *
1298 * @endcode
1299 *
1300 * We define the geometry here, this is called on each processor
1301 * and doesn't change in time. Once doing AMR, this won't need
1302 * to exist anymore.
1303 *
1304 * @code
1305 * template <int dim>
1306 * void HeatEquation<dim>::define()
1307 * {
1308 * const unsigned int initial_global_refinement = 6;
1309 *
1311 * triangulation.refine_global (initial_global_refinement);
1312 *
1313 * setup_system();
1314 *
1315 * tmp.reinit (dof_handler.n_dofs());
1316 * forcing_terms.reinit (dof_handler.n_dofs());
1317 * }
1318 *
1319 * @endcode
1320 *
1321 * Here we advance the solution forward in time. This is done
1322 * the same way as in the loop in @ref step_26 "step-26"'s run function.
1323 *
1324 * @code
1325 * template<int dim>
1326 * void HeatEquation<dim>::step(Vector<double>& braid_data,
1327 * double deltaT,
1328 * double a_time,
1329 * int a_time_idx)
1330 * {
1331 * a_time += deltaT;
1332 * ++a_time_idx;
1333 *
1334 * mass_matrix.vmult(system_rhs, braid_data);
1335 *
1336 * laplace_matrix.vmult(tmp, braid_data);
1337 *
1338 * system_rhs.add(-(1 - theta) * deltaT, tmp);
1339 *
1340 * #if DO_MFG
1341 * RightHandSideMFG<dim> rhs_function;
1342 * #else
1343 * RightHandSide<dim> rhs_function;
1344 * #endif
1345 * rhs_function.set_time(a_time);
1346 * VectorTools::create_right_hand_side(dof_handler,
1347 * QGauss<dim>(fe.degree+1),
1348 * rhs_function,
1349 * tmp);
1350 *
1351 * forcing_terms = tmp;
1352 * forcing_terms *= deltaT * theta;
1353 *
1354 * rhs_function.set_time(a_time - deltaT);
1355 * VectorTools::create_right_hand_side(dof_handler,
1356 * QGauss<dim>(fe.degree+1),
1357 * rhs_function,
1358 * tmp);
1359 *
1360 * forcing_terms.add(deltaT * (1 - theta), tmp);
1361 * system_rhs += forcing_terms;
1362 *
1363 * system_matrix.copy_from(mass_matrix);
1364 * system_matrix.add(theta * deltaT, laplace_matrix);
1365 *
1366 * constraints.condense (system_matrix, system_rhs);
1367 *
1368 * {
1369 * #if DO_MFG
1370 * @endcode
1371 *
1372 * If we are doing the method of manufactured solutions
1373 * then we set the boundary conditions to the exact solution.
1374 * Otherwise the boundary conditions are zero.
1375 *
1376 * @code
1377 * ExactValuesMFG<dim> boundary_values_function;
1378 * #else
1379 * BoundaryValues<dim> boundary_values_function;
1380 * #endif
1381 * boundary_values_function.set_time(a_time);
1382 *
1383 * std::map<types::global_dof_index, double> boundary_values;
1384 * VectorTools::interpolate_boundary_values(dof_handler,
1385 * 0,
1386 * boundary_values_function,
1387 * boundary_values);
1388 *
1389 * MatrixTools::apply_boundary_values(boundary_values,
1390 * system_matrix,
1391 * braid_data,
1392 * system_rhs);
1393 * }
1394 *
1395 * solve_time_step(braid_data);
1396 * }
1397 *
1398 * template<int dim>
1399 * int HeatEquation<dim>::size() const
1400 * {
1401 * return dof_handler.n_dofs();
1402 * }
1403 *
1404 * @endcode
1405 *
1406 * This function computes the error for the time step when doing
1407 * the method of manufactured solutions. First the exact values
1408 * is calculated, then the difference per cell is computed for
1409 * the various norms, and the error is computed. This is written
1410 * out to a pretty table.
1411 *
1412 * @code
1413 * template<int dim> void
1414 * HeatEquation<dim>::process_solution(double a_time,
1415 * int a_index,
1416 * const Vector<double>& a_vector)
1417 * {
1418 * @endcode
1419 *
1420 * Compute the exact value for the manufactured solution case
1421 *
1422 * @code
1423 * ExactValuesMFG<dim> exact_function;
1424 * exact_function.set_time(a_time);
1425 *
1426 * Vector<double> difference_per_cell (triangulation.n_active_cells());
1427 * VectorTools::integrate_difference(dof_handler,
1428 * a_vector,
1429 * exact_function,
1430 * difference_per_cell,
1431 * QGauss<dim>(fe.degree+1),
1432 * VectorTools::L2_norm);
1433 *
1434 * const double L2_error = VectorTools::compute_global_error(triangulation,
1435 * difference_per_cell,
1436 * VectorTools::L2_norm);
1437 *
1438 * VectorTools::integrate_difference(dof_handler,
1439 * a_vector,
1440 * exact_function,
1441 * difference_per_cell,
1442 * QGauss<dim>(fe.degree+1),
1443 * VectorTools::H1_seminorm);
1444 *
1445 * const double H1_error = VectorTools::compute_global_error(triangulation,
1446 * difference_per_cell,
1447 * VectorTools::H1_seminorm);
1448 *
1449 * const QTrapez<1> q_trapez;
1450 * const QIterated<dim> q_iterated (q_trapez, 5);
1451 * VectorTools::integrate_difference (dof_handler,
1452 * a_vector,
1453 * exact_function,
1454 * difference_per_cell,
1455 * q_iterated,
1456 * VectorTools::Linfty_norm);
1457 * const double Linfty_error = VectorTools::compute_global_error(triangulation,
1458 * difference_per_cell,
1459 * VectorTools::Linfty_norm);
1460 *
1461 * const unsigned int n_active_cells = triangulation.n_active_cells();
1462 * const unsigned int n_dofs = dof_handler.n_dofs();
1463 *
1464 * pout() << "Cycle " << a_index << ':'
1465 * << std::endl
1466 * << " Number of active cells: "
1467 * << n_active_cells
1468 * << std::endl
1469 * << " Number of degrees of freedom: "
1470 * << n_dofs
1471 * << std::endl;
1472 *
1473 * convergence_table.add_value("cycle", a_index);
1474 * convergence_table.add_value("cells", n_active_cells);
1475 * convergence_table.add_value("dofs", n_dofs);
1476 * convergence_table.add_value("L2", L2_error);
1477 * convergence_table.add_value("H1", H1_error);
1478 * convergence_table.add_value("Linfty", Linfty_error);
1479 *
1480 * convergence_table.set_precision("L2", 3);
1481 * convergence_table.set_precision("H1", 3);
1482 * convergence_table.set_precision("Linfty", 3);
1483 *
1484 * convergence_table.set_scientific("L2", true);
1485 * convergence_table.set_scientific("H1", true);
1486 * convergence_table.set_scientific("Linfty", true);
1487 *
1488 * convergence_table.set_tex_caption("cells", "\\# cells");
1489 * convergence_table.set_tex_caption("dofs", "\\# dofs");
1490 * convergence_table.set_tex_caption("L2", "@fL^2@f-error");
1491 * convergence_table.set_tex_caption("H1", "@fH^1@f-error");
1492 * convergence_table.set_tex_caption("Linfty", "@fL^\\infty@f-error");
1493 *
1494 * convergence_table.set_tex_format("cells", "r");
1495 * convergence_table.set_tex_format("dofs", "r");
1496 *
1497 * std::cout << std::endl;
1498 * convergence_table.write_text(std::cout);
1499 *
1500 * std::ofstream error_table_file("tex-conv-table.tex");
1501 * convergence_table.write_tex(error_table_file);
1502 * }
1503 * @endcode
1504
1505
1506<a name="ann-src/Utilities.cc"></a>
1507<h1>Annotated version of src/Utilities.cc</h1>
1508 *
1509 *
1510 *
1511 *
1512 * @code
1513 * #include "Utilities.hh"
1514 *
1515 * #include <string>
1516 * #include <fstream>
1517 *
1518 * #include <mpi.h>
1519 *
1520 * int procID = 0;
1521 *
1522 * @endcode
1523 *
1524 * The shared variables
1525 *
1526
1527 *
1528 *
1529 * @code
1530 * static std::string s_pout_filename ;
1531 * static std::string s_pout_basename ;
1532 * static std::ofstream s_pout ;
1533 *
1534 * static bool s_pout_init = false ;
1535 * static bool s_pout_open = false ;
1536 *
1537 * #ifdef USE_MPI
1538 * @endcode
1539 *
1540 * in parallel, compute the filename give the basename
1541 * [NOTE: dont call this before MPI is initialized.]
1542 *
1543 * @code
1544 * static void setFileName()
1545 * {
1546 * static const size_t ProcnumSize = 1 + 10 + 1 ; //'.' + 10digits + '\0'
1547 * char procnum[ProcnumSize] ;
1548 * snprintf( procnum ,ProcnumSize ,".%d" ,procID);
1549 * s_pout_filename = s_pout_basename + procnum ;
1550 * }
1551 *
1552 * @endcode
1553 *
1554 * in parallel, close the file if nec., open it and check for success
1555 *
1556 * @code
1557 * static void openFile()
1558 * {
1559 * if ( s_pout_open )
1560 * {
1561 * s_pout.close();
1562 * }
1563 * s_pout.open( s_pout_filename.c_str() );
1564 * @endcode
1565 *
1566 * if open() fails, we have problems, but it's better
1567 * to try again later than to make believe it succeeded
1568 *
1569 * @code
1570 * s_pout_open = (bool)s_pout ;
1571 * }
1572 *
1573 * #else
1574 * @endcode
1575 *
1576 * in serial, filename is always cout
1577 *
1578 * @code
1579 * static void setFileName()
1580 * {
1581 * s_pout_filename = "cout" ;
1582 * }
1583 *
1584 * @endcode
1585 *
1586 * in serial, this does absolutely nothing
1587 *
1588 * @code
1589 * static void openFile()
1590 * {
1591 * }
1592 * #endif
1593 *
1594 * std::ostream& pout()
1595 * {
1596 * #ifdef USE_MPI
1597 * @endcode
1598 *
1599 * the common case is _open == true, which just returns s_pout
1600 *
1601 * @code
1602 * if ( ! s_pout_open )
1603 * {
1604 * @endcode
1605 *
1606 * the uncommon cae: the file isn't opened, MPI may not be
1607 * initialized, and the basename may not have been set
1608 *
1609 * @code
1610 * int flag_i, flag_f;
1611 * MPI_Initialized(&flag_i);
1612 * MPI_Finalized(&flag_f);
1613 * @endcode
1614 *
1615 * app hasn't set a basename yet, so set the default
1616 *
1617 * @code
1618 * if ( ! s_pout_init )
1619 * {
1620 * s_pout_basename = "pout" ;
1621 * s_pout_init = true ;
1622 * }
1623 * @endcode
1624 *
1625 * if MPI not initialized, we cant open the file so return cout
1626 *
1627 * @code
1628 * if ( ! flag_i || flag_f)
1629 * {
1630 * return std::cout; // MPI hasn't been started yet, or has ended....
1631 * }
1632 * @endcode
1633 *
1634 * MPI is initialized, so file must not be, so open it
1635 *
1636 * @code
1637 * setFileName() ;
1638 * openFile() ;
1639 * @endcode
1640 *
1641 * finally, in case the open failed, return cout
1642 *
1643 * @code
1644 * if ( ! s_pout_open )
1645 * {
1646 * return std::cout ;
1647 * }
1648 * }
1649 * return s_pout ;
1650 * #else
1651 * return std::cout;
1652 * #endif
1653 * }
1654 * @endcode
1655
1656
1657<a name="ann-src/Utilities.hh"></a>
1658<h1>Annotated version of src/Utilities.hh</h1>
1659 *
1660 *
1661 *
1662 *
1663 * @code
1664 * #ifndef _UTILITIES_H_
1665 * #define _UTILITIES_H_
1666 *
1667 * #include <iostream>
1668 *
1669 * @endcode
1670 *
1671 * This preprocessor macro is used on function arguments
1672 * that are not used in the function. It is used to
1673 * suppress compiler warnings.
1674 *
1675 * @code
1676 * #define UNUSED(x) (void)(x)
1677 *
1678 * @endcode
1679 *
1680 * Contains the current MPI processor ID.
1681 *
1682 * @code
1683 * extern int procID;
1684 *
1685 * @endcode
1686 *
1687 * Function to return the ostream to write out to. In MPI
1688 * mode it returns a stream to a file named pout.<#> where
1689 * <#> is the procID. This allows the user to write output
1690 * from each processor to a separate file. In serial mode
1691 * (no MPI), it returns the standard output.
1692 *
1693 * @code
1694 * std::ostream& pout();
1695 * #endif
1696 * @endcode
1697
1698
1699<a name="ann-src/parallel_in_time.cc"></a>
1700<h1>Annotated version of src/parallel_in_time.cc</h1>
1701 *
1702 *
1703 *
1704 *
1705 * @code
1706 * /* ---------------------------------------------------------------------
1707 * *
1708 * * Copyright (C) 2013 - 2018 by the deal.II authors
1709 * *
1710 * * This file is part of the deal.II library.
1711 * *
1712 * * The deal.II library is free software; you can use it, redistribute
1713 * * it, and/or modify it under the terms of the GNU Lesser General
1714 * * Public License as published by the Free Software Foundation; either
1715 * * version 2.1 of the License, or (at your option) any later version.
1716 * * The full text of the license can be found in the file LICENSE at
1717 * * the top level of the deal.II distribution.
1718 * *
1719 * * ---------------------------------------------------------------------
1720 *
1721 * *
1722 * * Author: Joshua Christopher, Colorado State University, 2018
1723 * */
1724 *
1725 * #include "BraidFuncs.hh"
1726 * #include "HeatEquation.hh"
1727 * #include "Utilities.hh"
1728 *
1729 * #include <fstream>
1730 * #include <iostream>
1731 *
1732 * int main(int argc, char *argv[])
1733 * {
1734 * try
1735 * {
1736 * using namespace dealii;
1737 *
1738 * /* Initialize MPI */
1739 * MPI_Comm comm; //, comm_x, comm_t;
1740 * int rank;
1741 * MPI_Init(&argc, &argv);
1742 * comm = MPI_COMM_WORLD;
1743 * MPI_Comm_rank(comm, &rank);
1744 * procID = rank;
1745 *
1746 * @endcode
1747 *
1748 * Set up X-Braid
1749 *
1750 * @code
1751 * /* Initialize Braid */
1752 * braid_Core core;
1753 * double tstart = 0.0;
1754 * double tstop = 0.002;
1755 * int ntime = 10;
1756 * my_App *app = new(my_App);
1757 *
1758 * braid_Init(MPI_COMM_WORLD, comm, tstart, tstop, ntime, app,
1759 * my_Step, my_Init, my_Clone, my_Free, my_Sum, my_SpatialNorm,
1760 * my_Access, my_BufSize, my_BufPack, my_BufUnpack, &core);
1761 *
1762 * /* Define XBraid parameters
1763 * * See -help message forf descriptions */
1764 * int max_levels = 3;
1765 * @endcode
1766 *
1767 * int nrelax = 1;
1768 * int skip = 0;
1769 *
1770 * @code
1771 * double tol = 1.e-7;
1772 * @endcode
1773 *
1774 * int cfactor = 2;
1775 *
1776 * @code
1777 * int max_iter = 5;
1778 * @endcode
1779 *
1780 * int min_coarse = 10;
1781 * int fmg = 0;
1782 * int scoarsen = 0;
1783 * int res = 0;
1784 * int wrapper_tests = 0;
1785 *
1786 * @code
1787 * int print_level = 1;
1788 * int access_level = 1;
1789 * int use_sequential= 0;
1790 *
1791 * braid_SetPrintLevel( core, print_level);
1792 * braid_SetAccessLevel( core, access_level);
1793 * braid_SetMaxLevels(core, max_levels);
1794 * @endcode
1795 *
1796 * braid_SetMinCoarse( core, min_coarse );
1797 * braid_SetSkip(core, skip);
1798 * braid_SetNRelax(core, -1, nrelax);
1799 *
1800 * @code
1801 * braid_SetAbsTol(core, tol);
1802 * @endcode
1803 *
1804 * braid_SetCFactor(core, -1, cfactor);
1805 *
1806 * @code
1807 * braid_SetMaxIter(core, max_iter);
1808 * braid_SetSeqSoln(core, use_sequential);
1809 *
1810 * app->eq.define();
1811 * app->final_step = ntime;
1812 *
1813 * braid_Drive(core);
1814 *
1815 * @endcode
1816 *
1817 * Free the memory now that we are done
1818 *
1819 * @code
1820 * braid_Destroy(core);
1821 *
1822 * delete app;
1823 *
1824 * @endcode
1825 *
1826 * Clean up MPI
1827 * MPI_Comm_free(&comm);
1828 *
1829 * @code
1830 * MPI_Finalize();
1831 * }
1832 * catch (std::exception &exc)
1833 * {
1834 * std::cerr << std::endl << std::endl
1835 * << "----------------------------------------------------"
1836 * << std::endl;
1837 * std::cerr << "Exception on processing: " << std::endl << exc.what()
1838 * << std::endl << "Aborting!" << std::endl
1839 * << "----------------------------------------------------"
1840 * << std::endl;
1841 *
1842 * return 1;
1843 * }
1844 * catch (...)
1845 * {
1846 * std::cerr << std::endl << std::endl
1847 * << "----------------------------------------------------"
1848 * << std::endl;
1849 * std::cerr << "Unknown exception!" << std::endl << "Aborting!"
1850 * << std::endl
1851 * << "----------------------------------------------------"
1852 * << std::endl;
1853 * return 1;
1854 * }
1855 *
1856 * return 0;
1857 * }
1858 *
1859 * @endcode
1860
1861
1862<a name="ann-test/test_braid.cc"></a>
1863<h1>Annotated version of test/test_braid.cc</h1>
1864 *
1865 *
1866 *
1867 *
1868 * @code
1869 * #include "BraidFuncs.hh"
1870 *
1871 * #include <braid.h>
1872 * #include <braid_test.h>
1873 *
1874 * #include <iostream>
1875 *
1876 * int main(int argc, char** argv)
1877 * {
1878 * MPI_Comm comm;
1879 * int rank;
1880 * MPI_Init(&argc, &argv);
1881 * comm = MPI_COMM_WORLD;
1882 * MPI_Comm_rank(comm, &rank);
1883 *
1884 * my_App *app = new(my_App);
1885 * app->eq.define();
1886 *
1887 * double time = 0.2;
1888 *
1889 * braid_Int init_access_result = braid_TestInitAccess(app,
1890 * comm,
1891 * stdout,
1892 * time,
1893 * my_Init,
1894 * my_Access,
1895 * my_Free);
1896 * (void)init_access_result;
1897 *
1898 * braid_Int clone_result = braid_TestClone(app,
1899 * comm,
1900 * stdout,
1901 * time,
1902 * my_Init,
1903 * my_Access,
1904 * my_Free,
1905 * my_Clone);
1906 * (void)clone_result;
1907 *
1908 * braid_Int sum_result = braid_TestSum(app,
1909 * comm,
1910 * stdout,
1911 * time,
1912 * my_Init,
1913 * my_Access,
1914 * my_Free,
1915 * my_Clone,
1916 * my_Sum);
1917 * (void)sum_result;
1918 *
1919 * braid_Int norm_result = braid_TestSpatialNorm(app,
1920 * comm,
1921 * stdout,
1922 * time,
1923 * my_Init,
1924 * my_Free,
1925 * my_Clone,
1926 * my_Sum,
1927 * my_SpatialNorm);
1928 * (void)norm_result;
1929 *
1930 * braid_Int buf_result = braid_TestBuf(app,
1931 * comm,
1932 * stdout,
1933 * time,
1934 * my_Init,
1935 * my_Free,
1936 * my_Sum,
1937 * my_SpatialNorm,
1938 * my_BufSize,
1939 * my_BufPack,
1940 * my_BufUnpack);
1941 * (void)buf_result;
1942 *
1943 * @endcode
1944 *
1945 * /* Create spatial communicator for wrapper-tests */
1946 * braid_SplitCommworld(&comm, 1, &comm_x, &comm_t);
1947 *
1948 *
1949 * braid_TestAll(app, comm_x, stdout, 0.0, (tstop-tstart)/ntime,
1950 * 2*(tstop-tstart)/ntime, my_Init, my_Free, my_Clone,
1951 * my_Sum, my_SpatialNorm, my_BufSize, my_BufPack,
1952 * my_BufUnpack, my_Coarsen, my_Interp, my_Residual, my_Step);
1953 *
1954
1955 *
1956 *
1957 * @code
1958 * /* Finalize MPI */
1960 * }
1961 * @endcode
1962
1963
1964*/
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
Definition: vector.h:110
unsigned int level
Definition: grid_out.cc:4590
__global__ void set(Number *val, const Number s, const size_type N)
void set_flags(const FlagType &flags)
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())
pointer data()
void make_hanging_node_constraints(const DoFHandler< dim, spacedim > &dof_handler, AffineConstraints< number > &constraints)
void make_sparsity_pattern(const DoFHandler< dim, spacedim > &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)
void hyper_L(Triangulation< dim > &tria, const double left=-1., const double right=1., const bool colorize=false)
void mass_matrix(FullMatrix< double > &M, const FEValuesBase< dim > &fe, const double factor=1.)
Definition: l2.h:58
void create_laplace_matrix(const Mapping< dim, spacedim > &mapping, const DoFHandler< dim, spacedim > &dof, const Quadrature< dim > &q, SparseMatrix< double > &matrix, const Function< spacedim > *const a=nullptr, const AffineConstraints< double > &constraints=AffineConstraints< double >())
void create_mass_matrix(const Mapping< dim, spacedim > &mapping, const DoFHandler< dim, spacedim > &dof, const Quadrature< dim > &q, SparseMatrix< number > &matrix, const Function< spacedim, number > *const a=nullptr, const AffineConstraints< number > &constraints=AffineConstraints< number >())
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)
std::string int_to_string(const unsigned int value, const unsigned int digits=numbers::invalid_unsigned_int)
Definition: utilities.cc:473
const ::parallel::distributed::Triangulation< dim, spacedim > * triangulation
*** braid_TestAll(app, comm_x, stdout, 0.0,(tstop-tstart)/ntime, *2 *(tstop-tstart)/ntime, my_Init, my_Free, my_Clone, *my_Sum, my_SpatialNorm, my_BufSize, my_BufPack, *my_BufUnpack, my_Coarsen, my_Interp, my_Residual, my_Step)
****code ** MPI_Finalize()
*braid_SplitCommworld & comm
void advance(std::tuple< I1, I2 > &t, const unsigned int n)