From: Wolfgang Bangerth Date: Mon, 8 Mar 2004 22:24:59 +0000 (+0000) Subject: Detect and use METIS to partition triangulations. X-Git-Tag: v8.0.0~15665 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=1b012a2bec1f0e56c5674004f15d1db4806adeba;p=dealii.git Detect and use METIS to partition triangulations. git-svn-id: https://svn.dealii.org/trunk@8681 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/aclocal.m4 b/deal.II/aclocal.m4 index e63ca75f0f..44813481a5 100644 --- a/deal.II/aclocal.m4 +++ b/deal.II/aclocal.m4 @@ -3717,3 +3717,60 @@ AC_DEFUN(DEAL_II_CONFIGURE_PETSC, dnl ]) fi ]) + + + +dnl ------------------------------------------------------------ +dnl Check whether Metis is installed, and if so store the +dnl respective links +dnl +dnl Usage: DEAL_II_CONFIGURE_METIS +dnl +dnl ------------------------------------------------------------ +AC_DEFUN(DEAL_II_CONFIGURE_METIS, dnl +[ + dnl First check for the Metis directory + AC_MSG_CHECKING(for Metis library directory) + + AC_ARG_WITH(metis, + [ --with-metis=/path/to/metis Specify the path to the Metis installation, + of which the include and library directories + are subdirs; use this if you want to + override the METIS_DIR environment variable], + [ + USE_CONTRIB_METIS=yes + DEAL_II_METIS_DIR=$withval + AC_MSG_RESULT($DEAL_II_METIS_DIR) + + dnl Make sure that what was specified is actually correct + if test ! -d $DEAL_II_METIS_DIR/Lib ; then + AC_MSG_ERROR([Path to Metis specified with --with-metis does not + point to a complete Metis installation]) + fi + ], + [ + dnl Take something from the environment variables, if it is there + if test "x$METIS_DIR" != "x" ; then + USE_CONTRIB_METIS=yes + DEAL_II_METIS_DIR="$METIS_DIR" + AC_MSG_RESULT($DEAL_II_METIS_DIR) + + dnl Make sure that what this is actually correct + if test ! -d $DEAL_II_METIS_DIR/include \ + -o ! -d $DEAL_II_METIS_DIR/lib ; then + AC_MSG_ERROR([The path to Metis specified in the METIS_DIR + environment variable does not + point to a complete Metis installation]) + fi + else + USE_CONTRIB_METIS=no + DEAL_II_METIS_DIR="" + AC_MSG_RESULT(not found) + fi + ]) + if test "$USE_CONTRIB_METIS" = "yes" ; then + AC_DEFINE(DEAL_II_USE_METIS, 1, + [Defined if a Metis installation was found and is going + to be used]) + fi +]) diff --git a/deal.II/base/include/base/config.h.in b/deal.II/base/include/base/config.h.in index 865c05391f..cfcbb62c4d 100644 --- a/deal.II/base/include/base/config.h.in +++ b/deal.II/base/include/base/config.h.in @@ -164,6 +164,9 @@ provided by the compiler. */ #undef DEAL_II_USE_DIRECT_ERRNO_H +/* Defined if a Metis installation was found and is going to be used */ +#undef DEAL_II_USE_METIS + /* Flag indicating whether the library shall be compiled for multithreaded applications. If so, then it is set to one, otherwise to zero. */ #undef DEAL_II_USE_MT diff --git a/deal.II/configure b/deal.II/configure index cd87a4a653..0e4f390e0f 100755 --- a/deal.II/configure +++ b/deal.II/configure @@ -1,5 +1,5 @@ #! /bin/sh -# From configure.in Revision: 1.163 . +# From configure.in Revision. # Guess values for system-dependent variables and create Makefiles. # Generated by GNU Autoconf 2.57. # @@ -273,7 +273,7 @@ PACKAGE_BUGREPORT= ac_unique_file="deal.II" ac_subdirs_all="$ac_subdirs_all contrib tests" -ac_subst_vars='SHELL PATH_SEPARATOR PACKAGE_NAME PACKAGE_TARNAME PACKAGE_VERSION PACKAGE_STRING PACKAGE_BUGREPORT exec_prefix prefix program_transform_name bindir sbindir libexecdir datadir sysconfdir sharedstatedir localstatedir libdir includedir oldincludedir infodir mandir build_alias host_alias target_alias DEFS ECHO_C ECHO_N ECHO_T LIBS DEAL_II_VERSION DEAL_II_MAJOR DEAL_II_MINOR DEAL_II_PATH DEAL2_DIR build build_cpu build_vendor build_os host host_cpu host_vendor host_os target target_cpu target_vendor target_os CC CFLAGS LDFLAGS CPPFLAGS ac_ct_CC EXEEXT OBJEXT CXX CXXFLAGS ac_ct_CXX GXX_VERSION CXXFLAGSG CXXFLAGSO CXXFLAGSPIC SHLIBLD LDFLAGSPIC enablemultithreading withmultithreading F77 F77_VERSION F77FLAGSO F77FLAGSG F77FLAGSPIC F77LIBS enableshared lib_suffix AR RANLIB ac_ct_RANLIB enablemultigrid TECPLOT_LIBRARY_PATH TECPLOT_INCLUDE_PATH USE_CONTRIB_HSL USE_CONTRIB_PETSC DEAL_II_PETSC_DIR DEAL_II_PETSC_ARCH kdocdir kdocversion DOXYGEN PERL subdirs LIBOBJS LTLIBOBJS' +ac_subst_vars='SHELL PATH_SEPARATOR PACKAGE_NAME PACKAGE_TARNAME PACKAGE_VERSION PACKAGE_STRING PACKAGE_BUGREPORT exec_prefix prefix program_transform_name bindir sbindir libexecdir datadir sysconfdir sharedstatedir localstatedir libdir includedir oldincludedir infodir mandir build_alias host_alias target_alias DEFS ECHO_C ECHO_N ECHO_T LIBS DEAL_II_VERSION DEAL_II_MAJOR DEAL_II_MINOR DEAL_II_PATH DEAL2_DIR build build_cpu build_vendor build_os host host_cpu host_vendor host_os target target_cpu target_vendor target_os CC CFLAGS LDFLAGS CPPFLAGS ac_ct_CC EXEEXT OBJEXT CXX CXXFLAGS ac_ct_CXX GXX_VERSION CXXFLAGSG CXXFLAGSO CXXFLAGSPIC SHLIBLD LDFLAGSPIC enablemultithreading withmultithreading F77 F77_VERSION F77FLAGSO F77FLAGSG F77FLAGSPIC F77LIBS enableshared lib_suffix AR RANLIB ac_ct_RANLIB enablemultigrid TECPLOT_LIBRARY_PATH TECPLOT_INCLUDE_PATH USE_CONTRIB_HSL USE_CONTRIB_PETSC DEAL_II_PETSC_DIR DEAL_II_PETSC_ARCH USE_CONTRIB_METIS DEAL_II_METIS_DIR kdocdir kdocversion DOXYGEN PERL subdirs LIBOBJS LTLIBOBJS' ac_subst_files='' # Initialize some variables set by options. @@ -833,6 +833,10 @@ Optional Packages: installation; use this if you want to override the PETSC_ARCH environment variable + --with-metis=/path/to/metis Specify the path to the Metis installation, + of which the include and library directories + are subdirs; use this if you want to + override the METIS_DIR environment variable --with-cpu=cpu Optimize specifically for the given CPU type, rather than just generating code for this processor family @@ -7811,6 +7815,65 @@ fi; + echo "$as_me:$LINENO: checking for Metis library directory" >&5 +echo $ECHO_N "checking for Metis library directory... $ECHO_C" >&6 + + +# Check whether --with-metis or --without-metis was given. +if test "${with_metis+set}" = set; then + withval="$with_metis" + + USE_CONTRIB_METIS=yes + DEAL_II_METIS_DIR=$withval + echo "$as_me:$LINENO: result: $DEAL_II_METIS_DIR" >&5 +echo "${ECHO_T}$DEAL_II_METIS_DIR" >&6 + + if test ! -d $DEAL_II_METIS_DIR/Lib ; then + { { echo "$as_me:$LINENO: error: Path to Metis specified with --with-metis does not + point to a complete Metis installation" >&5 +echo "$as_me: error: Path to Metis specified with --with-metis does not + point to a complete Metis installation" >&2;} + { (exit 1); exit 1; }; } + fi + +else + + if test "x$METIS_DIR" != "x" ; then + USE_CONTRIB_METIS=yes + DEAL_II_METIS_DIR="$METIS_DIR" + echo "$as_me:$LINENO: result: $DEAL_II_METIS_DIR" >&5 +echo "${ECHO_T}$DEAL_II_METIS_DIR" >&6 + + if test ! -d $DEAL_II_METIS_DIR/include \ + -o ! -d $DEAL_II_METIS_DIR/lib ; then + { { echo "$as_me:$LINENO: error: The path to Metis specified in the METIS_DIR + environment variable does not + point to a complete Metis installation" >&5 +echo "$as_me: error: The path to Metis specified in the METIS_DIR + environment variable does not + point to a complete Metis installation" >&2;} + { (exit 1); exit 1; }; } + fi + else + USE_CONTRIB_METIS=no + DEAL_II_METIS_DIR="" + echo "$as_me:$LINENO: result: not found" >&5 +echo "${ECHO_T}not found" >&6 + fi + +fi; + if test "$USE_CONTRIB_METIS" = "yes" ; then + +cat >>confdefs.h <<\_ACEOF +#define DEAL_II_USE_METIS 1 +_ACEOF + + fi + + + + + echo "$as_me:$LINENO: result: " >&5 @@ -8765,6 +8828,8 @@ s,@USE_CONTRIB_HSL@,$USE_CONTRIB_HSL,;t t s,@USE_CONTRIB_PETSC@,$USE_CONTRIB_PETSC,;t t s,@DEAL_II_PETSC_DIR@,$DEAL_II_PETSC_DIR,;t t s,@DEAL_II_PETSC_ARCH@,$DEAL_II_PETSC_ARCH,;t t +s,@USE_CONTRIB_METIS@,$USE_CONTRIB_METIS,;t t +s,@DEAL_II_METIS_DIR@,$DEAL_II_METIS_DIR,;t t s,@kdocdir@,$kdocdir,;t t s,@kdocversion@,$kdocversion,;t t s,@DOXYGEN@,$DOXYGEN,;t t diff --git a/deal.II/configure.in b/deal.II/configure.in index 3ae1b2bb3b..6c479d0387 100644 --- a/deal.II/configure.in +++ b/deal.II/configure.in @@ -342,6 +342,10 @@ AC_SUBST(USE_CONTRIB_PETSC) AC_SUBST(DEAL_II_PETSC_DIR) AC_SUBST(DEAL_II_PETSC_ARCH) +DEAL_II_CONFIGURE_METIS +AC_SUBST(USE_CONTRIB_METIS) +AC_SUBST(DEAL_II_METIS_DIR) + dnl ------------------------------------------------------------- diff --git a/deal.II/deal.II/include/grid/grid_tools.h b/deal.II/deal.II/include/grid/grid_tools.h index 6623d7ae3f..8ffa04747c 100644 --- a/deal.II/deal.II/include/grid/grid_tools.h +++ b/deal.II/deal.II/include/grid/grid_tools.h @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 2001, 2002, 2003 by the deal.II authors +// Copyright (C) 2001, 2002, 2003, 2004 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -28,7 +28,7 @@ * cell that contains a given point. See the descriptions of the * individual functions for more information. * - * @author Wolfgang Bangerth, 2001, 2003 + * @author Wolfgang Bangerth, 2001, 2003, 2004 */ class GridTools { @@ -233,7 +233,38 @@ class GridTools typename Container::active_cell_iterator find_active_cell_around_point (const Container &container, const Point &p); - + + /** + * Use the METIS partitioner to generate + * a partitioning of the active cells + * making up the entire domain. After + * calling this function, the subdomain + * ids of all active cells will have + * values between zero and + * @p{n_partitions-1}. You can access the + * subdomain id of a cell by using + * @p{cell->subdomain_id()}. + * + * This function will generate an error + * if METIS is not installed unless + * @p{n_partitions} is one. I.e., you can + * write a program so that it runs in the + * single-processor single-partition case + * without METIS installed, and only + * requires METIS when multiple + * partitions are required. + */ + template + static + void partition_triangulation (Triangulation &triangulation, + const unsigned int n_partitions); + /** + * Exception + */ + DeclException1 (ExcInvalidNumberOfPartitions, + int, + << "The number of partitions you gave is " << arg1 + << ", but must be greater than zero."); /** * Exception */ diff --git a/deal.II/deal.II/source/grid/grid_tools.cc b/deal.II/deal.II/source/grid/grid_tools.cc index 874b445376..66b0ef18db 100644 --- a/deal.II/deal.II/source/grid/grid_tools.cc +++ b/deal.II/deal.II/source/grid/grid_tools.cc @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 2001, 2002, 2003 by the deal.II authors +// Copyright (C) 2001, 2002, 2003, 2004 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -17,8 +17,12 @@ #include #include #include +#include +#include #include #include +#include +#include #include #include @@ -238,6 +242,85 @@ GridTools::find_active_cell_around_point (const Container &container, +template +void +GridTools:: +partition_triangulation (Triangulation &triangulation, + const unsigned int n_partitions) +{ + Assert (n_partitions > 0, ExcInvalidNumberOfPartitions(n_partitions)); + + // check for an easy return + if (n_partitions == 1) + { + for (typename Triangulation::active_cell_iterator + cell = triangulation.begin_active(); + cell != triangulation.end(); ++cell) + cell->set_subdomain_id (0); + return; + } + + // we decompose the domain by first + // generating the connection graph of all + // cells with their neighbors, and then + // passing this graph off to METIS. To make + // things a little simpler and more + // general, we let the function + // DoFTools:make_flux_sparsity_pattern + // function generate the connection graph + // for us and reuse the SparsityPattern + // data structure for the connection + // graph. The connection structure of the + // mesh is obtained by using a fake + // piecewise constant finite element + // + // Since in 3d the generation of a + // sparsity pattern can be expensive, we + // take the detour of the compressed + // sparsity pattern, which is a little + // slower but more efficient in terms of + // memory + FE_DGQ fake_q0(0); + DoFHandler dof_handler (triangulation); + dof_handler.distribute_dofs (fake_q0); + Assert (dof_handler.n_dofs() == triangulation.n_active_cells(), + ExcInternalError()); + + CompressedSparsityPattern csp (dof_handler.n_dofs(), + dof_handler.n_dofs()); + DoFTools::make_flux_sparsity_pattern (dof_handler, csp); + + SparsityPattern sparsity_pattern; + sparsity_pattern.copy_from (csp); + + // partition this connection graph and get + // back a vector of indices, one per degree + // of freedom (which is associated with a + // cell) + std::vector partition_indices (triangulation.n_active_cells()); + sparsity_pattern.partition (n_partitions, partition_indices); + + // finally loop over all cells and set the + // subdomain ids. for this, get the DoF + // index of each cell and extract the + // subdomain id from the vector obtained + // above + std::vector dof_indices(1); + for (typename DoFHandler::active_cell_iterator + cell = dof_handler.begin_active(); + cell != dof_handler.end(); ++cell) + { + cell->get_dof_indices(dof_indices); + Assert (dof_indices[0] < triangulation.n_active_cells(), + ExcInternalError()); + Assert (partition_indices[dof_indices[0]] < n_partitions, + ExcInternalError()); + + cell->set_subdomain_id (partition_indices[dof_indices[0]]); + } +} + + #if deal_II_dimension != 1 template @@ -267,3 +350,8 @@ template MGDoFHandler::active_cell_iterator GridTools::find_active_cell_around_point (const MGDoFHandler &, const Point &p); + +template +void +GridTools::partition_triangulation (Triangulation &, + const unsigned int); diff --git a/deal.II/doc/development/Makefile b/deal.II/doc/development/Makefile index ff279d1f6a..1b86652764 100644 --- a/deal.II/doc/development/Makefile +++ b/deal.II/doc/development/Makefile @@ -27,6 +27,8 @@ makefiles.html: Makefile makefiles.1.html $D/common/Make.global_options makefile @echo '
  • USE_CONTRIB_PETSC=$(USE_CONTRIB_PETSC)' >> $@ @echo '
  • DEAL_II_PETSC_DIR=$(DEAL_II_PETSC_DIR)' >> $@ @echo '
  • DEAL_II_PETSC_ARCH=$(DEAL_II_PETSC_ARCH)' >> $@ + @echo '
  • USE_CONTRIB_METIS=$(USE_CONTRIB_METIS)' >> $@ + @echo '
  • DEAL_II_METIS_DIR=$(DEAL_II_METIS_DIR)' >> $@ @echo '
  • lib-path-base=$(lib-path-base)' >> $@ @echo '
  • lib-path-lac=$(lib-path-lac)' >> $@ @echo '
  • lib-path-deal2=$(lib-path-deal2)' >> $@ @@ -45,6 +47,7 @@ makefiles.html: Makefile makefiles.1.html $D/common/Make.global_options makefile @echo '
  • lib-contrib-hsl=$(lib-contrib-hsl)' >> $@ @echo '
  • lib-contrib-petsc.g=$(lib-contrib-petsc.g)' >> $@ @echo '
  • lib-contrib-petsc.o=$(lib-contrib-petsc.o)' >> $@ + @echo '
  • lib-contrib-metis=$(lib-contrib-metis)' >> $@ @echo '
  • include-path-base=$(include-path-base)' >> $@ @echo '
  • include-path-lac=$(include-path-lac)' >> $@ @echo '
  • include-path-deal2=$(include-path-deal2)' >> $@ @@ -52,6 +55,7 @@ makefiles.html: Makefile makefiles.1.html $D/common/Make.global_options makefile @echo '
  • include-path-contrib-petsc=$(include-path-contrib-petsc)' >> $@ @echo '
  • include-path-contrib-petsc-bmake=$(include-path-contrib-petsc-bmake)' >> $@ @echo '
  • include-path-contrib-petsc-mpi=$(include-path-contrib-petsc-mpi)' >> $@ + @echo '
  • include-path-contrib-metis=$(include-path-contrib-metis)' >> $@ @echo '
  • INCLUDE=$(INCLUDE)' >> $@ @echo '
  • CXXFLAGS.g=$(CXXFLAGS.g)' >> $@ @echo '
  • CXXFLAGS.o=$(CXXFLAGS.o)' >> $@ diff --git a/deal.II/doc/development/makefiles.1.html b/deal.II/doc/development/makefiles.1.html index 1dcb6c2b71..149b4f983f 100644 --- a/deal.II/doc/development/makefiles.1.html +++ b/deal.II/doc/development/makefiles.1.html @@ -210,6 +210,28 @@ upon configuration of the deal.II library.

    + + +
    USE_CONTRIB_METIS
    +

    + yes indicates whether the METIS library has been + detected upon configuration, and whether we shall interface with + it. no otherwise. (See the + ReadMe + file for more information on installation of METIS.) +

    +
    + + +
    DEAL_II_METIS_DIR
    +

    + If USE_CONTRIB_METIS is yes, then this + variable has the path to the METIS library, as either given by the + METIS_DIR environment variable, or the + --with-metis switch on the command line upon + configuration of the deal.II library. +

    +
    @@ -321,6 +343,16 @@ file for more information on installation of PETSc.)

    + +
    lib-contrib-metis
    +

    + Path and name of the METIS libraries, if + found and used at configuration time. + (See the + ReadMe + file for more information on installation of METIS.) +

    +
    @@ -349,6 +381,14 @@ +
    include-path-deal2
    +

    + Same for the /deal.II library. Include files are in + $(include-path-deal2)/*. +

    +
    + +
    include-path-contrib-hsl

    If some HSL files are @@ -377,11 +417,14 @@

    -
    include-path-deal2
    +
    include-path-contrib-metis +

    - Same for the /deal.II library. Include files are in - $(include-path-deal2)/*. -

    + If METIS is used, then this variable has the path to the METIS + include files. (See the + ReadMe + file for more information on installation of METIS.) +

    diff --git a/deal.II/lac/include/lac/sparsity_pattern.h b/deal.II/lac/include/lac/sparsity_pattern.h index b8ad6e0755..a92886edf3 100644 --- a/deal.II/lac/include/lac/sparsity_pattern.h +++ b/deal.II/lac/include/lac/sparsity_pattern.h @@ -731,6 +731,65 @@ class SparsityPattern : public Subscriptor */ bool optimize_diagonal () const; + /** + * Use the METIS partitioner to generate + * a partitioning of the degrees of + * freedom represented by this sparsity + * pattern. In effect, we view this + * sparsity pattern as a graph of + * connections between various degrees of + * freedom, where each nonzero entry in + * the sparsity pattern corresponds to an + * edge between two nodes in the + * connection graph. The goal is then to + * decompose this graph into groups of + * nodes so that a minimal number of + * edges are cut by the boundaries + * between node groups. This partitioning + * is done by METIS. Note that METIS can + * only partition symmetric sparsity + * patterns, and that of course the + * sparsity pattern has to be square. We + * do not check for symmetry of the + * sparsity pattern, since this is an + * expensive operation, but rather leave + * this as the responsibility of caller + * of this function. + * + * After calling this function, the + * output array will have values between + * zero and @p{n_partitions-1} for each + * node (i.e. row or column of the + * matrix). + * + * This function will generate an error + * if METIS is not installed unless + * @p{n_partitions} is one. I.e., you can + * write a program so that it runs in the + * single-processor single-partition case + * without METIS installed, and only + * requires METIS when multiple + * partitions are required. + * + * Note that the sparsity pattern itself + * is not changed by calling this + * function. However, you will likely use + * the information generated by calling + * this function to renumber degrees of + * freedom, after which you will of + * course have to regenerate the sparsity + * pattern. + * + * This function will rarely be called + * separately, since in finite element + * methods you will want to partition the + * mesh, not the matrix. This can be done + * by calling + * @p{GridTools::partition_triangulation}. + */ + void partition (const unsigned int n_partitions, + std::vector &partition_indices) const; + /** * Write the data of this object * en bloc to a file. This is @@ -928,7 +987,24 @@ class SparsityPattern : public Subscriptor int, int, << "The iterators denote a range of " << arg1 << " elements, but the given number of rows was " << arg2); - + /** + * Exception + */ + DeclException0 (ExcMETISNotInstalled); + /** + * Exception + */ + DeclException1 (ExcInvalidNumberOfPartitions, + int, + << "The number of partitions you gave is " << arg1 + << ", but must be greater than zero."); + /** + * Exception + */ + DeclException2 (ExcInvalidArraySize, + int, int, + << "The array has size " << arg1 << " but should have size " + << arg2); private: /** * Maximum number of rows that can diff --git a/deal.II/lac/source/sparsity_pattern.cc b/deal.II/lac/source/sparsity_pattern.cc index 74de399834..537617041e 100644 --- a/deal.II/lac/source/sparsity_pattern.cc +++ b/deal.II/lac/source/sparsity_pattern.cc @@ -23,6 +23,11 @@ #include #include +#ifdef DEAL_II_USE_METIS +extern "C" { +# include +} +#endif const unsigned int SparsityPattern::invalid_entry; @@ -899,6 +904,81 @@ SparsityPattern::memory_consumption () const +void +SparsityPattern:: +partition (const unsigned int n_partitions, + std::vector &partition_indices) const +{ + Assert (rows==cols, ExcNotQuadratic()); + Assert (is_compressed(), ExcNotCompressed()); + + Assert (n_partitions > 0, ExcInvalidNumberOfPartitions(n_partitions)); + Assert (partition_indices.size() == rows, + ExcInvalidArraySize (partition_indices.size(),rows)); + + // check for an easy return + if (n_partitions == 1) + { + std::fill_n (partition_indices.begin(), partition_indices.size(), 0U); + return; + } + + // Make sure that METIS is actually + // installed and detected +#ifndef DEAL_II_USE_METIS + AssertThrow (false, ExcMETISNotInstalled()); +#else + + // generate the data structures for + // METIS. Note that this is particularly + // simple, since METIS wants exactly our + // compressed row storage format. we only + // have to set up a few auxiliary arrays + int + n = static_cast(rows), + wgtflag = 0, // no weights on nodes or edges + numflag = 0, // C-style 0-based numbering + nparts = static_cast(n_partitions), // number of subdomains to create + dummy; // the numbers of edges cut by the + // resulting partition + + // use default options for METIS + int options[5] = { 0,0,0,0,0 }; + + // one more nuisance: we have to copy our + // own data to arrays that store signed + // integers :-( + std::vector int_rowstart (rowstart, rowstart+cols+1); + std::vector int_colnums (colnums, colnums+max_vec_len+1); + + std::vector int_partition_indices (rows); + + // Select which type of partitioning to create + + // Use recursive if the number of + // partitions is less than or equal to 8 + if (n_partitions <= 8) + METIS_PartGraphRecursive(&n, &int_rowstart[0], &int_colnums[0], + NULL, NULL, + &wgtflag, &numflag, &nparts, &options[0], + &dummy, &int_partition_indices[0]); + + // Otherwise use kway + else + METIS_PartGraphKway(&n, &int_rowstart[0], &int_colnums[0], + NULL, NULL, + &wgtflag, &numflag, &nparts, &options[0], + &dummy, &int_partition_indices[0]); + + // now copy back generated indices into the + // output array + std::copy (int_partition_indices.begin(), + int_partition_indices.end(), + partition_indices.begin()); +#endif +} + + // explicit instantiations template void SparsityPattern::copy_from (const FullMatrix &, bool); template void SparsityPattern::copy_from (const FullMatrix &, bool); diff --git a/tests/bits/metis_01.cc b/tests/bits/metis_01.cc new file mode 100644 index 0000000000..c2ab67ef42 --- /dev/null +++ b/tests/bits/metis_01.cc @@ -0,0 +1,92 @@ +//---------------------------- metis_01.cc --------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2004 by the deal.II authors +// +// This file is subject to QPL and may not be distributed +// without copyright and license information. Please refer +// to the file deal.II/doc/license.html for the text and +// further information on this license. +// +//---------------------------- metis_01.cc --------------------------- + + +// check GridTools::partition_triangulation + +#include "../tests.h" +#include +#include +#include +#include +#include +#include +#include + + +template +void test () +{ + Triangulation triangulation; + GridGenerator::hyper_cube (triangulation); + triangulation.refine_global (5-dim); + + // subdivide into 5 subdomains + deallog << "RECURSIVE" << std::endl; + GridTools::partition_triangulation (triangulation, 5); + for (typename Triangulation::active_cell_iterator + cell = triangulation.begin_active(); + cell != triangulation.end(); ++cell) + deallog << cell << ' ' << cell->subdomain_id() << std::endl; + + // subdivide into 9 subdomains. note that + // this uses the k-way subdivision (for + // more than 8 subdomains) rather than the + // recursive one + deallog << "K-WAY" << std::endl; + GridTools::partition_triangulation (triangulation, 9); + for (typename Triangulation::active_cell_iterator + cell = triangulation.begin_active(); + cell != triangulation.end(); ++cell) + deallog << cell << ' ' << cell->subdomain_id() << std::endl; +} + + + +int main () +{ + std::ofstream logfile("metis_01.output"); + deallog.attach(logfile); + deallog.depth_console(0); + + try + { + test<1> (); + test<2> (); + test<3> (); + } + 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; + }; +}