From: heister Date: Fri, 8 Feb 2013 23:21:39 +0000 (+0000) Subject: add test from Martin Kronbichler X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=ca54354a65aa72bb088a85c918182a51eb58b698;p=dealii-svn.git add test from Martin Kronbichler git-svn-id: https://svn.dealii.org/trunk@28287 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/tests/benchmarks/test_poisson/Makefile b/tests/benchmarks/test_poisson/Makefile new file mode 100644 index 0000000000..5958c739e0 --- /dev/null +++ b/tests/benchmarks/test_poisson/Makefile @@ -0,0 +1,144 @@ +# $Id: Makefile 25724 2012-07-24 23:35:36Z bangerth $ + + +# For the small projects Makefile, you basically need to fill in only +# four fields. +# +# The first is the name of the application. It is assumed that the +# application name is the same as the base file name of the single C++ +# file from which the application is generated. +target = poisson + +# The second field determines whether you want to run your program in +# debug or optimized mode. The latter is significantly faster, but no +# run-time checking of parameters and internal states is performed, so +# you should set this value to `on' while you develop your program, +# and to `off' when running production computations. +debug-mode = off + + +# As third field, we need to give the path to the top-level deal.II +# directory. You need to adjust this to your needs. Since this path is +# probably the most often needed one in the Makefile internals, it is +# designated by a single-character variable, since that can be +# reference using $D only, i.e. without the parentheses that are +# required for most other parameters, as e.g. in $(target). +D = ../deal.II + + +# The last field specifies the names of data and other files that +# shall be deleted when calling `make clean'. Object and backup files, +# executables and the like are removed anyway. Here, we give a list of +# files in the various output formats that deal.II supports. +clean-up-files = *gmv *gnuplot *gpl *eps *pov *vtk *ucd *.d2 + + + + +# +# +# Usually, you will not need to change anything beyond this point. +# +# +# The next statement tells the `make' program where to find the +# deal.II top level directory and to include the file with the global +# settings +include $D/common/Make.global_options + + +# Since the whole project consists of only one file, we need not +# consider difficult dependencies. We only have to declare the +# libraries which we want to link to the object file. deal.II has two +# libraries: one for the debug mode version of the +# application and one for optimized mode. +libs.g := $(lib-deal2.g) +libs.o := $(lib-deal2.o) + + +# We now use the variable defined above to switch between debug and +# optimized mode to select the set of libraries to link with. Included +# in the list of libraries is the name of the object file which we +# will produce from the single C++ file. Note that by default we use +# the extension .g.o for object files compiled in debug mode and .o for +# object files in optimized mode (or whatever local default on your +# system is instead of .o) +ifeq ($(debug-mode),on) + libraries = $(target).g.$(OBJEXT) $(libs.g) +else + libraries = $(target).$(OBJEXT) $(libs.o) +endif + + +# Now comes the first production rule: how to link the single object +# file produced from the single C++ file into the executable. Since +# this is the first rule in the Makefile, it is the one `make' selects +# if you call it without arguments. +all: $(target)$(EXEEXT) +$(target)$(EXEEXT) : $(libraries) + @echo ============================ Linking $@ + @$(CXX) -o $@ $^ $(LIBS) $(LDFLAGS) + + +# To make running the application somewhat independent of the actual +# program name, we usually declare a rule `run' which simply runs the +# program. You can then run it by typing `make run'. This is also +# useful if you want to call the executable with arguments which do +# not change frequently. You may then want to add them to the +# following rule: +run: $(target)$(EXEEXT) + @echo ============================ Running $< + @./$(target)$(EXEEXT) + + +# As a last rule to the `make' program, we define what to do when +# cleaning up a directory. This usually involves deleting object files +# and other automatically created files such as the executable itself, +# backup files, and data files. Since the latter are not usually quite +# diverse, you needed to declare them at the top of this file. +clean: + -rm -f *.$(OBJEXT) *~ Makefile.dep $(target)$(EXEEXT) $(clean-up-files) + + +# Since we have not yet stated how to make an object file from a C++ +# file, we should do so now. Since the many flags passed to the +# compiler are usually not of much interest, we suppress the actual +# command line using the `at' sign in the first column of the rules +# and write the string indicating what we do instead. +./%.g.$(OBJEXT) : + @echo "==============debug========= $( $@" + @$(CXX) $(CXXFLAGS.g) -c $< -o $@ +./%.$(OBJEXT) : + @echo "==============optimized===== $( $@" + @$(CXX) $(CXXFLAGS.o) -c $< -o $@ + + +# The following statement tells make that the rules `run' and `clean' +# are not expected to produce files of the same name as Makefile rules +# usually do. +.PHONY: all run clean + + +# Finally there is a rule which you normally need not care much about: +# since the executable depends on some include files from the library, +# besides the C++ application file of course, it is necessary to +# re-generate the executable when one of the files it depends on has +# changed. The following rule creates a dependency file +# `Makefile.dep', which `make' uses to determine when to regenerate +# the executable. This file is automagically remade whenever needed, +# i.e. whenever one of the cc-/h-files changed. Make detects whether +# to remake this file upon inclusion at the bottom of this file. +# +# If the creation of Makefile.dep fails, blow it away and fail +Makefile.dep: $(target).cc Makefile \ + $(shell echo $D/include/deal.II/*/*.h) + @echo ============================ Remaking $@ + @$D/common/scripts/make_dependencies $(INCLUDE) -B. $(target).cc \ + > $@ \ + || (rm -f $@ ; false) + @if test -s $@ ; then true ; else rm $@ ; false ; fi + +# To make the dependencies known to `make', we finally have to include +# them: +include Makefile.dep + + diff --git a/tests/benchmarks/test_poisson/poisson.cc b/tests/benchmarks/test_poisson/poisson.cc new file mode 100644 index 0000000000..90af0d8dc3 --- /dev/null +++ b/tests/benchmarks/test_poisson/poisson.cc @@ -0,0 +1,231 @@ + // Use a Helmholtz operator with constant coefficient on a hypercube mesh + // for which the assembly cost is strongly dominated by writing into the + // matrix (cell_similarity gives that essentially nothing needs to be done + // except for cell 0) + + +const unsigned int element_degree = 2; +const unsigned int dimension = 3; + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + + +using namespace dealii; + + +template +class HelmholtzProblem +{ +public: + HelmholtzProblem (const FiniteElement &fe); + void run (); + +private: + void setup_system (); + void assemble_system (); + void solve (); + + Triangulation triangulation; + const FiniteElement &fe; + DoFHandler dof_handler; + + ConstraintMatrix hanging_node_constraints; + + SparsityPattern sparsity_pattern; + SparseMatrix system_matrix; + Vector tri_sol, tri_rhs; + TimerOutput timer; +}; + + + +template +HelmholtzProblem::HelmholtzProblem (const FiniteElement &fe) : + fe (fe), + dof_handler (triangulation), + timer(std::cout, TimerOutput::summary, TimerOutput::wall_times) +{} + + + +template +void HelmholtzProblem::setup_system () +{ + timer.enter_subsection("setup mesh and matrix"); + + GridGenerator::hyper_cube (triangulation, -1, 1); + triangulation.refine_global(6); + dof_handler.distribute_dofs (fe); + std::cout << "Number of active cells: " + << triangulation.n_active_cells() + << std::endl + << "Number total degrees of freedom: " + << dof_handler.n_dofs() << std::endl; + + hanging_node_constraints.clear (); + IndexSet locally_relevant (dof_handler.locally_owned_dofs().size()); + DoFTools::extract_locally_relevant_dofs (dof_handler, locally_relevant); + hanging_node_constraints.reinit (locally_relevant); + DoFTools::make_hanging_node_constraints (dof_handler, + hanging_node_constraints); + VectorTools::interpolate_boundary_values (dof_handler, + 0, + ZeroFunction(), + hanging_node_constraints); + + { + CompressedSimpleSparsityPattern csp (dof_handler.n_dofs(), + dof_handler.n_dofs(), + locally_relevant); + DoFTools::make_sparsity_pattern (dof_handler, csp, + hanging_node_constraints, false); + sparsity_pattern.copy_from (csp); + } + system_matrix.reinit(sparsity_pattern); + tri_sol.reinit (dof_handler.n_dofs()); + tri_rhs.reinit (tri_sol); + + timer.leave_subsection(); +} + + +template +void HelmholtzProblem::assemble_system () +{ + timer.enter_subsection("write into matrix"); + + QGauss quadrature_formula(fe.degree+1); + + const unsigned int n_q_points = quadrature_formula.size(); + const unsigned int dofs_per_cell = fe.dofs_per_cell; + + FullMatrix cell_matrix (dofs_per_cell, dofs_per_cell); + Vector cell_rhs (dofs_per_cell); + + std::vector local_dof_indices (dofs_per_cell); + + FEValues fe_values (fe, quadrature_formula, + update_values | update_gradients | + update_JxW_values); + + typename DoFHandler::active_cell_iterator + cell = dof_handler.begin_active(), + endc = dof_handler.end(); + for (; cell!=endc; ++cell) + { + fe_values.reinit (cell); + if (fe_values.get_cell_similarity() != CellSimilarity::translation) + cell_matrix = 0; + + for (unsigned int q_point=0; q_pointget_dof_indices (local_dof_indices); + hanging_node_constraints.distribute_local_to_global (cell_matrix, + local_dof_indices, + system_matrix); + } + + timer.leave_subsection(); +} + + + +template +void HelmholtzProblem::solve () +{ + for (unsigned int i=0; i +void HelmholtzProblem::run () +{ + setup_system(); + assemble_system(); + solve(); +} + +int main (int argc, char **argv) +{ + + try + { + Utilities::System::MPI_InitFinalize mpi_initialization(argc, argv); + deallog.depth_console (0); + + FE_Q fe(element_degree); + HelmholtzProblem problem(fe); + problem.run(); + } + catch (std::exception &exc) + { + std::cerr << std::endl << std::endl + << "----------------------------------------------------" + << std::endl; + std::cerr << "Exception on processing: " << std::endl + << exc.what() << std::endl + << "Aborting!" << std::endl + << "----------------------------------------------------" + << std::endl; + return 1; + } + catch (...) + { + std::cerr << std::endl << std::endl + << "----------------------------------------------------" + << std::endl; + std::cerr << "Unknown exception!" << std::endl + << "Aborting!" << std::endl + << "----------------------------------------------------" + << std::endl; + return 1; + } + + return 0; +} + + + diff --git a/tests/benchmarks/testlist.sh b/tests/benchmarks/testlist.sh index 92675beef1..b8371151e1 100644 --- a/tests/benchmarks/testlist.sh +++ b/tests/benchmarks/testlist.sh @@ -1,2 +1,2 @@ #!/bin/bash -export TESTS="step-22 tablehandler test_assembly" +export TESTS="step-22 tablehandler test_assembly test_poisson"