Reference documentation for deal.II version Git e12e976e47 2020-11-24 09:52:11 +0100
\(\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\}}\)
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  *
1245  * MatrixCreator::create_mass_matrix(dof_handler,
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 */
1959  * MPI_Finalize();
1960  * }
1961  * @endcode
1962 
1963 
1964 */
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 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 initialize(const MatrixType &A, const typename BaseClass::AdditionalData &parameters=typename BaseClass::AdditionalData())
void set_flags(const FlagType &flags)
*** 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)
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
****code ** MPI_Finalize()
void attach_dof_handler(const DoFHandlerType &)
pointer data()
void write_vtk(std::ostream &out) const
void hyper_L(Triangulation< dim > &tria, const double left=-1., const double right=1., const bool colorize=false)
virtual void build_patches(const unsigned int n_subdivisions=0)
Definition: data_out.cc:1085
void mass_matrix(FullMatrix< double > &M, const FEValuesBase< dim > &fe, const double factor=1.)
Definition: l2.h:58
unsigned int level
Definition: grid_out.cc:4343
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 >())
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 >())
std::string int_to_string(const unsigned int value, const unsigned int digits=numbers::invalid_unsigned_int)
Definition: utilities.cc:474
Point< spacedim > point(const gp_Pnt &p, const double tolerance=1e-10)
Definition: utilities.cc:188
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
Definition: cuda.h:32
*braid_SplitCommworld & comm
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)
const ::parallel::distributed::Triangulation< dim, spacedim > * triangulation
void advance(std::tuple< I1, I2 > &t, const unsigned int n)
void make_hanging_node_constraints(const DoFHandler< dim, spacedim > &dof_handler, AffineConstraints< number > &constraints)