From: heltai Date: Thu, 5 Jun 2014 07:37:45 +0000 (+0000) Subject: Merged from trunk r32444 through r33019. X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=d05709ed11cc50dd6c53a0d15fbf870ece8e1e41;p=dealii-svn.git Merged from trunk r32444 through r33019. git-svn-id: https://svn.dealii.org/branches/branch_manifold_id_intermediate@33020 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/cmake/configure/configure_2_trilinos.cmake b/deal.II/cmake/configure/configure_2_trilinos.cmake index 7b680a4f5e..a79505384b 100644 --- a/deal.II/cmake/configure/configure_2_trilinos.cmake +++ b/deal.II/cmake/configure/configure_2_trilinos.cmake @@ -118,8 +118,7 @@ MACRO(FEATURE_TRILINOS_FIND_EXTERNAL var) ENDIF() # - # Trilinos has to be configured with 32bit indices if deal.II uses unsigned long - # long int. + # Trilinos has to be configured with 32bit indices if deal.II uses unsigned int. # IF(TRILINOS_WITH_NO_32BIT_INDICES AND NOT DEAL_II_WITH_64BIT_INDICES) MESSAGE(STATUS "deal.II was configured to use 32bit global indices but " diff --git a/deal.II/cmake/macros/macro_check_mpi_interface.cmake b/deal.II/cmake/macros/macro_check_mpi_interface.cmake index 0011282f69..7099bb2242 100644 --- a/deal.II/cmake/macros/macro_check_mpi_interface.cmake +++ b/deal.II/cmake/macros/macro_check_mpi_interface.cmake @@ -24,14 +24,26 @@ MACRO(CHECK_MPI_INTERFACE _feature _var) IF(DEAL_II_WITH_MPI) + SET(_nope FALSE) FOREACH(_library ${${_feature}_LIBRARIES}) - IF(_library MATCHES "/libmpi[^/]*\\.so") - LIST(FIND MPI_LIBRARIES "${_library}" _position) - IF("${_position}" STREQUAL "-1") + IF( _library MATCHES "/libmpi(|_cxx)\\.(a|so)[^/]*$") + + GET_FILENAME_COMPONENT(_file1 ${_library} REALPATH) + + SET(_not_found TRUE) + FOREACH(_mpi_library ${MPI_LIBRARIES}) + GET_FILENAME_COMPONENT(_file2 ${_mpi_library} REALPATH) + IF("${_file1}" STREQUAL "${_file2}") + SET(_not_found FALSE) + BREAK() + ENDIF() + ENDFOREACH() + + IF(_not_found) SET(_nope TRUE) - SET(_mpi_library ${_library}) + SET(_spurious_library ${_library}) BREAK() ENDIF() ENDIF() @@ -48,7 +60,7 @@ MACRO(CHECK_MPI_INTERFACE _feature _var) "Could not find a sufficient ${_feature} installation:\n" "${_feature} has to be compiled against the same MPI library as deal.II " "but the link line of ${_feature} contains:\n" - " ${_mpi_library}\n" + " ${_spurious_library}\n" "which is not listed in MPI_LIBRARIES:\n" " MPI_LIBRARIES = \"${_str}\"\n" ) diff --git a/deal.II/cmake/macros/macro_filter_system_libraries.cmake b/deal.II/cmake/macros/macro_filter_system_libraries.cmake index 49b4e75bf3..293edb030b 100644 --- a/deal.II/cmake/macros/macro_filter_system_libraries.cmake +++ b/deal.II/cmake/macros/macro_filter_system_libraries.cmake @@ -33,7 +33,7 @@ MACRO(FILTER_SYSTEM_LIBRARIES _feature) SET(_tmp_${_variable} ${${_variable}}) SET(${_variable} "") FOREACH(_lib ${_tmp_${_variable}}) - IF(_lib MATCHES "lib(c|quadmath|gfortran|m|rt|nsl|dl|pthread).so$") + IF(_lib MATCHES "lib(c|quadmath|gfortran|m|rt|nsl|dl|pthread)\\.(a|so)$") string(REGEX REPLACE ".*lib([a-z]+).so$" "\\1" _lib ${_lib}) ENDIF() LIST(APPEND ${_variable} ${_lib}) diff --git a/deal.II/cmake/modules/FindBOOST.cmake b/deal.II/cmake/modules/FindBOOST.cmake index b14b063f7c..33208996f5 100644 --- a/deal.II/cmake/modules/FindBOOST.cmake +++ b/deal.II/cmake/modules/FindBOOST.cmake @@ -42,7 +42,9 @@ IF(NOT BUILD_SHARED_LIBS) SET(Boost_USE_STATIC_LIBS TRUE) ENDIF() +LIST(REMOVE_ITEM CMAKE_MODULE_PATH ${CMAKE_SOURCE_DIR}/cmake/modules/) FIND_PACKAGE(Boost 1.44 COMPONENTS iostreams serialization system thread) +LIST(APPEND CMAKE_MODULE_PATH ${CMAKE_SOURCE_DIR}/cmake/modules/) # # Fall back to dynamic libraries if no static libraries could be found: diff --git a/deal.II/cmake/modules/FindMPI.cmake b/deal.II/cmake/modules/FindMPI.cmake index 111c75a58e..cf1671594c 100644 --- a/deal.II/cmake/modules/FindMPI.cmake +++ b/deal.II/cmake/modules/FindMPI.cmake @@ -152,12 +152,14 @@ DEAL_II_PACKAGE_HANDLE(MPI CXX_FLAGS OPTIONAL MPI_CXX_COMPILE_FLAGS LINKER_FLAGS OPTIONAL MPI_CXX_LINK_FLAGS CLEAR - MPI_MPI_H - MPI_HEADER_PATH - MPI_LIB - MPIEXEC - MPI_CXX_COMPILER MPI_C_COMPILER + MPI_CXX_COMPILER + MPIEXEC + MPI_EXTRA_LIBRARY MPI_Fortran_COMPILER + MPI_HEADER_PATH + MPI_LIB + MPI_LIBRARY + MPI_MPI_H ) diff --git a/deal.II/cmake/modules/FindPETSC.cmake b/deal.II/cmake/modules/FindPETSC.cmake index 4cc4aef519..fc0bb477d5 100644 --- a/deal.II/cmake/modules/FindPETSC.cmake +++ b/deal.II/cmake/modules/FindPETSC.cmake @@ -174,12 +174,13 @@ IF(NOT PETSC_PETSCVARIABLES MATCHES "-NOTFOUND") FILE(STRINGS "${PETSC_PETSCVARIABLES}" PETSC_EXTERNAL_LINK_LINE REGEX "^PETSC_WITH_EXTERNAL_LIB =.*") + SEPARATE_ARGUMENTS(PETSC_EXTERNAL_LINK_LINE) SET(_hints) SET(_petsc_libraries) SET(_cleanup_variables) - FOREACH(_token ${PETSC_EXTERNAL_LINK_LINE}}) + FOREACH(_token ${PETSC_EXTERNAL_LINK_LINE}) IF(_token MATCHES "^-L") # Build up hints with the help of all tokens passed with -L: STRING(REGEX REPLACE "^-L" "" _token "${_token}") @@ -188,21 +189,21 @@ IF(NOT PETSC_PETSCVARIABLES MATCHES "-NOTFOUND") # Search for every library that was specified with -l: STRING(REGEX REPLACE "^-l" "" _token "${_token}") - LIST(APPEND _cleanup_variables PETSC_LIBRARY_${_token}) + IF(NOT _token MATCHES "(petsc|stdc\\+\\+|gcc_s|clang_rt)") + LIST(APPEND _cleanup_variables PETSC_LIBRARY_${_token}) - IF(_token MATCHES "^(c|quadmath|gfortran|m|rt|nsl|dl|pthread)$") - FIND_SYSTEM_LIBRARY(PETSC_LIBRARY_${_token} NAMES ${_token}) - IF(NOT PETSC_LIBRARY_${_token} MATCHES "-NOTFOUND") - LIST(APPEND _petsc_libraries ${PETSC_LIBRARY_${_token}}) + IF(_token MATCHES "^(c|quadmath|gfortran|m|rt|nsl|dl|pthread)$") + FIND_SYSTEM_LIBRARY(PETSC_LIBRARY_${_token} NAMES ${_token}) + ELSE() + DEAL_II_FIND_LIBRARY(PETSC_LIBRARY_${_token} + NAMES ${_token} + HINTS ${_hints} + ) ENDIF() - ELSEIF(NOT _token MATCHES "(petsc|stdc\\+\\+|gcc_s)") - DEAL_II_FIND_LIBRARY(PETSC_LIBRARY_${_token} - NAMES ${_token} - HINTS ${_hints} - ) IF(NOT PETSC_LIBRARY_${_token} MATCHES "-NOTFOUND") LIST(APPEND _petsc_libraries ${PETSC_LIBRARY_${_token}}) ENDIF() + ENDIF() ENDIF() diff --git a/deal.II/cmake/scripts/run_testsuite.cmake b/deal.II/cmake/scripts/run_testsuite.cmake index 8d2f3a98ab..09d1ee4c0d 100644 --- a/deal.II/cmake/scripts/run_testsuite.cmake +++ b/deal.II/cmake/scripts/run_testsuite.cmake @@ -245,7 +245,7 @@ ENDIF() MESSAGE("-- CTEST_SITE: ${CTEST_SITE}") IF( "${TRACK}" STREQUAL "Regression Tests" - AND NOT CTEST_SITE MATCHES "tester" ) + AND NOT CTEST_SITE MATCHES "^(simserv04|tester)$" ) MESSAGE(FATAL_ERROR " I'm sorry ${CTEST_SITE}, I'm afraid I can't do that. The TRACK \"Regression Tests\" is not for you. @@ -352,25 +352,28 @@ ELSE() # Umm, no valid subversion info was found, try again with git-svn: FIND_PACKAGE(Git QUIET) IF(GIT_FOUND) - EXECUTE_PROCESS( - COMMAND ${GIT_EXECUTABLE} svn info - WORKING_DIRECTORY ${CTEST_SOURCE_DIRECTORY} - OUTPUT_VARIABLE _svn_WC_INFO - RESULT_VARIABLE _result - OUTPUT_STRIP_TRAILING_WHITESPACE - ) + GET_FILENAME_COMPONENT(_path "${CTEST_SOURCE_DIRECTORY}" PATH) + IF(EXISTS ${_path}/.git/svn) + EXECUTE_PROCESS( + COMMAND ${GIT_EXECUTABLE} svn info + WORKING_DIRECTORY ${CTEST_SOURCE_DIRECTORY} + OUTPUT_VARIABLE _svn_WC_INFO + RESULT_VARIABLE _result + OUTPUT_STRIP_TRAILING_WHITESPACE + ) - IF(${_result} EQUAL 0) - STRING(REGEX REPLACE "^(.*\n)?URL: ([^\n]+).*" - "\\2" _svn_WC_URL "${_svn_WC_INFO}") - STRING(REGEX REPLACE "^(.*\n)?Repository Root: ([^\n]+).*" - "\\2" _svn_WC_ROOT "${_svn_WC_INFO}") - STRING(REGEX REPLACE "^(.*\n)?Revision: ([^\n]+).*" - "\\2" _svn_WC_REVISION "${_svn_WC_INFO}") - STRING(REGEX REPLACE "^(.*\n)?Last Changed Date: ([^\n]+).*" - "\\2" _svn_WC_LAST_CHANGED_DATE "${_svn_WC_INFO}") - ELSE() - SET(_svn_WC_INFO) + IF(${_result} EQUAL 0) + STRING(REGEX REPLACE "^(.*\n)?URL: ([^\n]+).*" + "\\2" _svn_WC_URL "${_svn_WC_INFO}") + STRING(REGEX REPLACE "^(.*\n)?Repository Root: ([^\n]+).*" + "\\2" _svn_WC_ROOT "${_svn_WC_INFO}") + STRING(REGEX REPLACE "^(.*\n)?Revision: ([^\n]+).*" + "\\2" _svn_WC_REVISION "${_svn_WC_INFO}") + STRING(REGEX REPLACE "^(.*\n)?Last Changed Date: ([^\n]+).*" + "\\2" _svn_WC_LAST_CHANGED_DATE "${_svn_WC_INFO}") + ELSE() + SET(_svn_WC_INFO) + ENDIF() ENDIF() ENDIF() ENDIF() diff --git a/deal.II/doc/authors.html b/deal.II/doc/authors.html index 928ad2fa45..30c255ea90 100644 --- a/deal.II/doc/authors.html +++ b/deal.II/doc/authors.html @@ -4,7 +4,7 @@ deal.II Authors - + @@ -13,350 +13,10 @@

deal.II Authors

-

- deal.II is a global project. It is - administered by a group of principal - developers. Technical decisions are made by the principal - developers and a group of developers consisting of - long-term contributors with a global overview of the library. A - large number of contributors have provided substantial - patches over the years. -

-

- These three groups are listed below. -

- - -

Principal developers

- - - -

Developers

- - - - - -

Contributors

- -

Many people have contributed to deal.II over the years, some of them very - substantial parts of the library. Their work - is greatly appreciated: no open source project can survice without a - community. The following people contributed major parts - of the library (in alphabetical order), with many more that have - sent in fixes and small enhancements: -

- - - - -

- deal.II draws from some ideas which were first implemented in the - predecessor library, - DEAL, developed by Guido Kanschat, Franz-Theo - Suttmeier, and Roland Becker. It has also benefitted from - knowing about DiffPack, a library developed in Norway in the - 1990s. -

- -
-
- The deal.II Authors - $Date$ -
-
- - Valid HTML 4.01! - - Valid CSS! -
+ + This page is now part of the top-level project description. See + this page. + diff --git a/deal.II/doc/developers/cmake-internals.html b/deal.II/doc/developers/cmake-internals.html index dad7505cfb..30667346c4 100644 --- a/deal.II/doc/developers/cmake-internals.html +++ b/deal.II/doc/developers/cmake-internals.html @@ -706,7 +706,7 @@ EXPAND_INSTANTIATIONS(obj_dofs "${inst_in_files}")
- The deal.II Authors + The deal.II Authors $Date$
diff --git a/deal.II/doc/developers/indentation.html b/deal.II/doc/developers/indentation.html index 14ef04f532..4b58188a84 100644 --- a/deal.II/doc/developers/indentation.html +++ b/deal.II/doc/developers/indentation.html @@ -40,7 +40,7 @@
- The deal.II Authors + The deal.II Authors $Date$
diff --git a/deal.II/doc/developers/porting.html b/deal.II/doc/developers/porting.html index a2f4a9dc19..1819f1195b 100644 --- a/deal.II/doc/developers/porting.html +++ b/deal.II/doc/developers/porting.html @@ -157,7 +157,7 @@ $ cmake -DCMAKE_TOOLCHAIN_FILE=<...>/Toolchain.sample \
- The deal.II Authors + The deal.II Authors $Date$
diff --git a/deal.II/doc/developers/svn.html b/deal.II/doc/developers/svn.html index 7adc81d831..0576bcb4a3 100644 --- a/deal.II/doc/developers/svn.html +++ b/deal.II/doc/developers/svn.html @@ -85,7 +85,7 @@ $ svn merge ^/trunk
- The deal.II Authors + The deal.II Authors $Date$
diff --git a/deal.II/doc/developers/testsuite.html b/deal.II/doc/developers/testsuite.html index 815209d784..281f157a51 100644 --- a/deal.II/doc/developers/testsuite.html +++ b/deal.II/doc/developers/testsuite.html @@ -911,7 +911,7 @@ $ ctest -DCONFIG_FILE="[...]/Config.sample" [...]
- The deal.II Authors + The deal.II Authors $Date$
diff --git a/deal.II/doc/developers/toc.html b/deal.II/doc/developers/toc.html index 6ca72c7ed5..6fb2e9b667 100644 --- a/deal.II/doc/developers/toc.html +++ b/deal.II/doc/developers/toc.html @@ -96,7 +96,7 @@ have the following resources available:
- The deal.II Authors + The deal.II Authors $Date$
diff --git a/deal.II/doc/developers/writing-documentation.html b/deal.II/doc/developers/writing-documentation.html index 5a878e14d2..2f3e725428 100644 --- a/deal.II/doc/developers/writing-documentation.html +++ b/deal.II/doc/developers/writing-documentation.html @@ -456,7 +456,7 @@ void foobar ()
- The deal.II Authors + The deal.II Authors $Date$
diff --git a/deal.II/doc/documentation.html b/deal.II/doc/documentation.html index bede5f1563..d2e1a9e8a6 100644 --- a/deal.II/doc/documentation.html +++ b/deal.II/doc/documentation.html @@ -5,7 +5,7 @@ The deal.II Online Documentation - + @@ -59,7 +59,7 @@
- The deal.II Authors + The deal.II Authors $Date$
diff --git a/deal.II/doc/doxygen/tutorial/tutorial.h.in b/deal.II/doc/doxygen/tutorial/tutorial.h.in index 44fe64626e..ea3546ef73 100644 --- a/deal.II/doc/doxygen/tutorial/tutorial.h.in +++ b/deal.II/doc/doxygen/tutorial/tutorial.h.in @@ -1,7 +1,7 @@ // --------------------------------------------------------------------- // $Id$ // -// Copyright (C) 2005 - 2013 by the deal.II authors +// Copyright (C) 2005 - 2014 by the deal.II authors // // This file is part of the deal.II library. // @@ -69,12 +69,14 @@ @@MAP@@ * @enddot * - * Legend:
- * Green octagons: programs that show basic techniques;
- * Orange boxes: advanced techniques;
- * Yellow: applications in fluid dynamics;
- * Light blue: applications in solid mechanics;
- * Dark blue boxes: time dependent problems. + * Legend:
+ *
+ * Green: programs that show basic techniques;
+ * Orange: advanced techniques;
+ * Yellow: applications in fluid dynamics;
+ * Light blue: applications in solid mechanics;
+ * Dark blue: time dependent problems. + *
* * *

Tutorial programs listed by number

@@ -566,7 +568,6 @@ * Step-31, * Step-32, * Step-33, - * Step-40, * Step-42, * Step-43 * diff --git a/deal.II/doc/external-libs/arpack.html b/deal.II/doc/external-libs/arpack.html index 8d7182e6ac..34ea0caea7 100644 --- a/deal.II/doc/external-libs/arpack.html +++ b/deal.II/doc/external-libs/arpack.html @@ -102,7 +102,7 @@
- The deal.II Authors + The deal.II Authors $Date$
diff --git a/deal.II/doc/external-libs/mumps.html b/deal.II/doc/external-libs/mumps.html index 7bb7b236f2..9de62539d5 100644 --- a/deal.II/doc/external-libs/mumps.html +++ b/deal.II/doc/external-libs/mumps.html @@ -346,7 +346,7 @@
- The deal.II Authors + The deal.II Authors $Date$
diff --git a/deal.II/doc/external-libs/p4est.html b/deal.II/doc/external-libs/p4est.html index 767b54df50..26f1f315bf 100644 --- a/deal.II/doc/external-libs/p4est.html +++ b/deal.II/doc/external-libs/p4est.html @@ -57,7 +57,7 @@ cmake -DP4EST_DIR=/path/to/installation -DDEAL_II_WITH_P4EST=ON -DDEAL_II_WITH_M
- The deal.II Authors + The deal.II Authors $Date$
diff --git a/deal.II/doc/external-libs/petsc.html b/deal.II/doc/external-libs/petsc.html index 3230d590f5..73efed4803 100644 --- a/deal.II/doc/external-libs/petsc.html +++ b/deal.II/doc/external-libs/petsc.html @@ -149,7 +149,7 @@
- The deal.II Authors + The deal.II Authors $Date$
diff --git a/deal.II/doc/external-libs/slepc.html b/deal.II/doc/external-libs/slepc.html index 8de2413f8b..0ac53c2b5a 100644 --- a/deal.II/doc/external-libs/slepc.html +++ b/deal.II/doc/external-libs/slepc.html @@ -109,7 +109,7 @@
- The deal.II Authors + The deal.II Authors $Date$
diff --git a/deal.II/doc/external-libs/trilinos.html b/deal.II/doc/external-libs/trilinos.html index 1e3668e30a..5776eb809a 100644 --- a/deal.II/doc/external-libs/trilinos.html +++ b/deal.II/doc/external-libs/trilinos.html @@ -113,7 +113,7 @@
- The deal.II Authors + The deal.II Authors $Date$
diff --git a/deal.II/doc/license.html b/deal.II/doc/license.html index 13cab62b72..38e41ca003 100644 --- a/deal.II/doc/license.html +++ b/deal.II/doc/license.html @@ -4,7 +4,7 @@ The deal.II License - + @@ -24,7 +24,7 @@

However, deal.II is not in the public domain, it is property of and copyrighted by the - deal.II Authors, and + deal.II Authors, and there are restrictions on its use: The legally binding license is the GNU Lesser General Public License (LGPL) as published by the Free Software @@ -43,7 +43,7 @@ In this document, as well as in the program files, the terms authors of deal.II and deal.II authors generally refer to all - people mentioned at the deal.II Authors + people mentioned at the deal.II Authors page.

@@ -649,7 +649,7 @@
- The deal.II Authors + The deal.II Authors $Date$
diff --git a/deal.II/doc/mail.html b/deal.II/doc/mail.html index e41cfbe5fe..6bf3bf31e2 100644 --- a/deal.II/doc/mail.html +++ b/deal.II/doc/mail.html @@ -4,7 +4,7 @@ deal.II Communications - + @@ -87,7 +87,7 @@
- The deal.II Authors + The deal.II Authors $Date$
diff --git a/deal.II/doc/navbar.html b/deal.II/doc/navbar.html index af84e94087..875747a0d6 100644 --- a/deal.II/doc/navbar.html +++ b/deal.II/doc/navbar.html @@ -4,7 +4,7 @@ deal.II navigation bar - + @@ -38,7 +38,7 @@ Online Resources

FAQ
- News
+ News
Download
Mailing list
Wiki
@@ -51,7 +51,7 @@

Reports
Publications
- Authors
+ Authors
License
Contact

diff --git a/deal.II/doc/news/changes.h b/deal.II/doc/news/changes.h index 180f5f88d7..84a45349c6 100644 --- a/deal.II/doc/news/changes.h +++ b/deal.II/doc/news/changes.h @@ -148,16 +148,55 @@ inconvenience this causes.

Specific improvements

    -
  1. New: AnyData::try_read() is a function that allows users to check -whether an entry exists and get a pointer to it without throwing an -exception in case of failure. -
    -(Guido Kanschat, 2014/05/16) -
    2. - -
  2. New: The GMRES solver of deal.II can now write an estimate of +
  3. Bugfix: Filter libclang_rt* from the PETSc link line. +
    + (Matthias Maier, 2014/06/04) +
    4. + +
  4. Improved: ParameterHandler parsing added to DoFOutputOperator +
    + (Guido Kanschat, 2014/06/03) +
    5. + +
  5. Improved: ThetaTimestepping and other operators now use AnyData + to communicate with operators they rely on. This way, time step data + can be forwarded using the same mechanism, and no complicated back + access to member data is required + anymore. MeshWorker::IntegrationInfo uses AnyData as well. +
    + (Guido Kanschat, 2014/06/03) +
    6. + +
  6. Bugfix: CMake: Also clean CMAKE_MODULE_PATH prior to call to + FIND_PACKAGE(Boost) inside FindBOOST.cmake because apparently + "Boost" == "BOOST" for the Mac file system... +
    + (Matthias Maier, 2014/05/28) +
    7. + +
  7. Improved: CMake: Search results and error conditions for external + libraries are now much more verbose. Added an MPI sanity check. +
    + (Matthias Maier, 2014/05/28) +
    8. + +
  8. Fixed: The MatrixTools::apply_boundary_values() variant that works + on PETSc matrices could produce a deadlock in parallel if one processor + had no boundary values to apply. This is now fixed. +
    + (Michal Wichrowski, 2014/05/19) +
    9. + +
  9. New: AnyData::try_read() is a function that allows users to check + whether an entry exists and get a pointer to it without throwing an + exception in case of failure. +
    + (Guido Kanschat, 2014/05/16) +
    10. + +
  10. New: The GMRES solver of deal.II can now write an estimate of eigenvalues to the log file, in analogy to the CG solver. This is enabled - by the flag SolverGMRES<>::AdditionalData::compute_eigenvalues. + by the flag SolverGMRES::AdditionalData::compute_eigenvalues.
    (Martin Kronbichler, 2014/05/11)
    11. diff --git a/deal.II/doc/news/news.html b/deal.II/doc/news/news.html index 3caea173e7..06e1f09444 100755 --- a/deal.II/doc/news/news.html +++ b/deal.II/doc/news/news.html @@ -1,1004 +1,23 @@ + "http://www.w3.org/TR/html4/loose.dtd"> - - The deal.II news page + + The deal.II news page - + - - - - - + + + + + - - - - + - - - + - - - - -

    News on deal.II

    - -

    - This page lists news and announcements on - deal.II in reverse chronological order. If you - send us extensions to the library, please feel free to also send a - patch to this file such that we can announce your contributions in - your own words. -

    - -

    - Please refer to the archive of the deal.II - mailing list for additional news. -

    - - -

    -

    Recent changes
    -
    - Changes in the library since the last major release are - here. -
    -

    - - -
    -
    - 2013/12/24: Version 8.1 released -
    -
    - deal.II version 8.1 was released today. A full list of changes can be - found here. -
    - deal.II can be downloaded at - this site. -
    - -
    - 2013/12/18: step-26 for the heat equation -
    -
    - step-26 fills a long-standing gap: There was no tutorial - program solving the heat equation. There was also no tutorial - showing in relatively easy terms how to do adaptive meshes - in time dependent problems. This program fills both of these - needs. -
    - -
    - 2013/11/01: step-42 solves an elasto-plastic contact problem -
    -
    - step-42, written by Jörg Frohne, Timo Heister and Wolfgang - Bangerth, shows how to solve an elasto-plastic contact problem. - The program is an extension of step-41 to a much more complex - equation (nonlinear elasto-plasticity) and also demonstrates - how to compute solutions for such problems in parallel. -
    - -
    - 2013/10/01: step-51 implements an HDG method -
    -
    - step-51, contributed by Martin Kronbichler and Scott Miller, - implements a hybridized discontinuous Galerkin method for the scalar - convection-diffusion equation. The method uses standard DG terms for - describing the element interior but couples the local sub-problems - using face elements. The local contribution is eliminated from the - global linear system by a static condensation process. -
    - -
    - 2013/07/24: Version 8.0 released -
    -
    - deal.II version 8.0 was released today. This is a major new - release, featuring the conversion of the configuration and - build system to CMake, - support for parallel computations with more than 4 billion - unknowns, and many more. A full list of changes can be - found here. -
    - deal.II can be downloaded at - this site. -
    - -
    - 2013/06/05: 64-bit version of deal.II -
    -
    - Because the global dof index variables in deal.II were of type - unsigned int, it was not possible to create a dof_handler with more than - 4.3 billions of unknowns. The 64-bit version of deal.II uses unsigned long - long int and thus, allows for 18e18 unknowns. -
    - - -
    - 2013/03/07: Build system switched to CMake -
    -
    - After having used Autoconf as the configuration and build system - since the beginning of deal.II in 1997, we have now switched to - CMake. This implementation, - contributed by Matthias Maier of Heidelberg University, has taken - more than half a year, consists of more than 10,000 lines, and - touches every aspect of the installation process. It provides a more - modern interface to and more reliable way of interacting with other - packages. Furthermore, it no longer relies on shell scripts and - therefore integrates more easily with integrated development - environments as well as Microsoft Windows without Cygwin. -
    - -
    - 2013/02/27: Version 7.3 released -
    -
    - Version 7.3.0 was released today. This release provides one new - tutorial, work on parallel linear algebra, the new FE_Q_DG0 element, - and a large number of enhancements and bug fixes throughout the - entire library. A complete list of new features is - found here. -
    - -
    - 2012/09/03: Version 7.2 released -
    -
    - Version 7.2.0 was released today. This release provides four - new tutorial programs, a framework for matrix-free - computations, and a large number of enhancements and bug fixes - throughout the entire library. A complete list of new - features is found here. -
    - - -
    - 2012/06/06: step-43 — an efficient solver for - two-phase flow -
    -
    - The step-21 tutorial program was the first program in the - deal.II tutorial to demonstrate the solution of two-phase flow - problems, but it suffered from the fact that it was very - slow. step-43 demonstrates some of the techniques that can be - used to make solvers for such problems more efficient and - accurate. -
    - - -
    - 2012/06/03: New step-15 -
    -
    - Step-15 used to be a program that demonstrated techniques for solving - a rather quirky, nonlinear, 1d problem. It was an interesting problem, - but it served little purpose in teaching deal.II or in demonstrating - nonlinear solution techniques. It has now been replaced by a more - illustrative example, contributed in large part by Sven Wetterauer, in - which we solve the nonlinear elliptic minimal surface equation using a - Newton iteration. -
    - - -
    - 2012/05/05: Matrix-free implementation framework -
    -
    - A matrix-free framework has been contributed by Katharina - Kormann and Martin Kronbichler. The implementation is memory - efficient and includes parallelism via MPI for distributed - meshes, threads for shared-memory computations, and usage of - processors' vector units. The two new tutorial programs step-37 - (multigrid solver) and step-48 (nonlinear wave equation with - distributed grid) demonstrate the new features. -
    - - -
    - 2012/04/02: step-44 — a finite strain - elasticity solver -
    -
    - step-44, contributed by Andrew McBride and Jean-Paul - Pelteret, demonstrates one approach to modelling large - deformations of nearly-incompressible elastic solids. The - elastic response is governed by a non-linear, hyperelastic - free-energy function. The geometrical response is also - nonlinear, i.e., the program considers finite deformations. -
    - - -
    - 2012/04/02: step-41 solves the obstacle problem -
    -
    - step-41, a program contributed by Jörg Frohne, shows - how to solve the obstacle problem. This problem describes - the deformation of a body subject to an inequality - constraint that arises from the fact that the deformation is - constrained by contact of the body with some external - obstacle. The program uses the primal-dual active set method - to solve this problem. -
    - - -
    - 2011/10/09: Version 7.1 released -
    -
    - Version 7.1.0 was released today. This release provides two - new tutorial programs and many smaller enhancements - throughout the entire library. A complete list of new - features is found here. -
    - - -
    - 2011/10/06: step-32, a massively parallel thermal - convection solver -
    -
    - Long in the making, the step-32 tutorial program - demonstrates how to solve the Boussinesq equations - describing thermal convection on large numbers of - processors. -
    - - -
    - 2011/04/30: step-46 demonstrates coupling different models - on different parts of the domain -
    -
    - step-46 shows how to couple different models defined on subsets of - the domain, in this case Stokes flow around an elastic solid. The - trick here is that variables (here the flow velocity and pressure, - and the solid displacement) do not live on the entire domain, but - only on a part. The point of the program is how to represent this in - source code. -
    - - -
    - 2011/01/09: Version 7.0 released -
    -
    - Version 7.0.0 was released today. This release is a huge step - forward, with several significant new features such as support for - thousands or tens of thousands of processors. A complete list of new - features is found here. -
    - - -
    - 2011/01/02: step-38 shows how to solve partial differential - equations on curved manifolds -
    -
    - The new step-38 program, contributed by Andrea Bonito and M. Sebastian - Pauletti, shows how to discretize and solve partial differential - equations posed on curved manifolds embedded in higher dimensional - spaces. -
    - -
    - 2010/12/20: Graphical frontend to edit parameter files -
    -
    - Contributed by Martin Steigemann, deal.II now has Qt based graphical - user interface to edit input parameter files for - the ParameterHandler class. The executable is located in - lib/bin/dealii_parameter_gui. -
    - -
    - 2010/10/25: New directory and library structure -
    -
    - deal.II has had a somewhat quirky directory and library structure for - historical reasons, but this has now changed: All include and source - files are now under the top-level include/ - and source directories. Furthermore, we no longer build a - plethora of libraries but only libdeal_II.g.so (debug - version) and libdeal_II.so (optimized version). In - particular, we no longer build different versions of the library for - different space dimensions. -
    - -
    - 2010/10/23: Distributed mesh functionality available -
    -
    - After more than 2 years of work by Timo Heister, Martin Kronbichler, - and Wolfgang Bangerth, and with help by Carsten Burstedde, we have - merged a branch on which we have implemented the functionality - necessary to distribute meshes on hundreds, thousands, or more - processors. An overview of the framework - in which this is implemented is presented in a documentation module, - an accompanying paper, as well as in the new step-40 tutorial program. -
    - -
    - 2010/07/18: Version 6.3.1 released -
    -
    - Version 6.3.1 was released today. This release addresses a number of - problems in release 6.3.0 related to compatibility with certain - compilers, operating systems, and versions of external libraries. - A list of changes between - versions 6.3.0 and 6.3.1 can be found here. -
    - -
    - 2010/06/25: deal.II on Facebook -
    -
    - We built a facebook page titled deal.II today. Please, feel free - to join us there! -
    - -
    - 2010/06/25: Version 6.3 released -
    -
    - Version 6.3.0 was released today, containing several new - tutorial programs, a new threading framework and a framework - for assembling linear systems by just providing the local integrals. - A complete list of changes between - versions 6.2.0 and 6.3.0 can be found here. -
    - -
    - 2010/06/24: step-39 solves Poisson's equation with - the interior penalty discontinuous Galerkin method on - adaptively refined meshes with multigrid -
    -
    - Using the MeshWorker framework, this program highlights the - fact that loops over meshes are very generic and that - building matrices and level matrices and evaluating error - estimates are really the same task. -
    - -
    - 2010/06/24: The Meshworker framework has been released -
    -
    - Using its generic loop and the provided auxiliary classes, - this framework relieves the application programmer from coding - loops over cells, faces, and subfaces, such that focus can remain - on implementation of the local integrals. -
    - -
    - 2010/06/23: step-45 shows periodic boundary - conditions -
    -
    - Contributed by Markus Bürg, the new tutorial program step-45 - demonstrates an often requested feature: periodic boundary - conditions. -
    - - -
    - 2010/02/13: New multigrid, new step-16 -
    -
    - Over the last few months, the multigrid implementation has seen - significant rewrites, with much of the work done by Bärbel - Janssen. The goal — now achieved — was to finally fully - support multigrid also for continuous finite elements on adaptively - refined meshes (uniformly refined meshes and discontinuous elements - have worked for a long time). As part of this process, step-16 has - been rewritten and now solves the same problem step-6 solves, just - with a multigrid solver. -
    - - -
    - 2009/10/07: step-35 shows the solution of the - Navier-Stokes equations -
    -
    - Step-35, a tutorial program contributed by Abner - Salgado-Gonzalez, shows an implementation of a projection - scheme for the numerical solution of the Navier-Stokes - equations. -
    - - -
    - 2009/07/30: New tutorial program step-36 -
    -
    - The new step-36 tutorial program, contributed by Toby D. Young and - Wolfgang Bangerth, shows how to solve eigenvalue problems, using the - Schrödinger equation as an example. -
    - - -
    - 2009/06/05: Version 6.2.1 released -
    -
    - deal.II version 6.2.0 contained a couple of minor but - annoying bugs in its interaction with PETSc. The - documentation package also did not include any typeset - formulas. Release 6.2.1 corrects these mistakes. - A complete list of changes between - versions 6.2.0 and 6.2.1 can be found here. -
    - -
    - 2009/04/27: Version 6.2 released -
    -
    - Version 6.2.0 was released today! - A complete list of changes between - versions 6.1.0 and 6.2.0 can be found here. -
    - -
    - 2009/04/14: New tutorial program step-34 -
    -
    - step-34, written by Luca Heltai, demonstrates the ability of - deal.II to discretize and solve boundary integral equations - on curved manifolds. -
    - - -
    - 2008/12/05: Support for meshes embedded in - higher dimensions -
    -
    - Based on work by Cataldo Manigrasso and Luca Heltai - developed over the last few months on a branch, deal.II now - supports meshes that are embedded in higher space - dimensions. This allows, for example, to solve equations on - the surface of a sphere (e.g. to use boundary element - methods for exterior problems). -
    - - -
    - 2008/11/14: New tutorial program step-31 -
    -
    - The latest tutorial program, step-31, is contributed by - Martin Kronbichler and Wolfgang Bangerth and shows the - solution of the Boussinesq approximation to thermally driven - convection in an incompressible fluid. -
    - - -
    - 2008/09/10: New tutorial program step-28 -
    -
    - step-28, contributed by Yaqi Wang, illustrates the solution - of a coupled system of diffusion equations relevant to - nuclear reactor physics where we use different meshes for - different components of a vector-valued solution. -
    - - -
    - 2008/08/21: deal.II now exceeds - of 400,000 lines of code -
    -
    - Since today, deal.II has more than 400,000 - lines of code. The rate of growth continues to be close on - par with the 3000 lines of code per month that we have had - over many years already, with no signs of abating activity! -
    - - -
    - 2008/07/08: Support for anisotropic mesh refinement -
    -
    - Contributed by Tobias Leicht as well as Ralf Hartmann, deal.II - now supports anisotropic refinement of meshes. This functionality - is explained in the step-30 tutorial program. -
    - - -
    - 2008/05/28: Version 6.1 released -
    -
    - Version 6.1.0 was released today! - A complete list of changes between - versions 6.0.0 and 6.1.0 can be found here. -
    - -
    - 2008/05/23: step-33 illustrates the Euler equations -
    -
    - step-33, - a tutorial program written mostly by David Neckels, illustrates how - the Euler equations of compressible inviscid gas dynamics can be - solved. More generally, the program also provides a framework for - the solution of arbitrary hyperbolic, nonlinear conservation laws. -
    - -
    - 2008/03/28: step-22 demonstrates solving the Stokes - equations -
    -
    - The new step-22 - tutorial program, written by Martin Kronbichler and Wolfgang - Bangerth, shows how to solve, and in particular how to solve - efficiently, the Stokes equations that describe viscous - incompressible fluid flow. -
    - -
    - 2007/10/31: step-29 tutorial program added -
    -
    - Contributed by Moritz Allmaras, the step-29 tutorial program - explains how problems involving complex numbers can be treated in - deal.II. In the example, - the complex-valued solution to a Helmholtz equation is split into - real and imaginary parts and a vector-valued finite element field - is used to discretize them seperately. The program also demonstrates - how quantities derived from the finite element solution can be directly - included in the output by using the new DataPostprocessor - class. -
    - -
    - 2007/09/07: Version 6.0 released -
    -
    - Version 6.0.0 was released today! - A complete and long list of changes between - versions 5.2.0 and 6.0.0 can be found here. -
    - -
    - 2007/08/09: step-27 demonstrates hp technology -
    -
    - The new step-27 - tutorial program finally demonstrates how to use the - hp finite element technology that has been part of deal.II for some - 18 months now. It solves a Laplace equation with hp finite elements - and shows how to set finite element degrees, assemble matrices on - cells with varying polynomial degrees, and how to compute a simple - criterion for estimating the local Sobolev smoothness of a function. -
    - -
    - 2007/07/16: Wilkinson Prize -
    -
    - The principal authors of deal.II—Wolfgang - Bangerth, Ralf Hartmann and Guido Kanscht—have received the J. H. Wilkinson - Prize for Numerical Software for their creation of and work on - the deal.II library. The prize committee particularly noted the - reliability and usability of the software and the quality of the - documentation. -

    - We want to take the opportunity to thank our user community for their - support over the years that kept us working on this software. We - consider this award motivation to continue this work! -

    - -
    - 2007/02/09: Support for previously unorientable meshes -
    -
    - Up to now, in 3D only 'orientable' meshes could be used in - deal.II, where all lines are in standard orientation and the - faces can be either in standard orientation or with a - reversed normal. This has been generalized such that now all - hexahedral meshes can be used. -
    - -
    - 2006/10/30: step-21 is finished -
    -
    - Step-21 demonstrates the time dependent flow of a multiphase fluid - through a porous medium. It is an extension of step-20, using a mixed - formulation and block solvers. -
    - -
    - 2006/10/03: step-23/24/25 are finished -
    -
    - The step-23, step-24 and step-25 tutorial programs have been added, - demonstrating how to solve the time dependent wave equation - and a nonlinear variant of it, the sine-Gordon equation. -
    - -
    - 2006/08/24: deal.II is now part of the - SPEC CPU 2006 benchmark -
    -
    - A modified version of the deal.II library and the - step-14 tutorial program has been incorporated into the industry - standard SPEC CPU 2006 - benchmark testsuite that is used to determine the speed of - computers. It is part of the floating point testsuite CFP2006 as benchmark - 447.dealII. -

    - The inclusion of our software into an industry - standard benchmark testsuite will lead to better support of compilers - and machines for the kind of program we write, in particular for - unstructured hierarchical mesh computations. This will in turn - help make our codes run faster on future systems. -

    - -
    - 2006/2/12: step-20 is finished -
    -
    - The new step-20 tutorial program shows how to use Raviart-Thomas - elements to solve a mixed formulation of the Laplace equation. It - also demonstrates the use of block matrices and vectors, and how they - can be used to define more complex solvers and preconditioners - working on the Schur complement. -
    - -
    - 2006/1/3-5: First User Workshop -
    -
    - After 8 years of existence, it was time to have the first user - meeting. During a 3-day workshop held in Heidelberg, topics of - software development, use of deal.II for applications, and future - directions were discussed by about a dozen developers and users. -
    - -
    - 2005/12/20: hp functionality merged -
    -
    - After the code had been developed and maintained on a branch for 3 - years by Oliver Kayser-Herold and Wolfgang Bangerth, it has now been - merged into the main development, and will be available for the next - release. -
    - -
    - 2005/09/17: Version 5.2 released -
    -
    - Version 5.2.0 was released today! - A complete and long list of changes between - versions 5.1.0 and 5.2.0 can be found here. -
    - -
    - 2005/09/09: Two more tutorial programs -
    -
    - Two more tutorial programs, step-18 and step-19, are now - available. step-18 shows how to solve time-dependent but quasi-static - elasticity programs in parallel. step-19 demonstrates some data - handling and parameter file techniques. -
    - -
    - 2005/03/29: deal.II now has a Wiki -
    -
    - Hosted by Luca Heltai, there is now a wiki. See the link from the - deal.II front page. Users of - deal.II are encouraged to visit it, and, - as is usual with wikis, edit pages freely and post their - experiences and code snippets they want others to know about. -
    - -
    - 2005/03/07: Function parser library integrated -
    -
    - Luca Heltai has contributed code that integrates a third - party function - parser library) that takes a string describing a - mathematical function, and returns an object that can be used in - all places where function objects are used in deal.II. -
    - -
    - 2004/12/24: Version 5.1 released -
    -
    - Version 5.1.0 was released today! - A complete and long list of changes between - versions 5.0.0 and 5.1.0 can be found here. -
    - -
    - 2004/10/21: Qr elements of arbitrary order now - support hanging node constraints also in 3d -
    - -
    - Lagrange elements previously had hanging node constraint matrices - hard coded in 3d, which were only available for r=1,2. Thanks to - code written by Oliver Kayser-Herold, this information is now - computed on-the-fly, which allows the use of Lagrange elements of - arbitrary order also in 3d on meshes with hanging nodes. -
    - -
    - 2004/05/25: Version 5.0 released -
    -
    - Version 5.0.0 was released today! - A complete and long list of changes between - versions 4.0.0 and 5.0.0 can be found here. -
    - -
    - 2004/03/08: Interfaces to PETSc and METIS -
    -
    - In order to support parallel computations, we now have interfaces to - both the PETSc library for sequential and parallel linear algebra, - as well as to the METIS library for domain partitioning. -
    - -
    - 2004/03/08: MIPSpro compiler version 7.41 supported -
    -
    - deal.II now also supports compilation on SGIs by the MIPSpro - compiler version 7.41. -
    - -
    - 2003/06/06: Version 4.0 released -
    -
    - After one year, we made it and released Version 4.0.0 today! - A complete list of changes between - versions 3.4.0 and 4.0.0 can be found here. -
    - -
    - 2003/04/27: deal.II now exceeds - of 200,000 lines of code -
    -
    - Since today, deal.II has now more than 200,000 - lines of code. The rate of growth continues to be close to - 3000 lines of code per month. -
    - -
    - 2003/01/24: deal.II runs on AIX 5.1 -
    -
    - Using gcc3.2, deal.II can now also be compiled on AIX systems - with the PowerPC chip. -
    - -
    - 2002/12/18: deal.II runs on Mac OS X 10.2 -
    -
    - Some changes have been made to allow compiling and running deal.II on - Mac OS X 10.2 systems. -
    - -
    - 2002/12/05: Intel ICC 7.0 compiler supported -
    -
    - deal.II now also supports compilation by the recently released - Intel icc 7.0 compiler. -
    - -
    - 2002/06/28: deal.II 64-bit ready -
    -
    - deal.II also runs on 64-bit systems without changes, as - tests on a Sparc machine with Solaris 9 in 64-bit mode have - shown. We also have reports of people running deal.II on - IA64 systems. -
    - -
    - 2002/06/06: Version 3.4 released -
    -
    - Today, version 3.4 of the deal.II library - has been released. A complete list of changes between - versions 3.3.0 and 3.4.0 can be found here. -
    - -
    - 2002/05/22: Intel ECC Itanium compiler now supported -
    -
    - aclocal.m4 has been changed to detect the Intel ECC - compiler for Itanium systems. The ECC compiler can produce - optimized code for Intel's new IA-64 architecture. -
    - -
    - 2002/03/07: New web home www.dealii.org and new - email addresses -
    -
    - The deal.II homepage has been moved to www.dealii.org. Also the email addresses have been - changed. Though the old web home and the old addresses will - be supported for some time we recommend to adjust the - bookmarks and mail addresses accordingly. -
    - -
    - 2002/02/01: Version 3.3 released -
    -
    - Today, version 3.3 of the deal.II library - has been released. A complete list of changes between - versions 3.2.0 and 3.3.0 can be found here. -
    - -
    - 2002/01/31: deal.II now consists - of 160,000 lines of code -
    -
    - All .cc and .h files of - deal.II now sum up to more than 160,000(!) - lines of code. deal.II continues to grow - with a rate almost constant at 3000 lines of code per month. -
    - -
    - 2001/09/28: Version 3.2 released -
    -
    - After a long phase of implementation and testing, we today - released version 3.2 of the deal.II - library. There are many changes between version 3.1 and - release 3.2. A complete list of these changes can be found - here. -
    - -
    - 2001/07/06: Compaq's cxx compiler supported -
    -
    - The present development version now also supports - compilation by Compaq's cxx compiler. See the ReadMe - file for more information on using this compiler. -
    - -
    - 2001/06/28: deal.II 3.1.2 released -
    -
    - Version 3.1.2 is a bug-fix release based on 3.1. A list of - changes between the two versions can be found on - this - page. -
    - -
    - 2001/05/23: Intel ICC compiler supported -
    -
    - The present development version of the library is now fixed - so that it can also be compiled by - Intel's ICC C++ compiler for Linux. -
    - -
    - 2001/05/17: deal.II 3.1.1 released -
    -
    - Version 3.1.1 is a bug-fix release based on 3.1. A list of - changes between the two versions can be found on - this - page. -
    - -
    - 2001/02/26: Support for Subroutines from the Harwell Subroutine - Library (HSL) -
    -
    - The Harwell Subroutine Library (HSL) is a - library that, among much other functionality, offers some - sparse direct solvers. We have added support for an - interface to these solvers to the library. -
    - -
    - 2001/02/20: Version 3.1 released -
    -
    - After a long testing phase, we today released version 3.1 of - the deal.II library. - There are many changes between version 3.0 and release - 3.1. A complete list of these changes can be found - here. -
    - -
    - 2000/08/07: deal.II 3.0.1 -
    -
    - Version 3.0.1 is a bug-fix release based on 3.0. A list of - changes between the two versions can be found on - this - page. - (Announcement) -
    - -
    - 2000/04/13: deal.II 3.0 announced -
    -
    - Today, after 2½ years of development and testing, the - first public release of the deal.II - library was announced (full text). -
    -
    - - - - - - - + + This page is now part of the top-level project description. See + this page. + diff --git a/deal.II/doc/publications/index.html b/deal.II/doc/publications/index.html index f48418761d..4804c337e2 100644 --- a/deal.II/doc/publications/index.html +++ b/deal.II/doc/publications/index.html @@ -546,6 +546,13 @@ Lecture Notes in Computational Science and Engineering Volume 95, pp 57-73 , 2014. +
  11. J. Sheldon, S. Miller, J. Pitt +
    + Methodology for Comparing Coupling Algorithms for Fluid-Structure Interaction Problems +
    + World Journal of Mathematics, vol. 4, pp. 54-70, 2014. +
    12. +
  12. M. Stoll
    One-shot solution of a time-dependent time-periodic PDE-constrained optimization problem @@ -5458,7 +5465,7 @@ A Goal Oriented Software Library for Solving PDEs and Optimization Problems with
    - The deal.II Authors + The deal.II Authors
    diff --git a/deal.II/doc/readme.html b/deal.II/doc/readme.html index 28f2040c78..2ad389e13b 100644 --- a/deal.II/doc/readme.html +++ b/deal.II/doc/readme.html @@ -6,7 +6,7 @@ The deal.II Readme - + @@ -798,7 +798,8 @@ in order to determine the directions of future work on the library. We are also very much interested in incorporating any work and bug fixes by third parties into the library. Of course, - you will be credited for this on the authors page. + you will be credited for this on + the authors page.

    @@ -842,7 +843,7 @@ the license of the deal.II project granted to you.

  13. - We, the deal.II Authors, do not require + We, the deal.II Authors, do not require copyright assignments for contributions. This means that the copyright for code contributions in the deal.II project is held by its respective contributors who have each agreed to release their @@ -874,7 +875,7 @@
    - The deal.II Authors + The deal.II Authors $Date$
    diff --git a/deal.II/doc/reports/index.html b/deal.II/doc/reports/index.html index 40cfbda80c..bed84536cc 100644 --- a/deal.II/doc/reports/index.html +++ b/deal.II/doc/reports/index.html @@ -80,7 +80,7 @@
    - The deal.II Authors + The deal.II Authors $Date$
    diff --git a/deal.II/doc/users/cmake.html b/deal.II/doc/users/cmake.html index 7bc712f4f9..97f706bfca 100644 --- a/deal.II/doc/users/cmake.html +++ b/deal.II/doc/users/cmake.html @@ -1033,7 +1033,7 @@ make deal_II.g VERBOSE=ON
    - The deal.II Authors + The deal.II Authors $Date$
    diff --git a/deal.II/doc/users/cmakelists.html b/deal.II/doc/users/cmakelists.html index 42623597ed..c423b104f8 100644 --- a/deal.II/doc/users/cmakelists.html +++ b/deal.II/doc/users/cmakelists.html @@ -838,7 +838,7 @@ D=/path/install/dir
    - The deal.II Authors + The deal.II Authors $Date$
    diff --git a/deal.II/doc/users/doxygen.html b/deal.II/doc/users/doxygen.html index b346d84320..a7c2b4e04d 100644 --- a/deal.II/doc/users/doxygen.html +++ b/deal.II/doc/users/doxygen.html @@ -52,7 +52,7 @@ TAGFILES = deal.tag=http://www.dealii.org/X.Y.Z/doxygen/deal.II
    - The deal.II Authors + The deal.II Authors $Date$
    diff --git a/deal.II/doc/users/navbar.html b/deal.II/doc/users/navbar.html index 822aefa4d8..7df7379a6b 100644 --- a/deal.II/doc/users/navbar.html +++ b/deal.II/doc/users/navbar.html @@ -6,7 +6,7 @@ deal.II Development Navigation Bar - + @@ -51,7 +51,7 @@ Online Resources

    FAQ
    - News
    + News
    Download
    Mailing list
    Wiki
    diff --git a/deal.II/doc/users/toc.html b/deal.II/doc/users/toc.html index bcb950a371..a8ab44758d 100644 --- a/deal.II/doc/users/toc.html +++ b/deal.II/doc/users/toc.html @@ -119,7 +119,7 @@ this.

    - The deal.II Authors + The deal.II Authors $Date$
    diff --git a/deal.II/examples/CMakeLists.txt b/deal.II/examples/CMakeLists.txt index ba2e146d3b..422dbc6608 100644 --- a/deal.II/examples/CMakeLists.txt +++ b/deal.II/examples/CMakeLists.txt @@ -26,7 +26,6 @@ IF(DEAL_II_COMPONENT_EXAMPLES) # PATTERN "CMakeFiles*" EXCLUDE PATTERN "doc*" EXCLUDE - PATTERN "doxygen*" EXCLUDE PATTERN ".svn*" EXCLUDE # # Glob Includes: @@ -35,6 +34,7 @@ IF(DEAL_II_COMPONENT_EXAMPLES) PATTERN "*.prm" PATTERN "*.inp" PATTERN "step*/CMakeLists.txt" + PATTERN "doxygen/CMakeLists.txt" # # Special files: # diff --git a/deal.II/examples/doxygen/CMakeLists.txt b/deal.II/examples/doxygen/CMakeLists.txt new file mode 100644 index 0000000000..5621f8129f --- /dev/null +++ b/deal.II/examples/doxygen/CMakeLists.txt @@ -0,0 +1,29 @@ +## +# CMake script for compiling the doxygen examples +## + +# Usually, you will not need to modify anything beyond this point... + +FILE(GLOB sources ${CMAKE_CURRENT_SOURCE_DIR}/*.cc) + +CMAKE_MINIMUM_REQUIRED(VERSION 2.8.8) + +FIND_PACKAGE(deal.II 8.0 QUIET + HINTS ${deal.II_DIR} ${DEAL_II_DIR} ../ ../../ $ENV{DEAL_II_DIR} + ) +IF(NOT ${deal.II_FOUND}) + MESSAGE(FATAL_ERROR "\n" + "*** Could not locate deal.II. ***\n\n" + "You may want to either pass a flag -DDEAL_II_DIR=/path/to/deal.II to cmake\n" + "or set an environment variable \"DEAL_II_DIR\" that contains this path." + ) +ENDIF() + +DEAL_II_INITIALIZE_CACHED_VARIABLES() +PROJECT(doxygen_examples) + +FOREACH(ccfile ${sources}) + GET_FILENAME_COMPONENT(file ${ccfile} NAME_WE) + ADD_EXECUTABLE(${file} ${ccfile}) + DEAL_II_SETUP_TARGET(${file}) +ENDFOREACH() diff --git a/deal.II/examples/doxygen/theta_timestepping.cc b/deal.II/examples/doxygen/theta_timestepping.cc index 94871bd102..caf390eba9 100644 --- a/deal.II/examples/doxygen/theta_timestepping.cc +++ b/deal.II/examples/doxygen/theta_timestepping.cc @@ -30,33 +30,25 @@ using namespace dealii; using namespace Algorithms; -class Explicit - : public Operator > +class Explicit : public Operator > { public: Explicit(const FullMatrix &matrix); - void operator() (NamedData*> &out, - const NamedData*> &in); + void operator() (AnyData &out, const AnyData &in); - void initialize_timestep_data(const TimestepData &); private: - const TimestepData *timestep_data; SmartPointer, Explicit> matrix; FullMatrix m; }; -class Implicit - : public Operator > +class Implicit : public Operator > { public: Implicit(const FullMatrix &matrix); - void operator() (NamedData*> &out, - const NamedData*> &in); - - void initialize_timestep_data(const TimestepData &); -private: - const TimestepData *timestep_data; + void operator() (AnyData &out, const AnyData &in); + + private: SmartPointer, Implicit> matrix; FullMatrix m; }; @@ -66,10 +58,10 @@ private: int main() { FullMatrix matrix(2); - matrix(0,0) = 1.; - matrix(1,1) = 1.; - matrix(0,1) = 31.4; - matrix(1,0) = -31.4; + matrix(0,0) = 0.; + matrix(1,1) = 0.; + matrix(0,1) = 3.14; + matrix(1,0) = -3.14; OutputOperator > out; out.initialize_stream(std::cout); @@ -77,16 +69,13 @@ int main() Explicit op_explicit(matrix); Implicit op_implicit(matrix); ThetaTimestepping > solver(op_explicit, op_implicit); - op_explicit.initialize_timestep_data(solver.explicit_data()); - op_implicit.initialize_timestep_data(solver.implicit_data()); solver.set_output(out); Vector value(2); value(0) = 1.; - NamedData*> indata; - NamedData*> outdata; - Vector *p = &value; - outdata.add(p, "value"); + AnyData indata; + AnyData outdata; + outdata.add(&value, "value"); solver.notify(Events::initial); solver(outdata, indata); @@ -102,24 +91,18 @@ Explicit::Explicit(const FullMatrix &M) void -Explicit::initialize_timestep_data(const TimestepData &t) -{ - timestep_data = &t; -} - - -void -Explicit::operator() (NamedData*> &out, const NamedData*> &in) +Explicit::operator() (AnyData &out, const AnyData &in) { + const double timestep = *in.read_ptr("Timestep"); if (this->notifications.test(Events::initial) || this->notifications.test(Events::new_timestep_size)) { - m.equ(-timestep_data->step, *matrix); + m.equ(-timestep, *matrix); for (unsigned int i=0; inotifications.clear(); - unsigned int i = in.find("Previous iterate"); - m.vmult(*out(0), *in(i)); + m.vmult(*out.entry*>(0), + *in.read_ptr >("Previous iterate")); } @@ -132,26 +115,19 @@ Implicit::Implicit(const FullMatrix &M) void -Implicit::initialize_timestep_data(const TimestepData &t) -{ - timestep_data = &t; -} - - -void -Implicit::operator() (NamedData*> &out, const NamedData*> &in) +Implicit::operator() (AnyData &out, const AnyData &in) { + const double timestep = *in.read_ptr("Timestep"); if (this->notifications.test(Events::initial) || this->notifications.test(Events::new_timestep_size)) { - m.equ(timestep_data->step, *matrix); + m.equ(timestep, *matrix); for (unsigned int i=0; inotifications.clear(); - - unsigned int i = in.find("Previous time"); - m.vmult(*out(0), *in(i)); + m.vmult(*out.entry*>(0), + *in.read_ptr >("Previous time")); } diff --git a/deal.II/examples/step-16/step-16.cc b/deal.II/examples/step-16/step-16.cc index e8a2035f0b..32aa7d2bde 100644 --- a/deal.II/examples/step-16/step-16.cc +++ b/deal.II/examples/step-16/step-16.cc @@ -1,7 +1,7 @@ /* --------------------------------------------------------------------- * $Id$ * - * Copyright (C) 2003 - 2013 by the deal.II authors + * Copyright (C) 2003 - 2014 by the deal.II authors * * This file is part of the deal.II library. * @@ -245,21 +245,32 @@ namespace Step16 // The following function extends what the corresponding one in step-6 // did, with the exception of also distributing the degrees of freedom on // each level of the mesh which is needed for the multigrid hierarchy. + // We then output the number of degrees (globally and on each level) and + // reset the solution and right hand side vectors to the correct size. + // As in other programs where we use a sequence of meshes that become + // finer and finer (starting with step-5), the code as is here sets + // the solution vector to zero, rather than trying to interpolate the + // solution obtained on the previous mesh to the current one. This + // follows from the empirical observation that interpolating the solution + // and using it as a starting guess does not significant decrease the number + // of linear CG iterations one has to perform to solve a linear system. + // (The situation is very different for nonlinear problems, where using + // a previous time step's solution, or a solution from a coarser grid, + // typically drastically makes the solution of a nonlinear system cheaper.) template void LaplaceProblem::setup_system () { dof_handler.distribute_dofs (fe); dof_handler.distribute_mg_dofs (fe); - // Here we output not only the degrees of freedom on the finest level, but - // also in the multilevel structure - deallog << "Number of degrees of freedom: " - << dof_handler.n_dofs(); - - for (unsigned int l=0; lREADME file. -There is an additional caveat: PETSc appears not to co-operate well when using -threads and the program crashes when deal.II is compiled in its usual mode -supporting multithreading. The program can therefore only be compiled and run -if multithreading is disabled upon configuration of deal.II using the ---disable-threads switch to ./configure. - Now, for the details: as mentioned, the program does not compute anything new, so the use of finite element classes etc. is exactly the same as before. The diff --git a/deal.II/examples/step-18/doc/intro.dox b/deal.II/examples/step-18/doc/intro.dox index 5e494a505a..a750e9e10d 100644 --- a/deal.II/examples/step-18/doc/intro.dox +++ b/deal.II/examples/step-18/doc/intro.dox @@ -534,9 +534,4 @@ the results section will show what displacement pattern comes out of this simulation. -@note PETSc appears not to co-operate well when using -threads and the program crashes when deal.II is compiled in its usual mode -supporting multithreading. The program can therefore only be compiled and run -if multithreading is disabled upon configuration of deal.II using the ---disable-threads switch to ./configure. diff --git a/deal.II/examples/step-32/doc/intro.dox b/deal.II/examples/step-32/doc/intro.dox index 6a20022310..46a21c1b61 100644 --- a/deal.II/examples/step-32/doc/intro.dox +++ b/deal.II/examples/step-32/doc/intro.dox @@ -827,10 +827,8 @@ processors busy and you won't get any additional speedup from using threads. For example, if your cluster nodes have 8 cores as they often have at the time of writing this, and if your batch scheduler puts 8 MPI processes on each node, then using threads doesn't make the program any -faster. Consequently, you probably want to use the ---disable-threads flag when configuring your deal.II installation -for this machine. On the other hand, if you want to run this program on a -single multicore machine, then it may make sense to use threads. That said, at +faster. Consequently, you probably want to either configure your deal.II without threads, or set the number of threads in MPI_InitFinalize to 1 (third argument), or "export DEAL_II_NUM_THREADS=1" before running. +That said, at the time of writing this, we only use the WorkStream class for assembling (parts of) linear systems, while 75% or more of the run time of the program is spent in the linear solvers that are not parallelized — in other words, diff --git a/deal.II/examples/step-40/doc/intro.dox b/deal.II/examples/step-40/doc/intro.dox index 0fec405ec8..dcc07c9b98 100644 --- a/deal.II/examples/step-40/doc/intro.dox +++ b/deal.II/examples/step-40/doc/intro.dox @@ -20,12 +20,7 @@ to work properly, this program needs access to the Hypre preconditioner package implementing algebraic multigrid; it can be installed as part of PETSc but has to be explicitly enabled during PETSc configuration; see the page linked to from the installation -instructions for PETSc. Finally, note that -PETSc appears not to co-operate well when using -threads and the program crashes when deal.II is compiled in its usual mode -supporting multithreading. The program can therefore only be compiled and run -if multithreading is disabled upon configuration of deal.II using the ---disable-threads switch to ./configure. +instructions for PETSc. diff --git a/deal.II/examples/step-47/doc/kind b/deal.II/examples/step-47/doc/kind index 15a13db451..c1d9154931 100644 --- a/deal.II/examples/step-47/doc/kind +++ b/deal.II/examples/step-47/doc/kind @@ -1 +1 @@ -basic +techniques diff --git a/deal.II/examples/step-52/doc/intro.dox b/deal.II/examples/step-52/doc/intro.dox index 4ae88a593d..926c81cdb4 100644 --- a/deal.II/examples/step-52/doc/intro.dox +++ b/deal.II/examples/step-52/doc/intro.dox @@ -8,9 +8,9 @@ This program shows how to use Runge-Kutta methods to solve a time-dependent problem. -

    Problem statement

    +

    Problem statement

    -In this example, we solve the energy-integrated time-dependent diffusion +In this example, we solve the one-group time-dependent diffusion approximation of the neutron transport equation (see step-28 for the time-independent multigroup diffusion). We assume that the medium is not fissible and therefore, the neutron flux satisfies the following equation: @@ -19,30 +19,28 @@ fissible and therefore, the neutron flux satisfies the following equation: - \Sigma_a(x) \phi(x,t) + S(x,t) @f} augmented by appropriate boundary conditions. Here, $v$ is the velocity of -neutrons, $D$ is the diffusion coefficient, $\Sigma_a$ is the absorption -cross section, and $S$ is a source. Because we are only interested in the -time dependence, we assume that $D$ and $\Sigma_a$ are constant. In this -example, we are only interested in the error in time and thus, we are looking -for a solution of the form: +neutrons (for simplicity we assume it is equal to 1), $D$ is the diffusion coefficient, +$\Sigma_a$ is the absorption cross section, and $S$ is a source. Because we are +only interested in the time dependence, we assume that $D$ and $\Sigma_a$ are +constant. We are looking for a solution on a square domain $[0,b]\times[0,b]$ of +the form: @f{eqnarray*} \phi(x,t) = A\sin(\omega t)(bx-x^2). @f} -By using quadratic finite elements, we will not have any spatial error. We -impose the following boundary conditions: homogeneous Dirichlet fo $x=0$ and +By using quadratic finite elements, all the error will be due to the time discretization. We +impose the following boundary conditions: homogeneous Dirichlet for $x=0$ and $x=b$ and homogeneous Neumann conditions for $y=0$ and $y=b$. The source is given by: @f{eqnarray*} S=A\left(\frac{1}{v}\omega \cos(\omega t)(bx -x^2) + \sin(\omega t) \left(\Sigma_a (bx-x^2)+2D\right) \right). @f} -Because the solution is a sine, we know that -\f$\phi\left(x,\frac{\pi}{\omega}\right) = 0$. Therefore, we can easily -compute the error at this time since it is simply the norm of the solution -found. +Because the solution is a sine, we know that $\phi\left(x,\pi\right) = 0$. +Therefore, the error at this time is simply the norm of the numerical solution.

    Runge-Kutta

    -The Runke-Kutta methods implemented in deal.II assume that the equation to be +The Runge-Kutta methods implemented in deal.II assume that the equation to be solved can be written as: @f{eqnarray*} \frac{dy}{dt} = f(t,y). @@ -61,11 +59,11 @@ y_{n+1} = y_n + \sum_{i=1}^s b_i k_i @f} where @f{eqnarray*} -k_i = h M^{-1} f(t_n+c_ih,y_n+\sum_{j=1}^sa_{ij}k_j) +k_i = h M^{-1} f\left(t_n+c_ih,y_n+\sum_{j=1}^sa_{ij}k_j\right) @f} -with $a_{ij}$, $b_i$, and $c_i$ are known coefficient and $h$ is the time step -used. The methods currently implemented in deal.II can be divided in three -categories: +where $a_{ij}$, $b_i$, and $c_i$ are known coefficients and $h$ is the time step +used. At the time of the writing of this tutorial, the methods implemented in +deal.II can be divided in three categories:
    1. explicit Runge-Kutta
    2. embedded (or adaptive) Runge-Kutta @@ -73,28 +71,39 @@ categories:

    Explicit Runge-Kutta

    -These methods that include for forward Euler, third order Runge-Kutta, and -fourth order Runge-Kutta, require a function to evaluate $M^{-1}f(t,y). These +These methods that include forward Euler, third order Runge-Kutta, and +fourth order Runge-Kutta, require a function to evaluate $M^{-1}f(t,y)$. These methods become unstable when the time step chosen is too large.

    Embedded Runge-Kutta

    These methods include Heun-Euler, Bogacki-Shampine, Dormand-Prince (ode45 in Matlab), Fehlberg, and Cash-Karp. These methods use a low order method to -estimate the error and decide if the time step needs to be refined or it can be -coarsen. Only embedded explicit methods have been implemented so far. +estimate the error and decide if the time step needs to be refined or coarsen. +Only embedded explicit methods have been implemented at the time of the writing.

    Implicit Runge-Kutta

    -These methods include backward Euler, implicit midpoint, Crank-Nicolson, and the +These methods include backward Euler, implicit midpoint, Crank-Nicolson, and a two stages SDIRK. These methods require to evaluate $M^{-1}f(t,y)$ and -$\left(I-\Delta t M^{-1} \frac{\partial f}{\partial Y}\right) = \left(M - \Delta -t \frac{\partial f}{\partial y}\right)^{-1} M$. These methods are always stable. +$\left(I-\Delta t M^{-1} \frac{\partial f}{\partial y}\right)^{-1}$ or equivalently +$\left(M - \Delta t \frac{\partial f}{\partial y}\right)^{-1} M$. These methods are +always stable. + +

    Weak form

    +To use the Runge-Kutta methods, we need to be able to evaluate: +@f{eqnarray*} +f = \oint D b_i \frac{\partial b_j}{\partial n} d\boldsymbol{r} - \int D \nabla +b_i \nabla b_j \phi_j d\boldsymbol{r} -\int \Sigma_a b_i b_j \phi_j +d\boldsymbol{r} + \int b_j S d\boldsymbol{r} +@f} +and $\frac{\partial f}{\partial y}$. Because $f$ is linear in $y$ (or $\phi$ in +this case) $\frac{\partial f}{\partial y} y = f$.

    Remarks

    -To simplify the problem, we solve the domain in two dimensional and the mesh is +To simplify the problem, the domain is two dimensional and the mesh is uniform (there is no need to adapt the mesh since we use quadratic finite elements and the exact solution is quadratic). Going from a two dimensional domain to a three dimensional domain is not very challenging. However if the -mesh must be adapted, we cannot forget to: +mesh must be adapted, it is important to remember to:
    1. project the solution to the new mesh when the mesh is changed. The mesh used should be the same at the beginning and at the end of the time step. diff --git a/deal.II/examples/step-52/doc/results.dox b/deal.II/examples/step-52/doc/results.dox index 4a70bbe4bc..cc3aeb9f40 100644 --- a/deal.II/examples/step-52/doc/results.dox +++ b/deal.II/examples/step-52/doc/results.dox @@ -1,6 +1,6 @@

      Results

      -The output of this program consist of the console output and solutions given in +The output of this program consists of the console output and the solutions given in vtu format. The console output is: @@ -24,7 +24,7 @@ Cash-Karp error: 2.5579e-10 Number of steps done: 200 @endcode -Like expected the high-order methods give a more accurate solutions. We see that +Like expected the high-order methods give a more accurate solution. We see that the Heun-Euler method adapted the number of time steps in order to satisfy the tolerance. The others embedded methods did not need to change the number of time steps. diff --git a/deal.II/examples/step-52/step-52.cc b/deal.II/examples/step-52/step-52.cc index 1031e49cdd..f75fa19b4c 100644 --- a/deal.II/examples/step-52/step-52.cc +++ b/deal.II/examples/step-52/step-52.cc @@ -20,7 +20,7 @@ // @sect3{Include files} -// The first task as usal is to include the functionality of these well-known +// The first task as usual is to include the functionality of these well-known // deal.II library files and some C++ header files. #include #include @@ -70,6 +70,7 @@ namespace Step52 public: Diffusion(); + // This function is the driver that will call the other functions. void run(); private: @@ -80,33 +81,36 @@ namespace Step52 // the time. void assemble_system(); - // Compute the intensity of the source at the given point. + // Compute the intensity of the source at a given time for a given point. double get_source(double time,const Point<2> &point) const; - // Evaluate the diffusion equation \f$M^{-1}(f(t,y))\f$ + // Evaluate the diffusion equation $M^{-1}(f(t,y))$ at a given time and + // for a given y. Vector evaluate_diffusion(const double time, const Vector &y) const; - // Evaluate \f$\left(I-\tau M^{-1} \frac{\partial f(t,y)}{\partial y}\right)^{-1} = - // \left(M-\tau \frac{\partial f}{\partial y}\right)^{-1} M \f$ - Vector id_minus_tau_J_inverse(const double time, const double tau, - const Vector &y); + // Evaluate $\left(I-\tau M^{-1} \frac{\partial f(t,y)}{\partial y}\right)^{-1}$ or + // equivalently $\left(M-\tau \frac{\partial f}{\partial y}\right)^{-1} M$ at a given + // time, for a given $\tau$ and y. + Vector id_minus_tau_J_inverse(const double time, + const double tau, + const Vector &y); - // Output the results as vtu + // Output the results as vtu files. void output_results(unsigned int time_step,TimeStepping::runge_kutta_method method) const; - // Driver for the explicit methods + // Driver for the explicit methods. void explicit_method(TimeStepping::runge_kutta_method method, const unsigned int n_time_steps, const double initial_time, const double final_time); - // Driver for the implicit methods + // Driver for the implicit methods. void implicit_method(TimeStepping::runge_kutta_method method, const unsigned int n_time_steps, const double initial_time, const double final_time); - // Driver for the embedded explicit methods. Returns the number of steps + // Driver for the embedded explicit methods. This function returns the number of steps // executed. unsigned int embedded_explicit_method(TimeStepping::runge_kutta_method method, const unsigned int n_time_steps, @@ -114,6 +118,7 @@ namespace Step52 const double final_time); + // The next parameters are self-explanatory. unsigned int fe_degree; double diffusion_coefficient; @@ -140,7 +145,7 @@ namespace Step52 - // We choose quadratic finite elements so that there are no spatial error. + // We choose quadratic finite elements and we initialize the parameters. Diffusion::Diffusion() : fe_degree(2), @@ -152,6 +157,7 @@ namespace Step52 + // @sect5{Diffusion::setup_system} void Diffusion::setup_system() { dof_handler.distribute_dofs(fe); @@ -165,6 +171,7 @@ namespace Step52 DoFTools::make_sparsity_pattern(dof_handler,c_sparsity,constraint_matrix); sparsity_pattern.copy_from(c_sparsity); + // Initialize the matrices and the solution vector. system_matrix.reinit(sparsity_pattern); mass_matrix.reinit(sparsity_pattern); mass_minus_tau_Jacobian.reinit(sparsity_pattern); @@ -173,6 +180,7 @@ namespace Step52 + // @sect5{Diffusion::assemble_system} void Diffusion::assemble_system() { system_matrix = 0.; @@ -196,7 +204,8 @@ namespace Step52 cell = dof_handler.begin_active(), endc = dof_handler.end(); - // Compute \f$-\int D \nabla b \cdot \nabla b - \int \Sigma_a b b\f$ and \f$\int b b\f$ + // Compute $-\int D \nabla b_i \cdot \nabla b_j d\boldsymbol{r} - \int \Sigma_a b_i b_j d\boldsymbol{r}$ + // and the mass matrix $\int b_i b_j d\boldsymbol{r}$. for (; cell!=endc; ++cell) { cell_matrix = 0.; @@ -229,6 +238,9 @@ namespace Step52 + // @sect5{Diffusion::get_source} + // + // Compute the source for a given time and a given point. double Diffusion::get_source(double time,const Point<2> &point) const { const double pi = 3.14159265358979323846; @@ -246,15 +258,16 @@ namespace Step52 + // @sect5{Diffusion:evaluate_diffusion} + // + // Evaluate the weak form of the diffusion equation at a given time t and for a given vector y. Vector Diffusion::evaluate_diffusion(const double time, const Vector &y) const { Vector tmp(dof_handler.n_dofs()); tmp = 0.; - // Compute system_matrix*y + // Compute $tmp=system\_matrix\cdot y$. system_matrix.vmult(tmp,y); - - // Compute the source term const QGauss<2> quadrature_formula(fe_degree+1); FEValues<2> fe_values(fe, quadrature_formula, @@ -300,7 +313,9 @@ namespace Step52 } - + // @sect5{Diffusion::id_minus_tau_J_inverse} + // + // We compute $\left(M-\tau \frac{\partial f}{\partial y}\right)^{-1} M$. Vector Diffusion::id_minus_tau_J_inverse(const double time, const double tau, const Vector &y) { @@ -308,10 +323,17 @@ namespace Step52 Vector result(y); SparseDirectUMFPACK inverse_mass_minus_tau_Jacobian; + // Compute $M-\tau \frac{\partial f}{\partial y}$. mass_minus_tau_Jacobian.copy_from(mass_matrix); mass_minus_tau_Jacobian.add(-tau,system_matrix); + + // Inverse the matrix to get $\left(M-\tau \frac{\partial f}{\partial y}\right)^{-1}$. inverse_mass_minus_tau_Jacobian.initialize(mass_minus_tau_Jacobian); + + // Compute $tmp=My$. mass_matrix.vmult(tmp,y); + + // Compute $\left(M-\tau \frac{\partial f}{\partial y}\right)^{-1} tmp = \left(M-\tau \frac{\partial f}{\partial y}\right)^{-1} My$. inverse_mass_minus_tau_Jacobian.vmult(result,tmp); return result; @@ -319,6 +341,9 @@ namespace Step52 + // @sect5{Diffusion::output_results} + // + // We output the solution in vtu files. void Diffusion::output_results(unsigned int time_step,TimeStepping::runge_kutta_method method) const { std::string method_name; @@ -401,7 +426,7 @@ namespace Step52 } - + // @sect5{Diffusion::explicit_method} void Diffusion::explicit_method(TimeStepping::runge_kutta_method method, const unsigned int n_time_steps, const double initial_time, @@ -415,7 +440,7 @@ namespace Step52 output_results(0,method); for (unsigned int i=0; iDiffusion::implicit_method} void Diffusion::implicit_method(TimeStepping::runge_kutta_method method, const unsigned int n_time_steps, const double initial_time, @@ -442,7 +468,7 @@ namespace Step52 output_results(0,method); for (unsigned int i=0; iDiffusion::embedded_explicit_method} unsigned int Diffusion::embedded_explicit_method(TimeStepping::runge_kutta_method method, const unsigned int n_time_steps, const double initial_time, @@ -487,11 +514,11 @@ namespace Step52 unsigned int n_steps=0; while (timefinal_time) time_step = final_time-time; - // Because we use a member function, we need to bind this to the + // Because we use a member function, we need to bind $this$ to the // function. time = embedded_explicit_runge_kutta.evolve_one_time_step( std_cxx1x::bind(&Diffusion::evaluate_diffusion,this,std_cxx1x::_1,std_cxx1x::_2), @@ -501,7 +528,7 @@ namespace Step52 if ((n_steps+1)%10==0) output_results(n_steps+1,method); - // Update the time step + // Update the time step. time_step = embedded_explicit_runge_kutta.get_status().delta_t_guess; ++n_steps; } @@ -511,14 +538,15 @@ namespace Step52 + // @sect5{Diffusion::run} void Diffusion::run() { - // Create the grid (a square [0,5]x[0,5]) and refine the mesh four times. - // The final gird has 16 times 16 cells, for a total of 256. + // Create the grid (a [0,5]x[0,5] square) and refine the mesh four times. + // The final grid has 16 by 16 cells, for a total of 256. GridGenerator::hyper_cube(triangulation, 0., 5.); triangulation.refine_global(4); - // Set the boundary indicator for x=0 and x=5 to 1 + // Set the boundary indicator for x=0 and x=5 to 1. typename Triangulation<2>::active_cell_iterator cell = triangulation.begin_active(), endc = triangulation.end(); @@ -542,56 +570,58 @@ namespace Step52 const double initial_time = 0.; const double final_time = 10.; - // Use forward Euler + // Use forward Euler. explicit_method(TimeStepping::FORWARD_EULER,n_time_steps,initial_time,final_time); std::cout<<"Forward Euler error: "<main() function} +// // The following main function is similar to previous examples as // well, and need not be commented on. int main () diff --git a/deal.II/include/deal.II/algorithms/any_data.h b/deal.II/include/deal.II/algorithms/any_data.h index 8845a5071e..9d67972c95 100644 --- a/deal.II/include/deal.II/algorithms/any_data.h +++ b/deal.II/include/deal.II/algorithms/any_data.h @@ -36,116 +36,171 @@ DEAL_II_NAMESPACE_OPEN * for long enough, then change to a map. */ class AnyData : - public Subscriptor + public Subscriptor { - public: - /// Default constructor for empty object - AnyData(); - - /// Number of stored data objects. - unsigned int size() const; - - /// Add a new data object - template - void add(type entry, const std::string& name); - - /** - * @brief Merge the data of another NamedData to the end of this object. - */ - void merge(const AnyData& other); - - /** - * @brief Access to stored data object by name. - * - * Find the object with given name, try to convert it to - * type and return it. This function throws an exception - * if either the name does not exist or if the conversion - * fails. If such an exception is not desired, use try_read() - * instead. - */ - template - type entry (const std::string& name); - - /** - * @brief Read-only access to stored data object by name. - * - * Find the object with given name, try to convert it to - * type and return it. This function throws an exception - * if either the name does not exist or if the conversion - * fails. If such an exception is not desired, use try_read() - * instead. - */ - template - const type entry (const std::string& name) const; - - /** - * @brief Dedicated read only access by name. - * - * For a constant object, this function equals entry(). For a - * non-const object, it forces read only access to the data. In - * particular, it throws an exception if the object is not found - * or cannot be converted to type. If such an exception is not - * desired, use try_read() instead. - */ - template - const type read (const std::string& name) const; - - /** - * @brief Dedicated read only access by name without exceptions. - * - * This function tries to find the name in the list and return a - * pointer to the associated object. If either the name is not - * found or the object cannot be converted to the return type, a - * null pointer is returned. - */ - template - const type * try_read (const std::string& name) const; - +public: + /// Default constructor for empty object + AnyData(); + + /// Number of stored data objects. + unsigned int size() const; + + /// Add a new data object + template + void add(type entry, const std::string &name); + + /** + * @brief Merge the data of another NamedData to the end of this object. + */ + void merge(const AnyData &other); + + /** + * @brief Access to stored data object by name. + * + * Find the object with given name, try to convert it to + * type and return it. This function throws an exception + * if either the name does not exist or if the conversion + * fails. If such an exception is not desired, use try_read() + * instead. + */ + template + type entry (const std::string &name); + + /** + * @brief Read-only access to stored data object by name. + * + * Find the object with given name, try to convert it to + * type and return it. This function throws an exception + * if either the name does not exist or if the conversion + * fails. If such an exception is not desired, use try_read() + * instead. + */ + template + const type entry (const std::string &name) const; + + /** + * @brief Dedicated read only access by name. + * + * For a constant object, this function equals entry(). For a + * non-const object, it forces read only access to the data. In + * particular, it throws an exception if the object is not found + * or cannot be converted to type. If such an exception is not + * desired, use try_read() instead. + * + * @warning Do not use this function for stored objects which are + * pointers. Use read_ptr() instead! + */ + template + const type read (const std::string &name) const; + + /** + * @brief Dedicated read only access by name for pointer data. + * + * If the stored data object is a pointer to a constant object, the logic + * of access becomes fairly complicated. Namely, the standard read + * function may fail, depending on whether it was a const pointer + * or a regular pointer. This function fixes the logic and + * ascertains that the object does not become mutable by accident. + */ + template + const type *read_ptr (const std::string &name) const; + + /** + * Perform the same action as read_ptr(), but do not throw an + * exception if the pointer does not exist. Return a null pointer + * instead. + */ + template + const type *try_read_ptr (const std::string &name) const; + + /** + * @brief Dedicated read only access by name without exceptions. + * + * This function tries to find the name in the list and return a + * pointer to the associated object. If either the name is not + * found or the object cannot be converted to the return type, a + * null pointer is returned. + */ + template + const type *try_read (const std::string &name) const; + /** * Access to stored data object by index. */ - template - type entry (const unsigned int i); - - /// Read-only access to stored data object by index. - template - const type entry (const unsigned int i) const; - - /// Dedicated read only access. - template - const type read (const unsigned int i) const; - - /// Name of object at index. - const std::string &name(const unsigned int i) const; - - /** - * @brief Find index of a named object - * - * Try to find the objecty and return its index in the list. Throw - * an exception if the object has not been found. - */ - unsigned int find(const std::string &name) const; - - /// Find out if object is of a certain type - template - bool is_type(const unsigned int i) const; - - /// Conversion from old NamedData - template - AnyData(const NamedData&); - - /// The requested type and the stored type are different + template + type entry (const unsigned int i); + + /// Read-only access to stored data object by index. + template + const type entry (const unsigned int i) const; + + /// Dedicated read only access. + template + const type read (const unsigned int i) const; + + /// Dedicated read only access to pointer object. + template + const type *read_ptr (const unsigned int i) const; + + /// Dedicated read only access to pointer object without exception. + template + const type *try_read_ptr (const unsigned int i) const; + + /// Dedicated read only access without exception. + template + const type *try_read (const unsigned int i) const; + + /// Name of object at index. + const std::string &name(const unsigned int i) const; + + /** + * @brief Find index of a named object + * + * Try to find the objecty and return its index in the list. Throw + * an exception if the object has not been found. + */ + unsigned int find(const std::string &name) const; + + /** + * @brief Try to find index of a named object + * + * Try to find the objecty and return its index in the + * list. returns numbers::invalid_unsigned_int if the name was not + * found. + */ + unsigned int try_find(const std::string &name) const; + + /// Find out if object is of a certain type + template + bool is_type(const unsigned int i) const; + + /// Conversion from old NamedData + template + AnyData(const NamedData &); + + /// An entry with this name does not exist in the AnyData object. + DeclException1(ExcNameNotFound, std::string &, + << "No entry with the name " << arg1 << " exists"); + + /// The requested type and the stored type are different DeclException2(ExcTypeMismatch, - char*, char*, - << "The requested type " << arg1 - << " and the stored type " << arg2 - << " must coincide"); - - private: - /// The stored data - std::vector data; - /// The names of the stored data - std::vector names; + char *, char *, + << "The requested type " << arg1 + << " and the stored type " << arg2 + << " must coincide"); + + /** + * Exception indicating that a function expected a vector to have a + * certain name, but NamedData had a different name in that + * position. + */ + DeclException2(ExcNameMismatch, int, std::string, + << "Name at position " << arg1 << " is not equal to " << arg2); +private: + /// The stored data + std::vector data; + /// The names of the stored data + std::vector names; }; @@ -156,9 +211,9 @@ AnyData::AnyData() template inline -AnyData::AnyData(const NamedData& other) +AnyData::AnyData(const NamedData &other) { - for (unsigned int i=0;i(&data[i]); + type *p = boost::any_cast(&data[i]); Assert(p != 0, - ExcTypeMismatch(typeid(type).name(),data[i].type().name())); + ExcTypeMismatch(typeid(type).name(),data[i].type().name())); return *p; } @@ -191,9 +246,11 @@ const type AnyData::entry (const unsigned int i) const { AssertIndexRange(i, size()); - const type* p = boost::any_cast(&data[i]); + const type *p = boost::any_cast(&data[i]); + if (p==0 ) + p = boost::any_cast(&data[i]); Assert(p != 0, - ExcTypeMismatch(typeid(type).name(),data[i].type().name())); + ExcTypeMismatch(typeid(type).name(),data[i].type().name())); return *p; } @@ -204,15 +261,58 @@ const type AnyData::read(const unsigned int i) const { AssertIndexRange(i, size()); - const type* p = boost::any_cast(&data[i]); + const type *p = boost::any_cast(&data[i]); + if (p==0) + p = boost::any_cast(&data[i]); Assert(p != 0, - ExcTypeMismatch(typeid(type).name(),data[i].type().name())); + ExcTypeMismatch(typeid(type).name(),data[i].type().name())); return *p; } +template inline -const std::string& +const type * +AnyData::read_ptr(const unsigned int i) const +{ + AssertIndexRange(i, size()); + const type *const *p = boost::any_cast(&data[i]); + if (p==0) + p = boost::any_cast(&data[i]); + Assert(p != 0, + ExcTypeMismatch(typeid(type *).name(),data[i].type().name())); + return *p; +} + + +template +inline +const type * +AnyData::try_read_ptr(const unsigned int i) const +{ + AssertIndexRange(i, size()); + const type *const *p = boost::any_cast(&data[i]); + if (p==0) + p = boost::any_cast(&data[i]); + return *p; +} + + +template +inline +const type * +AnyData::try_read(const unsigned int i) const +{ + AssertIndexRange(i, size()); + const type *p = boost::any_cast(&data[i]); + if (p==0) + p = boost::any_cast(&data[i]); + return p; +} + + +inline +const std::string & AnyData::name(const unsigned int i) const { AssertIndexRange(i, size()); @@ -222,17 +322,29 @@ AnyData::name(const unsigned int i) const inline unsigned int -AnyData::find(const std::string& n) const +AnyData::try_find(const std::string &n) const { std::vector::const_iterator it = std::find(names.begin(), names.end(), n); - - Assert(it != names.end(), ExcMessage("An entry with this name does not exist")); - + + if (it == names.end()) + return numbers::invalid_unsigned_int; + return it - names.begin(); } +inline +unsigned int +AnyData::find(const std::string &n) const +{ + const unsigned int i = try_find(n); + Assert(i != numbers::invalid_unsigned_int, ExcNameNotFound(n)); + + return i; +} + + template inline bool @@ -245,12 +357,12 @@ AnyData::is_type(const unsigned int i) const template inline type -AnyData::entry (const std::string& n) +AnyData::entry (const std::string &n) { const unsigned int i = find(n); - type* p = boost::any_cast(&data[i]); + type *p = boost::any_cast(&data[i]); Assert(p != 0, - ExcTypeMismatch(typeid(type).name(),data[i].type().name())); + ExcTypeMismatch(typeid(type).name(),data[i].type().name())); return *p; } @@ -258,12 +370,12 @@ AnyData::entry (const std::string& n) template inline const type -AnyData::entry (const std::string& n) const +AnyData::entry (const std::string &n) const { const unsigned int i = find(n); - const type* p = boost::any_cast(&data[i]); + const type *p = boost::any_cast(&data[i]); Assert(p != 0, - ExcTypeMismatch(typeid(type).name(),data[i].type().name())); + ExcTypeMismatch(typeid(type).name(),data[i].type().name())); return *p; } @@ -271,20 +383,51 @@ AnyData::entry (const std::string& n) const template inline const type -AnyData::read(const std::string& n) const +AnyData::read(const std::string &n) const +{ + const unsigned int i = find(n); + const type *p = boost::any_cast(&data[i]); + Assert(p != 0, + ExcTypeMismatch(typeid(type).name(),data[i].type().name())); + return *p; +} + + +template +inline +const type * +AnyData::read_ptr(const std::string &n) const { const unsigned int i = find(n); - const type* p = boost::any_cast(&data[i]); + const type *const *p = boost::any_cast(&data[i]); + if (p==0) + p = boost::any_cast(&data[i]); Assert(p != 0, - ExcTypeMismatch(typeid(type).name(),data[i].type().name())); + ExcTypeMismatch(typeid(type).name(),data[i].type().name())); + return *p; +} + + +template +inline +const type * +AnyData::try_read_ptr(const std::string &n) const +{ + const unsigned int i = try_find(n); + if (i == numbers::invalid_unsigned_int) + return 0; + + const type *const *p = boost::any_cast(&data[i]); + if (p==0) + p = boost::any_cast(&data[i]); return *p; } template inline -const type* -AnyData::try_read(const std::string& n) const +const type * +AnyData::try_read(const std::string &n) const { // Try to find name std::vector::const_iterator it = @@ -295,7 +438,7 @@ AnyData::try_read(const std::string& n) const // Compute index and return casted pointer unsigned int i=it-names.begin(); - const type* p = boost::any_cast(&data[i]); + const type *p = boost::any_cast(&data[i]); return p; } @@ -303,7 +446,7 @@ AnyData::try_read(const std::string& n) const template inline void -AnyData::add(type ent, const std::string& n) +AnyData::add(type ent, const std::string &n) { boost::any e = ent; data.push_back(e); @@ -313,7 +456,7 @@ AnyData::add(type ent, const std::string& n) inline void -AnyData::merge(const AnyData& other) +AnyData::merge(const AnyData &other) { for (unsigned int i=0; i &residual, Operator &inverse_derivative); /** - * Declare the parameters - * applicable to Newton's method. + * Declare the parameters applicable to Newton's method. */ - void declare_parameters (ParameterHandler ¶m); + static void declare_parameters (ParameterHandler ¶m); /** - * Read the parameters. + * Read the parameters in the ParameterHandler. */ - void initialize (ParameterHandler ¶m); + void parse_parameters (ParameterHandler ¶m); + + /** + * @deprecated Use parse_parameters(). + */ + void initialize (ParameterHandler ¶m) DEAL_II_DEPRECATED; /** * Initialize the pointer @@ -111,7 +115,7 @@ namespace Algorithms /** * @deprecated Use the function using AnyData */ - virtual void operator() (NamedData &out, const NamedData &in); + virtual void operator() (NamedData &out, const NamedData &in); virtual void notify(const Event &); diff --git a/deal.II/include/deal.II/algorithms/newton.templates.h b/deal.II/include/deal.II/algorithms/newton.templates.h index acb61f4712..9e3d62466d 100644 --- a/deal.II/include/deal.II/algorithms/newton.templates.h +++ b/deal.II/include/deal.II/algorithms/newton.templates.h @@ -55,7 +55,7 @@ namespace Algorithms template void - Newton::initialize (ParameterHandler ¶m) + Newton::parse_parameters (ParameterHandler ¶m) { param.enter_subsection("Newton"); control.parse_parameters (param); @@ -65,6 +65,13 @@ namespace Algorithms param.leave_subsection (); } + template + void + Newton::initialize (ParameterHandler ¶m) + { + parse_parameters(param); + } + template void Newton::initialize (OutputOperator &output) @@ -80,7 +87,7 @@ namespace Algorithms inverse_derivative->notify(e); } - + template double Newton::threshold(const double thr) @@ -89,8 +96,8 @@ namespace Algorithms assemble_threshold = thr; return t; } - - + + template void Newton::operator() (NamedData &out, const NamedData &in) @@ -105,7 +112,7 @@ namespace Algorithms Assert (out.size() == 1, ExcNotImplemented()); deallog.push ("Newton"); - VECTOR &u = *out.entry(0); + VECTOR &u = *out.entry(0); if (debug>2) deallog << "u: " << u.l2_norm() << std::endl; @@ -117,14 +124,14 @@ namespace Algorithms res->reinit(u); AnyData src1; AnyData src2; - src1.add(&u, "Newton iterate"); + src1.add(&u, "Newton iterate"); src1.merge(in); - src2.add(res, "Newton residual"); + src2.add(res, "Newton residual"); src2.merge(src1); AnyData out1; - out1.add(res, "Residual"); + out1.add(res, "Residual"); AnyData out2; - out2.add(Du, "Update"); + out2.add(Du, "Update"); unsigned int step = 0; // fill res with (f(u), v) diff --git a/deal.II/include/deal.II/algorithms/operator.h b/deal.II/include/deal.II/algorithms/operator.h index 6254d35019..2f5e935127 100644 --- a/deal.II/include/deal.II/algorithms/operator.h +++ b/deal.II/include/deal.II/algorithms/operator.h @@ -103,7 +103,7 @@ namespace Algorithms * of reassembling the matrix. */ Event notifications; - + }; /** @@ -141,7 +141,7 @@ namespace Algorithms { public: Operator(); - + /** * Implementation of the function in the base class in order to do * compatibility conversions between the old and the new @@ -157,6 +157,11 @@ namespace Algorithms */ virtual void operator() (NamedData &out, const NamedData &in); + /** + * Set this true to avoid compatibility warnings. + */ + bool silent_compatibility; + private: /** * While we are providing compatibility functions to the old @@ -198,12 +203,12 @@ namespace Algorithms /** * Output all the vectors in AnyData. */ - virtual OutputOperator &operator<< (const AnyData& vectors); - + virtual OutputOperator &operator<< (const AnyData &vectors); + /** * @deprecated Output all the vectors in NamedData. */ - virtual OutputOperator &operator<< (const NamedData &vectors); + OutputOperator &operator<< (const NamedData &vectors); protected: unsigned int step; private: diff --git a/deal.II/include/deal.II/algorithms/operator.templates.h b/deal.II/include/deal.II/algorithms/operator.templates.h index e90d771c54..0c1ee9d6cd 100644 --- a/deal.II/include/deal.II/algorithms/operator.templates.h +++ b/deal.II/include/deal.II/algorithms/operator.templates.h @@ -24,56 +24,56 @@ namespace Algorithms { template Operator::Operator() - : compatibility_flag(false) + : silent_compatibility(false), compatibility_flag(false) {} template void - Operator::operator() (AnyData& out, const AnyData& in) + Operator::operator() (AnyData &out, const AnyData &in) { // Had this function been overloaded in a derived clas, it would // not have been called. Therefore, we have to start the // compatibility engine. But before, we have to avoid an endless loop. Assert(!compatibility_flag, ExcMessage("Compatibility resolution of Operator generates and endless loop\n" - "Please provide an operator() in a derived class")); + "Please provide an operator() in a derived class")); compatibility_flag = true; - - NamedData new_out; - for (unsigned int i=0;i new_out; + for (unsigned int i=0; i(i)) - new_out.add(out.entry(i), out.name(i)); - else - deallog << "Cannot convert AnyData argument " << out.name(i) << " to NamedData" - << std::endl; + if (out.is_type(i)) + new_out.add(out.entry(i), out.name(i)); + else if (!silent_compatibility) + deallog << "Cannot convert AnyData argument " << out.name(i) << " to NamedData" + << std::endl; } - - NamedData new_in; - for (unsigned int i=0;i new_in; + for (unsigned int i=0; i(i)) - { - // This const cast is due to the wrong constness handling - // in NamedData. As soon as NamedData is gone, this code - // will not be necessary anymore. And deprecating begins - // now. - VECTOR* p = const_cast (in.entry(i)); - new_in.add(p, in.name(i)); - } - else if (in.is_type(i)) - { - // This const cast is due to the wrong constness handling - // in NamedData. As soon as NamedData is gone, this code - // will not be necessary anymore. And deprecating begins - // now. - VECTOR* p = const_cast (in.entry(i)); - new_in.add(p, in.name(i)); - } - else - deallog << "Cannot convert AnyData argument " << in.name(i) - << " to NamedData" << std::endl; + // deallog << "Convert " << in.name(i) << std::endl; + if (in.is_type(i)) + { + // This const cast is due to the wrong constness handling + // in NamedData. As soon as NamedData is gone, this code + // will not be necessary anymore. And deprecating begins + // now. + VECTOR *p = const_cast (in.entry(i)); + new_in.add(p, in.name(i)); + } + else if (in.is_type(i)) + { + // This const cast is due to the wrong constness handling + // in NamedData. As soon as NamedData is gone, this code + // will not be necessary anymore. And deprecating begins + // now. + VECTOR *p = const_cast (in.entry(i)); + new_in.add(p, in.name(i)); + } + else if (!silent_compatibility) + deallog << "Cannot convert AnyData argument " << in.name(i) + << " to NamedData" << std::endl; } this->operator() (new_out, new_in); compatibility_flag = false; @@ -82,21 +82,21 @@ namespace Algorithms template void - Operator::operator() (NamedData& out, const NamedData& in) + Operator::operator() (NamedData &out, const NamedData &in) { // Had this function been overloaded in a derived clas, it would // not have been called. Therefore, we have to start the // compatibility engine. But before, we have to avoid an endless loop. Assert(!compatibility_flag, ExcMessage("Compatibility resolution of Operator generates and endless loop\n" - "Please provide an operator() in a derived class")); + "Please provide an operator() in a derived class")); compatibility_flag = true; - + AnyData new_out; - for (unsigned int i=0;ioperator() (new_out, new_in); @@ -122,60 +122,45 @@ namespace Algorithms template OutputOperator & - OutputOperator::operator<< (const AnyData& vectors) + OutputOperator::operator<< (const AnyData &vectors) { if (os == 0) { - //TODO: make this possible - //deallog << ' ' << step; - //for (unsigned int i=0;iprint(deallog); - //deallog << std::endl; + deallog << "Step " << step << std::endl; + for (unsigned int i=0; i(i); + if (v == 0) continue; + deallog << vectors.name(i); + for (unsigned int j=0; jsize(); ++j) + deallog << ' ' << (*v)(j); + deallog << std::endl; + } + deallog << std::endl; } else { (*os) << ' ' << step; for (unsigned int i=0; i(i)) - { - const VECTOR& v = *vectors.entry(i); - for (unsigned int j=0; j(i)) - { - const VECTOR& v = *vectors.entry(i); - for (unsigned int j=0; j(i); + if (v == 0) continue; + for (unsigned int j=0; jsize(); ++j) + (*os) << ' ' << (*v)(j); + (*os) << std::endl; + } (*os) << std::endl; } return *this; } - + template OutputOperator & OutputOperator::operator<< (const NamedData &vectors) { - if (os == 0) - { - //TODO: make this possible - //deallog << ' ' << step; - //for (unsigned int i=0;iprint(deallog); - //deallog << std::endl; - } - else - { - (*os) << ' ' << step; - for (unsigned int i=0; isize(); ++j) - (*os) << ' ' << (*vectors(i))(j); - (*os) << std::endl; - } + const AnyData newdata = vectors; + (*this) << newdata; return *this; } } diff --git a/deal.II/include/deal.II/algorithms/theta_timestepping.h b/deal.II/include/deal.II/algorithms/theta_timestepping.h index e6bb52410a..ff6262cbab 100644 --- a/deal.II/include/deal.II/algorithms/theta_timestepping.h +++ b/deal.II/include/deal.II/algorithms/theta_timestepping.h @@ -1,7 +1,7 @@ // --------------------------------------------------------------------- // $Id$ // -// Copyright (C) 2010 - 2013 by the deal.II authors +// Copyright (C) 2010 - 2014 by the deal.II authors // // This file is part of the deal.II library. // @@ -114,7 +114,7 @@ namespace Algorithms * communicate the current time step information through AnyData as * well. Therefore, the AnyData objects handed as input to * #op_explicit and #op_implicit contain two entries of type - * const double reference named "Time" and "Timestep". Note + * const double* named "Time" and "Timestep". Note * that "Time" refers to the time at the beginning of the current * step for #op_explicit and at the end for #op_implicit, * respectively. @@ -140,7 +140,7 @@ namespace Algorithms * @until End of declarations * * These operators will be implemented after the main program. But let - * us look at how they get used. First, let us define a matrix to be + * us look first at how they get used. First, let us define a matrix to be * used for our system and also an OutputOperator in order to write * the data of each timestep to a file. * @@ -148,10 +148,9 @@ namespace Algorithms * @until out.initialize * * Now we create objects for the implicit and explicit parts of the - * steps as well as the ThetaTimestepping itself. Notice how the - * TimestepData of ThetaTimestepping gets forwarded to the inner - * operators. There are two different data objects, because the - * timestep size is modified by #theta. + * steps as well as the ThetaTimestepping itself. We initialize the + * timestepping with the output operator in order to be able to see + * the output in every step. * * @until set_output * @@ -159,29 +158,47 @@ namespace Algorithms * is filled with the initial value and is also the vector where the * solution at each timestep will be. Because the interface of * Operator has to be able to handle several vectors, we need to store - * it in a NamedData object. Notice, that we need to create the - * intermediate pointer p. If we would use - * &value directly in the add function, the - * resulting object would be constant. + * it in an AnyData object. Since our problem has no additional + * parameters, the input AnyData object remains empty. * * @until add * - * Finally, we are ready to tell the solver, that we are looknig at + * Finally, we are ready to tell the solver, that we are starting at * the initial timestep and run it. * - * @until outdata - * @skip Explicit::initialize + * @until } + * + * First the constructor, which simply copies the system matrix into + * the member pointer for later use. + * + * @skip Explicit:: + * @until } * * Now we need to study the application of the implicit and explicit * operator. We assume that the pointer matrix points to - * the matrix created in the main program, and that - * timestep_data points to the correct data object of - * ThetaTimestepping. + * the matrix created in the main program (the constructor did this + * for us). Here, we first get the time step size from the AnyData + * object that was provided as input. Then, if we are in the first + * step or if the timestep has changed, we fill the local matrix + * m, such that with the given matrix \f$M\f$, it becomes + * \f[ + * m = I - \Delta t M. + * \f] + * After we have worked off the notifications, we clear them, such + * that the matrix is only generated when necessary. * * @skipline void - * @until vmult + * @until clear + * + * Now we multiply the input vector with the new matrix and store on output. + * * @until } + * The code for the implicit operator is almost the same, except + * that we change the sign in front of the timestep and use the + * inverse of t he matrix. * + * @until vmult + * @until } * @author Guido Kanschat * @date 2010 */ @@ -223,21 +240,38 @@ namespace Algorithms */ void set_output(OutputOperator &output); - void declare_parameters (ParameterHandler ¶m); - void initialize (ParameterHandler ¶m); + /** + * Declare parameters in a parameter handler. + */ + static void declare_parameters (ParameterHandler ¶m); + + /** + * Read the parameters in the ParameterHandler. + */ + void parse_parameters (ParameterHandler ¶m); + + /** + * @deprecated Use parse_parameters(). + */ + void initialize (ParameterHandler ¶m) DEAL_II_DEPRECATED; /** * The current time in the * timestepping scheme. */ - const double ¤t_time() const; + double current_time() const; /** * The current step size. */ - const double &step_size() const; + double step_size() const; /** * The weight between implicit and explicit part. */ - const double &theta() const; + double theta() const; + + /** + * Set a new weight and return the old + */ + double theta(double new_theta); /** * The data handed to the #op_explicit time stepping operator. @@ -354,6 +388,32 @@ namespace Algorithms { output = &out; } + + + template + inline + double ThetaTimestepping::theta () const + { + return vtheta; + } + + + template + inline + double ThetaTimestepping::theta (double new_theta) + { + const double tmp = vtheta; + vtheta = new_theta; + return tmp; + } + + + template + inline + double ThetaTimestepping::current_time () const + { + return control.now(); + } } DEAL_II_NAMESPACE_CLOSE diff --git a/deal.II/include/deal.II/algorithms/theta_timestepping.templates.h b/deal.II/include/deal.II/algorithms/theta_timestepping.templates.h index 348ab379cb..a3402f6dd0 100644 --- a/deal.II/include/deal.II/algorithms/theta_timestepping.templates.h +++ b/deal.II/include/deal.II/algorithms/theta_timestepping.templates.h @@ -51,7 +51,7 @@ namespace Algorithms template void - ThetaTimestepping::initialize (ParameterHandler ¶m) + ThetaTimestepping::parse_parameters (ParameterHandler ¶m) { param.enter_subsection("ThetaTimestepping"); control.parse_parameters (param); @@ -61,14 +61,22 @@ namespace Algorithms } + template + void + ThetaTimestepping::initialize (ParameterHandler ¶m) + { + parse_parameters(param); + } + + template void ThetaTimestepping::operator() (NamedData &out, const NamedData &in) { Operator::operator() (out, in); } - - + + template void ThetaTimestepping::operator() (AnyData &out, const AnyData &in) @@ -77,7 +85,7 @@ namespace Algorithms deallog.push ("Theta"); - VECTOR& solution = *out.entry(0); + VECTOR &solution = *out.entry(0); GrowingVectorMemory mem; typename VectorMemory::Pointer aux(mem); aux->reinit(solution); @@ -90,25 +98,33 @@ namespace Algorithms // vector associated with the old // timestep AnyData src1; - src1.add(&solution, "Previous iterate"); - src1.add(d_explicit.time, "Time"); - src1.add(d_explicit.step, "Timestep"); + src1.add(&solution, "Previous iterate"); + src1.add(&d_explicit.time, "Time"); + src1.add(&d_explicit.step, "Timestep"); + src1.add(&vtheta, "Theta"); src1.merge(in); AnyData src2; - src2.add(&solution, "Previous iterate"); AnyData out1; - out1.add(aux, "Solution"); + out1.add(aux, "Solution"); // The data provided to the inner solver - src2.add(aux, "Previous time"); - src2.add(d_implicit.time, "Time"); - src2.add(d_explicit.step, "Timestep"); + src2.add(aux, "Previous time"); + src2.add(&solution, "Previous iterate"); + src2.add(&d_implicit.time, "Time"); + src2.add(&d_implicit.step, "Timestep"); + src2.add(&vtheta, "Theta"); src2.merge(in); if (output != 0) (*output) << 0U << out; + // Avoid warnings because time and timestep cannot be converted to VECTOR* + const bool explicit_silent = op_explicit->silent_compatibility; + const bool implicit_silent = op_implicit->silent_compatibility; + op_explicit->silent_compatibility = true; + op_implicit->silent_compatibility = true; + for (unsigned int count = 1; d_explicit.time < control.final(); ++count) { const bool step_change = control.advance(); @@ -135,6 +151,8 @@ namespace Algorithms d_explicit.time = control.now(); } + op_explicit->silent_compatibility = explicit_silent; + op_implicit->silent_compatibility = implicit_silent; deallog.pop(); } } diff --git a/deal.II/include/deal.II/fe/fe_tools.h b/deal.II/include/deal.II/fe/fe_tools.h index 5b9e220f74..854f26041f 100644 --- a/deal.II/include/deal.II/fe/fe_tools.h +++ b/deal.II/include/deal.II/fe/fe_tools.h @@ -1,7 +1,7 @@ // --------------------------------------------------------------------- // $Id$ // -// Copyright (C) 2000 - 2013 by the deal.II authors +// Copyright (C) 2000 - 2014 by the deal.II authors // // This file is part of the deal.II library. // @@ -697,21 +697,64 @@ namespace FETools /** * Gives the patchwise extrapolation of a @p dof1 function @p z1 to a @p - * dof2 function @p z2. @p dof1 and @p dof2 need to be DoFHandler based on - * the same triangulation. - * - * This function is interesting for e.g. extrapolating patchwise a piecewise - * linear solution to a piecewise quadratic solution. - * - * Note that the resulting field does not satisfy continuity requirements of - * the given finite elements. - * + * dof2 function @p z2. @p dof1 and @p dof2 need to be DoFHandler objects + * based on the same triangulation. This function is used, for example, for + * extrapolating patchwise a piecewise linear solution to a piecewise + * quadratic solution. + * + * The function's name is historical and probably not particularly well + * chosen. The function performs the following operations, one after the + * other: + * + * - It interpolates directly from every cell of @p dof1 to the + * corresponding cell of @dof2 using the interpolation matrix of the + * finite element spaces used on these cells and provided by + * the finite element objects involved. This step is done using the + * FETools::interpolate() function. + * - It then performs a loop over all non-active cells of @dof2. If + * such a non-active cell has at least one active child, then we + * call the children of this cell a "patch". We then interpolate + * from the children of this patch to the patch, using the finite + * element space associated with @p dof2 and immediately interpolate + * back to the children. In essence, this information throws away + * all information in the solution vector that lives on a scale + * smaller than the patch cell. + * - Since we traverse non-active cells from the coarsest to the + * finest levels, we may find patches that correspond to child + * cells of previously treated patches if the mesh had been + * refined adaptively (this cannot happen if the mesh has been + * refined globally because there the children of a patch are + * all active). We also perform the operation described above + * on these patches, but it is easy to see that on patches that + * are children of previously treated patches, the operation is + * now the identity operation (since it interpolates from the + * children of the current patch a function that had previously + * been interpolated to these children from an even coarser patch). + * Consequently, this does not alter the solution vector any more. + * + * The name of the function originates from the fact that it can be + * used to construct a representation of a function of higher polynomial + * degree on a once coarser mesh. For example, if you imagine that you + * start with a $Q_1$ function on globally refined mesh, and that @p dof2 + * is associated with a $Q_2$ element, then this function computes the + * equivalent of the operator $I_{2h}^{(2)}$ interpolating the original + * piecewise linear function onto a quadratic function on a once coarser + * mesh with mesh size $2h$ (but representing this function on the original + * mesh with size $h$). If the exact solution is sufficiently smooth, + * then $u^\ast=I_{2h}^{(2)}u_h$ is typically a better approximation to + * the exact solution $u$ of the PDE than $u_h$ is. In other words, this + * function provides a postprocessing step that improves the solution in + * a similar way one often obtains by extrapolating a sequence of solutions, + * explaining the origin of the function's name. + * + * @note The resulting field does not satisfy continuity requirements of + * the given finite elements if the algorithm outlined above is used. * When you use continuous elements on grids with hanging nodes, please use * the @p extrapolate function with an additional ConstraintMatrix argument, * see below. * - * Since this function operates on patches of cells, it is required that the - * underlying grid is refined at least once for every coarse grid cell. If + * @note Since this function operates on patches of cells, it requires that + * the underlying grid is refined at least once for every coarse grid cell. If * this is not the case, an exception will be raised. */ template @@ -722,13 +765,15 @@ namespace FETools /** * Gives the patchwise extrapolation of a @p dof1 function @p z1 to a @p - * dof2 function @p z2. @p dof1 and @p dof2 need to be DoFHandler based on + * dof2 function @p z2. @p dof1 and @p dof2 need to be DoFHandler objects based on * the same triangulation. @p constraints is a hanging node constraints - * object corresponding to @p dof2. This object is particular important when + * object corresponding to @p dof2. This object is necessary when * interpolating onto continuous elements on grids with hanging nodes * (locally refined grids). * - * Otherwise, the same holds as for the other @p extrapolate function. + * Otherwise, the function does the same as the other @p extrapolate + * function above (for which the documentation provides an extensive + * description of its operation). */ template void extrapolate (const DoFHandler &dof1, @@ -803,17 +848,15 @@ namespace FETools * accordingly. * * The name must be in the form which is returned by the - * FiniteElement::get_name function, where a few modifications are allowed: - * - *
        - *
      • Dimension template parameters <2> etc. can be + * FiniteElement::get_name function, where dimension template parameters <2> etc. can be * omitted. Alternatively, the explicit number can be replaced by * dim or d. If a number is given, it must match * the template parameter of this function. * - *
      • The powers used for FESystem may either be numbers or can be + * The names of FESystem elements follow the pattern FESystem[FE_Base1^p1-FE_Base2^p2] + * The powers p1 etc. may either be numbers or can be * replaced by dim or d. - *
      + * * * If no finite element can be reconstructed from this string, an exception * of type @p FETools::ExcInvalidFEName is thrown. diff --git a/deal.II/include/deal.II/grid/tria_accessor.h b/deal.II/include/deal.II/grid/tria_accessor.h index 6d693744fa..13b3a76d56 100644 --- a/deal.II/include/deal.II/grid/tria_accessor.h +++ b/deal.II/include/deal.II/grid/tria_accessor.h @@ -3136,12 +3136,11 @@ public: * may be different as to whether * a point is in- or outside the * cell in real space. - - In case of codim>0, the point is first projected - to the manifold where the cell is embedded and - then check if this projection is inside the cell. - - */ + * + * In case of codim>0, the point is first projected + * to the manifold where the cell is embedded and + * then check if this projection is inside the cell. + */ bool point_inside (const Point &p) const; /** diff --git a/deal.II/include/deal.II/lac/petsc_solver.h b/deal.II/include/deal.II/lac/petsc_solver.h index b3c9ac5c84..b87797cd4f 100644 --- a/deal.II/include/deal.II/lac/petsc_solver.h +++ b/deal.II/include/deal.II/lac/petsc_solver.h @@ -1190,7 +1190,7 @@ namespace PETScWrappers /** * The method allows to take advantage * if the system matrix is symmetric by - * using LDL^T decomposition unstead of + * using LDL^T decomposition instead of * more expensive LU. The argument * indicates whether the matrix is * symmetric or not. diff --git a/deal.II/include/deal.II/lac/trilinos_sparse_matrix.h b/deal.II/include/deal.II/lac/trilinos_sparse_matrix.h index c5ab727b66..21b2f10545 100644 --- a/deal.II/include/deal.II/lac/trilinos_sparse_matrix.h +++ b/deal.II/include/deal.II/lac/trilinos_sparse_matrix.h @@ -1325,7 +1325,10 @@ namespace TrilinosWrappers const TrilinosScalar new_diag_value = 0); /** - * Make an in-place transpose of a matrix. + * Sets an internal flag so that all operations performed by the matrix, + * i.e., multiplications, are done in transposed order. However, this does + * not reshape the matrix to transposed form directly, so care should be + * taken when using this flag. */ void transpose (); diff --git a/deal.II/include/deal.II/meshworker/functional.h b/deal.II/include/deal.II/meshworker/functional.h index a651ba1acc..0205951f21 100644 --- a/deal.II/include/deal.II/meshworker/functional.h +++ b/deal.II/include/deal.II/meshworker/functional.h @@ -19,6 +19,7 @@ #define __deal2__mesh_worker_functional_h #include +#include #include #include #include @@ -107,49 +108,28 @@ namespace MeshWorker { public: /** - * The data type for - * communicating the cell and - * face vectors. - */ - typedef NamedData*> DataVectors; - - /** - * The initialization - * function, specifying the - * @p results vectors and - * whether face data should - * be collected separately. + * The initialization function, specifying the @p results + * vectors and whether face data should be collected separately. * - * @p results should contain - * two block vectors named - * "cells" and "faces" (the - * latter only if - * @p separate_faces is - * true). In each of the two, - * each block should have - * equal size and be large - * enough to accommodate all - * user indices set in the - * cells and faces covered by - * the loop it is used - * in. Typically, for - * estimators, this is - * Triangulation::n_active_cells() - * and - * Triangulation::n_faces(), - * respectively. + * @p results should contain two block vectors named "cells" and + * "faces" (the latter only if @p separate_faces is true). In + * each of the two, each block should have equal size and be + * large enough to accommodate all user indices set in the cells + * and faces covered by the loop it is used in. Typically, for + * estimators, this is Triangulation::n_active_cells() and + * Triangulation::n_faces(), respectively. * - * The use of BlockVector may - * seem cumbersome, but it - * allows us to assemble - * several functionals at the - * same time, one in each - * block. The typical - * situation for error - * estimate is just having a - * single block in each vector. + * The use of BlockVector may seem cumbersome, but it allows us + * to assemble several functionals at the same time, one in each + * block. The typical situation for error estimate is just + * having a single block in each vector. */ - void initialize(DataVectors &results, + void initialize(AnyData &results, bool separate_faces = true); + + /** + * @deprecated + */ + void initialize(NamedData*> &results, bool separate_faces = true); /** * Initialize the local data @@ -188,7 +168,7 @@ namespace MeshWorker */ number operator() (const unsigned int i) const; private: - DataVectors results; + AnyData results; bool separate_faces; }; //----------------------------------------------------------------------// @@ -247,13 +227,29 @@ namespace MeshWorker template inline void - CellsAndFaces::initialize(DataVectors &r, bool sep) + CellsAndFaces::initialize(AnyData &r, bool sep) + { + Assert(r.name(0) == "cells", AnyData::ExcNameMismatch(0, "cells")); + if (sep) + { + Assert(r.name(1) == "faces", AnyData::ExcNameMismatch(1, "faces")); + AssertDimension(r.entry*>(0)->n_blocks(), + r.entry*>(1)->n_blocks()); + } + + results = r; + separate_faces = sep; + } + + template + inline void + CellsAndFaces::initialize(NamedData*> &r, bool sep) { - Assert(r.name(0) == "cells", typename DataVectors::ExcNameMismatch(0, "cells")); + Assert(r.name(0) == "cells", AnyData::ExcNameMismatch(0, "cells")); if (sep) { - Assert(r.name(1) == "faces", typename DataVectors::ExcNameMismatch(1, "faces")); - AssertDimension(r(0)->n_blocks(), r(1)->n_blocks()); + Assert(r.name(1) == "faces", AnyData::ExcNameMismatch(1, "faces")); + AssertDimension(r(0)->n_blocks(),r(1)->n_blocks()); } results = r; @@ -266,7 +262,7 @@ namespace MeshWorker inline void CellsAndFaces::initialize_info(DOFINFO &info, bool) const { - info.initialize_numbers(results(0)->n_blocks()); + info.initialize_numbers(results.entry*>(0)->n_blocks()); } @@ -275,14 +271,15 @@ namespace MeshWorker inline void CellsAndFaces::assemble(const DOFINFO &info) { + BlockVector* v; + if (separate_faces && + info.face_number != deal_II_numbers::invalid_unsigned_int) + v = results.entry*>(1); + else + v = results.entry*>(0); + for (unsigned int i=0; iblock(i)(info.face->user_index()) += info.value(i); - else - results(0)->block(i)(info.cell->user_index()) += info.value(i); - } + v->block(i)(info.cell->user_index()) += info.value(i); } @@ -296,15 +293,17 @@ namespace MeshWorker { if (separate_faces) { + BlockVector* v1 = results.entry*>(1); const double J = info1.value(i) + info2.value(i); - results(1)->block(i)(info1.face->user_index()) += J; + v1->block(i)(info1.face->user_index()) += J; if (info2.face != info1.face) - results(1)->block(i)(info2.face->user_index()) += J; + v1->block(i)(info2.face->user_index()) += J; } else { - results(0)->block(i)(info1.cell->user_index()) += .5*info1.value(i); - results(0)->block(i)(info2.cell->user_index()) += .5*info2.value(i); + BlockVector* v0 = results.entry*>(0); + v0->block(i)(info1.cell->user_index()) += .5*info1.value(i); + v0->block(i)(info2.cell->user_index()) += .5*info2.value(i); } } } diff --git a/deal.II/include/deal.II/meshworker/integration_info.h b/deal.II/include/deal.II/meshworker/integration_info.h index 4eb299711e..a120eced0f 100644 --- a/deal.II/include/deal.II/meshworker/integration_info.h +++ b/deal.II/include/deal.II/meshworker/integration_info.h @@ -248,7 +248,7 @@ namespace MeshWorker * undertaken to use this class. * * First, you should consider if you need values from any vectors in a - * NamedData object. If so, fill the VectorSelector objects + * AnyData object. If so, fill the VectorSelector objects * #cell_selector, #boundary_selector and #face_selector with their names * and the data type (value, gradient, Hessian) to be extracted. * @@ -313,17 +313,22 @@ namespace MeshWorker * formulas, which integrate polynomial bilinear forms exactly. */ template + void initialize(const FiniteElement &el, + const Mapping &mapping, + const AnyData &data, + const VECTOR& dummy, + const BlockInfo *block_info = 0); + /** + * @deprecated Use AnyData instead of NamedData. + */ + template void initialize(const FiniteElement &el, const Mapping &mapping, const NamedData &data, const BlockInfo *block_info = 0); /** - * Initialize the IntegrationInfo objects contained. - * - * Before doing so, add update flags necessary to produce the data - * needed and also set uninitialized quadrature rules to Gauss - * formulas, which integrate polynomial bilinear forms exactly. + * @deprecated Use AnyData instead of NamedData. */ template void initialize(const FiniteElement &el, @@ -767,6 +772,44 @@ namespace MeshWorker } + template + template + void + IntegrationInfoBox::initialize( + const FiniteElement &el, + const Mapping &mapping, + const AnyData &data, + const VECTOR& dummy, + const BlockInfo *block_info) + { + initialize(el, mapping, block_info); + std_cxx1x::shared_ptr > p; + VectorDataBase* pp; + + p = std_cxx1x::shared_ptr >(new VectorData (cell_selector)); + // Public member function of parent class was not found without + // explicit cast + pp = &*p; + pp->initialize(data); + cell_data = p; + cell.initialize_data(p); + + p = std_cxx1x::shared_ptr >(new VectorData (boundary_selector)); + pp = &*p; + pp->initialize(data); + boundary_data = p; + boundary.initialize_data(p); + + p = std_cxx1x::shared_ptr >(new VectorData (face_selector)); + pp = &*p; + pp->initialize(data); + face_data = p; + face.initialize_data(p); + subface.initialize_data(p); + neighbor.initialize_data(p); + } + + template template void diff --git a/deal.II/include/deal.II/meshworker/simple.h b/deal.II/include/deal.II/meshworker/simple.h index ed692c57f3..297e2421a6 100644 --- a/deal.II/include/deal.II/meshworker/simple.h +++ b/deal.II/include/deal.II/meshworker/simple.h @@ -19,6 +19,7 @@ #define __deal2__mesh_worker_simple_h #include +#include #include #include #include @@ -60,8 +61,18 @@ namespace MeshWorker class ResidualSimple { public: - void initialize(NamedData &results); + /** + * Initialize with an AnyData object holding the result of + * assembling. + * + * Assembling currently writes into the first vector of results. + */ + void initialize(AnyData& results); + /** + * @deprecated Use initialize(AnyData&) instead. + */ + void initialize(NamedData &results); /** * Initialize the constraints. */ @@ -114,7 +125,7 @@ namespace MeshWorker /** * The global residal vectors filled by assemble(). */ - NamedData > > residuals; + AnyData residuals; /** * A pointer to the object containing constraints. */ @@ -500,6 +511,13 @@ namespace MeshWorker //----------------------------------------------------------------------// + template + inline void + ResidualSimple::initialize(AnyData &results) + { + residuals = results; + } + template inline void ResidualSimple::initialize(NamedData &results) @@ -537,18 +555,19 @@ namespace MeshWorker { for (unsigned int k=0; k(k); if (constraints == 0) { for (unsigned int i=0; idistribute_local_to_global(info.vector(k).block(0), info.indices, (*residuals(k))); + constraints->distribute_local_to_global(info.vector(k).block(0), info.indices, *v); else for (unsigned int i=0; i != info.vector(k).n_blocks(); ++i) - constraints->distribute_local_to_global(info.vector(k).block(i), info.indices_by_block[i], (*residuals(k))); + constraints->distribute_local_to_global(info.vector(k).block(i), info.indices_by_block[i], *v); } } } @@ -561,30 +580,31 @@ namespace MeshWorker { for (unsigned int k=0; k(k); if (constraints == 0) { for (unsigned int i=0; idistribute_local_to_global - (info1.vector(k).block(0), info1.indices, (*residuals(k))); + (info1.vector(k).block(0), info1.indices, *v); constraints->distribute_local_to_global - (info2.vector(k).block(0), info2.indices, (*residuals(k))); + (info2.vector(k).block(0), info2.indices, *v); } else if (info1.indices_by_block.size() != 0 && info2.indices_by_block.size() != 0) { for (unsigned int i=0; idistribute_local_to_global - (info1.vector(k).block(i), info1.indices_by_block[i], (*residuals(k))); + (info1.vector(k).block(i), info1.indices_by_block[i], *v); constraints->distribute_local_to_global - (info2.vector(k).block(i), info2.indices_by_block[i], (*residuals(k))); + (info2.vector(k).block(i), info2.indices_by_block[i], *v); } } else @@ -1132,7 +1152,7 @@ namespace MeshWorker inline void SystemSimple::initialize(MATRIX &m, VECTOR &rhs) { - NamedData data; + AnyData data; VECTOR *p = &rhs; data.add(p, "right hand side"); diff --git a/deal.II/include/deal.II/meshworker/vector_selector.h b/deal.II/include/deal.II/meshworker/vector_selector.h index 5976df4dbd..e35975a1bb 100644 --- a/deal.II/include/deal.II/meshworker/vector_selector.h +++ b/deal.II/include/deal.II/meshworker/vector_selector.h @@ -18,6 +18,7 @@ #define __deal2__mesh_worker_vector_selector_h #include +#include #include #include #include @@ -32,7 +33,7 @@ namespace MeshWorker /** * A class that selects vectors from a list of named vectors. * - * Since the number of vectors in a NamedData object may grow with every + * Since the number of vectors in an AnyData object may grow with every * nesting of applications or loops, it is important to be able to * select those, which are actually used in computing residuals etc. * This class organizes the selection. @@ -49,23 +50,13 @@ namespace MeshWorker { public: /** - * Add a vector to the - * selection of finite element - * functions. The arguments are - * the name of the vector and - * indicators, which - * information is to be - * extracted from the - * vector. The name refers to - * an entry in a NamedData - * object, which will be - * identified by initialize(). - * The three bool parameters - * indicate, whether values, - * greadients and Hessians of - * the finite element function - * are to be computed on each - * cell or face. + * Add a vector to the selection of finite element functions. The + * arguments are the name of the vector and indicators, which + * information is to be extracted from the vector. The name refers + * to an entry in a AnyData object, which will be identified by + * initialize(). The three bool parameters indicate, whether + * values, greadients and Hessians of the finite element function + * are to be computed on each cell or face. */ void add(const std::string &name, const bool values = true, @@ -73,10 +64,8 @@ namespace MeshWorker const bool hessians = false); /** - * Does the same as the function - * above but it is possible to - * select a block of the global - * vector. + * Does the same as the function above but it is possible to + * select a block of the global vector. */ // void add(const std::string& name, // const unsigned int selected_block, @@ -85,21 +74,19 @@ namespace MeshWorker // bool hessians = false); /** - * Initialize the selection - * field with a data - * vector. While add() only - * enters the names of vectors, - * which will be used in the integration - * loop over cells and faces, - * this function links the - * names to actual vectos in a - * NamedData object. + * Initialize the selection field with a data vector. While add() + * only enters the names of vectors, which will be used in the + * integration loop over cells and faces, this function links the + * names to actual vectos in a AnyData object. * - * @note This function caches - * the index associated with a - * name. Therefore, it must be - * called every time after the - * NamedData object has changed. + * @note This function caches the index associated with a + * name. Therefore, it must be called every time after the AnyData + * object has changed. + */ + void initialize(const AnyData &); + + /** + * @deprecated Use AnyData instead of NamedData. */ template void initialize(const NamedData &); @@ -110,20 +97,17 @@ namespace MeshWorker bool empty () const; /** - * Returns true if values are - * selected for any vector. + * Returns true if values are selected for any vector. */ bool has_values () const; /** - * Returns true if gradients are - * selected for any vector. + * Returns true if gradients are selected for any vector. */ bool has_gradients () const; /** - * Returns true if hessians are - * selected for any vector. + * Returns true if hessians are selected for any vector. */ bool has_hessians () const; @@ -162,6 +146,12 @@ namespace MeshWorker * selection to the stream. */ template + void print (STREAM &s, const AnyData &v) const; + + /** + * @deprecated Use AnyData instead of NamedData. + */ + template void print (STREAM &s, const NamedData &v) const; /** @@ -183,14 +173,12 @@ namespace MeshWorker NamedSelection value_selection; /** - * Selection of the vectors - * used to compute gradients. + * Selection of the vectors used to compute gradients. */ NamedSelection gradient_selection; /** - * Selection of the vectors - * used to compute hessians. + * Selection of the vectors used to compute hessians. */ NamedSelection hessian_selection; }; @@ -213,81 +201,61 @@ namespace MeshWorker VectorDataBase(); /** - * Constructor from a base - * class object + * Constructor from a base class object */ VectorDataBase(const VectorSelector &); /** - * Virtual, but empty - * destructor. + * Initialize with a AnyData object and cache the indices in the + * VectorSelector base class. + * + * @note Make sure the VectorSelector base class was filled with + * reasonable data before calling this function. + */ + void initialize(const AnyData&); + + /** + * Virtual, but empty destructor. */ virtual ~VectorDataBase(); /** - * The only function added to - * VectorSelector is an - * abstract virtual function - * implemented in the derived - * class template and called by - * IntegrationInfo. + * The only function added to VectorSelector is an abstract + * virtual function implemented in the derived class template and + * called by IntegrationInfo. * - * Depending on the selections - * made in our base class, this - * fills the first three - * arguments with the local - * data of the finite element - * functions. It is usually - * called either for the whole - * FESystem, or for each base - * element separately. + * Depending on the selections made in our base class, this fills + * the first three arguments with the local data of the finite + * element functions. It is usually called either for the whole + * FESystem, or for each base element separately. * - * @param values is the vector - * filled with the values of - * the finite element function - * in the quadrature points. + * @param values is the vector filled with the values of the + * finite element function in the quadrature points. * - * @param gradients is the vector - * filled with the derivatives of - * the finite element function - * in the quadrature points. + * @param gradients is the vector filled with the derivatives of + * the finite element function in the quadrature points. * - * @param hessians is the - * vector filled with the - * second derivatives of the - * finite element function in - * the quadrature points. + * @param hessians is the vector filled with the second + * derivatives of the finite element function in the quadrature + * points. * - * @param fe is the - * FEValuesBase object which is - * used to compute the function - * values. Its UpdateFlags must - * have been set appropriately. + * @param fe is the FEValuesBase object which is used to compute + * the function values. Its UpdateFlags must have been set + * appropriately. * - * @param index is the local - * index vector. If @p fe - * refers to base elements of - * the system, this vector - * should be sorted by block - * and the arguments @p start - * and @p size below specify - * the subset of @p indices - * used. + * @param index is the local index vector. If @p fe refers to base + * elements of the system, this vector should be sorted by block + * and the arguments @p start and @p size below specify the subset + * of @p indices used. * - * @param component is the first - * index in @p values, @p - * gradients and @p hessians - * entered in this function. + * @param component is the first index in @p values, @p gradients + * and @p hessians entered in this function. * - * @param n_comp is the number - * of components to be filled. + * @param n_comp is the number of components to be filled. * - * @param start is the first - * index of this block in @p - * indices, or zero if - * no base elements are used. + * @param start is the first index of this block in @p indices, or + * zero if no base elements are used. * - * @param size is the number of - * dofs per cell of the current + * @param size is the number of dofs per cell of the current * element or base element. */ virtual void fill( @@ -302,15 +270,10 @@ namespace MeshWorker const unsigned int size) const; /** - * Fill the local data vector - * from level vectors. Performs - * exactly what the other - * fill() does, but uses the - * cell level to access a - * single level out of a - * hierarchy of level vectors, - * instead of a global data - * vector on the active cells. + * Fill the local data vector from level vectors. Performs exactly + * what the other fill() does, but uses the cell level to access a + * single level out of a hierarchy of level vectors, instead of a + * global data vector on the active cells. */ virtual void mg_fill( std::vector > > &values, @@ -328,13 +291,15 @@ namespace MeshWorker * The memory used by this object. */ std::size_t memory_consumption () const; + protected: + AnyData data; }; /** * Based on VectorSelector, this is the class that implements the * function VectorDataBase::fill() for a certain type of vector, using - * NamedData to identify vectors by name. + * AnyData to identify vectors by name. * * @ingroup MeshWorker * @author Guido Kanschat, 2009 @@ -355,30 +320,16 @@ namespace MeshWorker VectorData(const VectorSelector &); /** - * Initialize with a NamedData - * object and cache the indices - * in the VectorSelector base - * class. - * - * @note Make sure the - * VectorSelector base class - * was filled with reasonable - * data before calling this - * function. + * @deprecated Use AnyData instead of NamedData. */ void initialize(const NamedData &); /** - * Initialize with a single - * vector and cache the indices - * in the VectorSelector base - * class. + * Initialize with a single vector and cache the indices in the + * VectorSelector base class. * - * @note Make sure the - * VectorSelector base class - * was filled with reasonable - * data before calling this - * function. + * @note Make sure the VectorSelector base class was filled with + * reasonable data before calling this function. */ void initialize(const VECTOR *, const std::string &name); @@ -393,27 +344,36 @@ namespace MeshWorker const unsigned int start, const unsigned int size) const; - /** + virtual void mg_fill( + std::vector > > &values, + std::vector > > > &gradients, + std::vector > > > &hessians, + const FEValuesBase &fe, + const unsigned int level, + const std::vector &index, + const unsigned int component, + const unsigned int n_comp, + const unsigned int start, + const unsigned int size) const; + + /** * The memory used by this object. */ std::size_t memory_consumption () const; - - private: - NamedData > > data; }; /** * Based on VectorSelector, this is the class that implements the * function VectorDataBase::fill() for a certain type of multilevel vectors, using - * NamedData to identify vectors by name. + * AnyData to identify vectors by name. * * @ingroup MeshWorker * @author Guido Kanschat, 2010 */ template class MGVectorData : - public VectorDataBase + public VectorData { public: /** @@ -421,58 +381,29 @@ namespace MeshWorker */ MGVectorData(); /** - * Constructor using a - * prefilled VectorSelector + * Constructor using a prefilled VectorSelector */ MGVectorData(const VectorSelector &); /** - * Initialize with a NamedData - * object and cache the indices - * in the VectorSelector base - * class. - * - * @note Make sure the - * VectorSelector base class - * was filled with reasonable - * data before calling this - * function. + * @deprecated Use AnyData instead of NamedData. */ void initialize(const NamedData*> &); /** - * Initialize with a single - * vector and cache the indices - * in the VectorSelector base - * class. + * Initialize with a single vector and cache the indices in the + * VectorSelector base class. * - * @note Make sure the - * VectorSelector base class - * was filled with reasonable - * data before calling this - * function. + * @note Make sure the VectorSelector base class was filled with + * reasonable data before calling this function. */ void initialize(const MGLevelObject *, const std::string &name); - virtual void mg_fill( - std::vector > > &values, - std::vector > > > &gradients, - std::vector > > > &hessians, - const FEValuesBase &fe, - const unsigned int level, - const std::vector &index, - const unsigned int component, - const unsigned int n_comp, - const unsigned int start, - const unsigned int size) const; /** * The memory used by this object. */ std::size_t memory_consumption () const; - - private: - NamedData,MGVectorData > > data; }; @@ -504,6 +435,14 @@ namespace MeshWorker //} + inline void + VectorSelector::initialize(const AnyData &src) + { + value_selection.initialize(src); + gradient_selection.initialize(src); + hessian_selection.initialize(src); + } + template inline void VectorSelector::initialize(const NamedData &src) @@ -596,6 +535,23 @@ namespace MeshWorker } + template + inline void + VectorSelector::print(STREAM &s, const AnyData &v) const + { + s << "values: "; + for (unsigned int i=0; i inline void VectorSelector::print(STREAM &s, const NamedData &v) const diff --git a/deal.II/include/deal.II/meshworker/vector_selector.templates.h b/deal.II/include/deal.II/meshworker/vector_selector.templates.h index dc0e6fbbff..f367c3d06d 100644 --- a/deal.II/include/deal.II/meshworker/vector_selector.templates.h +++ b/deal.II/include/deal.II/meshworker/vector_selector.templates.h @@ -40,6 +40,15 @@ namespace MeshWorker {} + template + void + VectorDataBase::initialize(const AnyData &d) + { + this->data = d; + VectorSelector::initialize(d); + } + + template void VectorDataBase::fill( @@ -93,7 +102,7 @@ namespace MeshWorker void VectorData::initialize(const NamedData &d) { - data = d; + this->data = d; VectorSelector::initialize(d); } @@ -103,8 +112,8 @@ namespace MeshWorker VectorData::initialize(const VECTOR *v, const std::string &name) { SmartPointer > p = v; - data.add(p, name); - VectorSelector::initialize(data); + this->data.add(p, name); + VectorSelector::initialize(this->data); } @@ -125,25 +134,26 @@ namespace MeshWorker AssertDimension(gradients.size(), this->n_gradients()); AssertDimension(hessians.size(), this->n_hessians()); + const AnyData& data = this->data; for (unsigned int i=0; in_values(); ++i) { - const VECTOR &src = *data(this->value_index(i)); + const VECTOR* src = data.read_ptr(this->value_index(i)); VectorSlice > > dst(values[i], component, n_comp); - fe.get_function_values(src, make_slice(index, start, size), dst, true); + fe.get_function_values(*src, make_slice(index, start, size), dst, true); } for (unsigned int i=0; in_gradients(); ++i) { - const VECTOR &src = *data(this->gradient_index(i)); + const VECTOR* src = data.read_ptr(this->gradient_index(i)); VectorSlice > > > dst(gradients[i], component, n_comp); - fe.get_function_gradients(src, make_slice(index, start, size), dst, true); + fe.get_function_gradients(*src, make_slice(index, start, size), dst, true); } for (unsigned int i=0; in_hessians(); ++i) { - const VECTOR &src = *data(this->hessian_index(i)); + const VECTOR* src = data.read_ptr(this->hessian_index(i)); VectorSlice > > > dst(hessians[i], component, n_comp); - fe.get_function_hessians(src, make_slice(index, start, size), dst, true); + fe.get_function_hessians(*src, make_slice(index, start, size), dst, true); } } @@ -153,7 +163,7 @@ namespace MeshWorker VectorData::memory_consumption () const { std::size_t mem = VectorSelector::memory_consumption(); - mem += sizeof (data); + mem += sizeof (this->data); return mem; } @@ -166,8 +176,8 @@ namespace MeshWorker template MGVectorData::MGVectorData(const VectorSelector &s) - : - VectorDataBase(s) + : + VectorData(s) {} @@ -175,7 +185,7 @@ namespace MeshWorker void MGVectorData::initialize(const NamedData*> &d) { - data = d; + this->data = d; VectorSelector::initialize(d); } @@ -186,14 +196,14 @@ namespace MeshWorker { SmartPointer, MGVectorData > p = v; - data.add(p, name); - VectorSelector::initialize(data); + this->data.add(p, name); + VectorSelector::initialize(this->data); } template void - MGVectorData::mg_fill( + VectorData::mg_fill( std::vector > > &values, std::vector > > > &gradients, std::vector > > > &hessians, @@ -209,37 +219,28 @@ namespace MeshWorker AssertDimension(gradients.size(), this->n_gradients()); AssertDimension(hessians.size(), this->n_hessians()); + const AnyData& data = this->data; for (unsigned int i=0; in_values(); ++i) { - const VECTOR &src = (*data(this->value_index(i)))[level]; + const MGLevelObject* src = data.read_ptr >(this->value_index(i)); VectorSlice > > dst(values[i], component, n_comp); - fe.get_function_values(src, make_slice(index, start, size), dst, true); + fe.get_function_values((*src)[level], make_slice(index, start, size), dst, true); } for (unsigned int i=0; in_gradients(); ++i) { - const VECTOR &src = (*data(this->value_index(i)))[level]; + const MGLevelObject* src = data.read_ptr >(this->value_index(i)); VectorSlice > > > dst(gradients[i], component, n_comp); - fe.get_function_gradients(src, make_slice(index, start, size), dst, true); + fe.get_function_gradients((*src)[level], make_slice(index, start, size), dst, true); } for (unsigned int i=0; in_hessians(); ++i) { - const VECTOR &src = (*data(this->value_index(i)))[level]; + const MGLevelObject* src = data.read_ptr >(this->value_index(i)); VectorSlice > > > dst(hessians[i], component, n_comp); - fe.get_function_hessians(src, make_slice(index, start, size), dst, true); + fe.get_function_hessians((*src)[level], make_slice(index, start, size), dst, true); } } - - - template - std::size_t - MGVectorData::memory_consumption () const - { - std::size_t mem = VectorSelector::memory_consumption(); - mem += sizeof (data); - return mem; - } } DEAL_II_NAMESPACE_CLOSE diff --git a/deal.II/include/deal.II/numerics/data_out.h b/deal.II/include/deal.II/numerics/data_out.h index 64bb199e4c..0bc554148c 100644 --- a/deal.II/include/deal.II/numerics/data_out.h +++ b/deal.II/include/deal.II/numerics/data_out.h @@ -121,9 +121,20 @@ namespace internal * the sequence of cells to be converted into patches itself, but relies * on the two functions first_cell() and next_cell(). By default, they * return the first active cell, and the next active cell, respectively. - * Since they are @p virtual functions, you may overload them to select other - * cells for output. If cells are not active, interpolated values are taken - * for output instead of the exact values on active cells. + * Since they are @p virtual functions, you can write your own class + * derived from DataOut in which you overload these two functions to select other + * cells for output. This may, for example, include only cells that are + * in parts of a domain (e.g., if you don't care about the solution elsewhere, + * think for example a buffer region in which you attenuate outgoing waves + * in the Perfectly Matched Layer method) or if you don't want output to + * be generated at all levels of an adaptively refined mesh because this + * creates too much data (in this case, the set of cells returned by + * your implementations of first_cell() and next_cell() will include + * non-active cells, and DataOut::build_patches() will simply take + * interpolated values of the solution instead of the exact values on these + * cells children for output). Once you derive your own class, you would + * just create an object of this type instead of an object of type DataOut, + * and everything else will remain the same. * * The two functions are not constant, so you may store information within * your derived class about the last accessed cell. This is useful if the diff --git a/deal.II/include/deal.II/numerics/dof_output_operator.h b/deal.II/include/deal.II/numerics/dof_output_operator.h index 70f89b5467..ed4bf51c14 100644 --- a/deal.II/include/deal.II/numerics/dof_output_operator.h +++ b/deal.II/include/deal.II/numerics/dof_output_operator.h @@ -19,10 +19,13 @@ #define __deal2__dof_output_operator_h #include +#include #include +#include #include #include #include +#include #include @@ -30,18 +33,40 @@ DEAL_II_NAMESPACE_OPEN namespace Algorithms { + /** + * An output operator writing a separate file in each step and + * writing the vectors as finite element functions with respect to a + * given DoFHandler. + */ template class DoFOutputOperator : public OutputOperator { public: + /* + * Constructor. The filename is the common base name of + * all files and the argument digits should be the number + * of digits of the highest number in the sequence. File names by + * default have the form "outputNN" with NNN the number set by the + * last step command. Numbers with less digits are filled with + * zeros from the left. + */ + DoFOutputOperator (const std::string filename_base = std::string("output"), + const unsigned int digits = 3); + + void parse_parameters(ParameterHandler ¶m); void initialize (DoFHandler &dof_handler); virtual OutputOperator & - operator<<(const NamedData &vectors); + operator << (const AnyData &vectors); private: SmartPointer, DoFOutputOperator > dof; + + const std::string filename_base; + const unsigned int digits; + + DataOut out; }; template diff --git a/deal.II/include/deal.II/numerics/dof_output_operator.templates.h b/deal.II/include/deal.II/numerics/dof_output_operator.templates.h index b9529a01c4..5bb629a8e1 100644 --- a/deal.II/include/deal.II/numerics/dof_output_operator.templates.h +++ b/deal.II/include/deal.II/numerics/dof_output_operator.templates.h @@ -16,28 +16,53 @@ #include -#include DEAL_II_NAMESPACE_OPEN namespace Algorithms { + template + DoFOutputOperator::DoFOutputOperator ( + const std::string filename_base, + const unsigned int digits) + : + filename_base(filename_base), + digits(digits) + { + out.set_default_format(DataOutBase::gnuplot); + } + + + template + void + DoFOutputOperator::parse_parameters(ParameterHandler ¶m) + { + out.parse_parameters(param); + } + template OutputOperator & DoFOutputOperator::operator<<( - const NamedData &vectors) + const AnyData &data) { Assert ((dof!=0), ExcNotInitialized()); - DataOut out; out.attach_dof_handler (*dof); - out.add_data_vector (*vectors(vectors.find("solution")), "solution"); - out.add_data_vector (*vectors(vectors.find("update")), "update"); - out.add_data_vector (*vectors(vectors.find("residual")), "residual"); + for (unsigned int i=0; i(i); + if (p!=0) + { + out.add_data_vector (*p, data.name(i)); + } + } std::ostringstream streamOut; - streamOut << "Newton_" << std::setw(3) << std::setfill('0') << this->step; + streamOut << filename_base + << std::setw(digits) << std::setfill('0') << this->step + << out.default_suffix(); std::ofstream out_filename (streamOut.str().c_str()); - out.build_patches (2); - out.write_gnuplot (out_filename); + out.build_patches (); + out.write (out_filename); + out.clear (); return *this; } } diff --git a/deal.II/source/base/parameter_handler.cc b/deal.II/source/base/parameter_handler.cc index 29446174a1..a81e36349b 100644 --- a/deal.II/source/base/parameter_handler.cc +++ b/deal.II/source/base/parameter_handler.cc @@ -2066,7 +2066,12 @@ ParameterHandler::print_parameters_section (std::ostream &out, const std::string value = p->second.get("value"); // print name - out << "\\item {\\it Parameter name:} {\\tt " << demangle(p->first) << "}\n\n" + out << "\\item {\\it Parameter name:} {\\tt " << demangle(p->first) << "}\n" + << "\\phantomsection\\label{parameters:"; + for (unsigned int i=0; ifirst); + out << "}\n\n" << std::endl; out << "\\index[prmindex]{" diff --git a/deal.II/source/lac/petsc_solver.cc b/deal.II/source/lac/petsc_solver.cc index 4524e91c37..d9c95f30a4 100644 --- a/deal.II/source/lac/petsc_solver.cc +++ b/deal.II/source/lac/petsc_solver.cc @@ -1,7 +1,7 @@ // --------------------------------------------------------------------- // $Id$ // -// Copyright (C) 2004 - 2013 by the deal.II authors +// Copyright (C) 2004 - 2014 by the deal.II authors // // This file is part of the deal.II library. // @@ -736,7 +736,11 @@ namespace PETScWrappers * factorization here we start to see differences with the base * class solve function */ +#if DEAL_II_PETSC_VERSION_GTE(3,2,0) ierr = PCFactorSetMatSolverPackage (solver_data->pc, MATSOLVERMUMPS); +#else + ierr = PCFactorSetMatSolverPackage (solver_data->pc, MAT_SOLVER_MUMPS); +#endif AssertThrow (ierr == 0, ExcPETScError (ierr)); /** diff --git a/deal.II/source/numerics/matrix_tools.cc b/deal.II/source/numerics/matrix_tools.cc index 6a726729c1..dbb7d854f6 100644 --- a/deal.II/source/numerics/matrix_tools.cc +++ b/deal.II/source/numerics/matrix_tools.cc @@ -1,7 +1,7 @@ // --------------------------------------------------------------------- // $Id$ // -// Copyright (C) 1998 - 2013 by the deal.II authors +// Copyright (C) 1998 - 2014 by the deal.II authors // // This file is part of the deal.II library. // @@ -2334,75 +2334,83 @@ namespace MatrixTools Assert (matrix.n() == solution.size(), ExcDimensionMismatch(matrix.n(), solution.size())); - // if no boundary values are to be applied - // simply return - if (boundary_values.size() == 0) - return; - - const std::pair local_range - = matrix.local_range(); - Assert (local_range == right_hand_side.local_range(), - ExcInternalError()); - Assert (local_range == solution.local_range(), - ExcInternalError()); - - // determine the first nonzero diagonal - // entry from within the part of the - // matrix that we can see. if we can't - // find such an entry, take one - PetscScalar average_nonzero_diagonal_entry = 1; - for (types::global_dof_index i=local_range.first; i constrained_rows; - for (std::map::const_iterator - dof = boundary_values.begin(); - dof != boundary_values.end(); - ++dof) - if ((dof->first >= local_range.first) && - (dof->first < local_range.second)) - constrained_rows.push_back (dof->first); - - // then eliminate these rows and set - // their diagonal entry to what we have - // determined above. note that for petsc - // matrices interleaving read with write - // operations is very expensive. thus, we - // here always replace the diagonal - // element, rather than first checking - // whether it is nonzero and in that case - // preserving it. this is different from - // the case of deal.II sparse matrices - // treated in the other functions. - matrix.clear_rows (constrained_rows, average_nonzero_diagonal_entry); + // if no boundary values are to be applied, then + // jump straight to the compress() calls that we still have + // to perform because they are collective operations + if (boundary_values.size() > 0) + { + const std::pair local_range + = matrix.local_range(); + Assert (local_range == right_hand_side.local_range(), + ExcInternalError()); + Assert (local_range == solution.local_range(), + ExcInternalError()); + + // determine the first nonzero diagonal + // entry from within the part of the + // matrix that we can see. if we can't + // find such an entry, take one + PetscScalar average_nonzero_diagonal_entry = 1; + for (types::global_dof_index i=local_range.first; i indices; - std::vector solution_values; - for (std::map::const_iterator - dof = boundary_values.begin(); - dof != boundary_values.end(); - ++dof) - if ((dof->first >= local_range.first) && - (dof->first < local_range.second)) - { - indices.push_back (dof->first); - solution_values.push_back (dof->second); - } - solution.set (indices, solution_values); + // figure out which rows of the matrix we + // have to eliminate on this processor + std::vector constrained_rows; + for (std::map::const_iterator + dof = boundary_values.begin(); + dof != boundary_values.end(); + ++dof) + if ((dof->first >= local_range.first) && + (dof->first < local_range.second)) + constrained_rows.push_back (dof->first); - // now also set appropriate values for - // the rhs - for (unsigned int i=0; i indices; + std::vector solution_values; + for (std::map::const_iterator + dof = boundary_values.begin(); + dof != boundary_values.end(); + ++dof) + if ((dof->first >= local_range.first) && + (dof->first < local_range.second)) + { + indices.push_back (dof->first); + solution_values.push_back (dof->second); + } + solution.set (indices, solution_values); - right_hand_side.set (indices, solution_values); + // now also set appropriate values for + // the rhs + for (unsigned int i=0; i constrained_rows; + matrix.clear_rows (constrained_rows, 1.); + } + // clean up matrix.compress (); solution.compress (VectorOperation::insert); diff --git a/tests/algorithms/newton_01.output b/tests/algorithms/newton_01.output index 38c1c23723..48e97a1986 100644 --- a/tests/algorithms/newton_01.output +++ b/tests/algorithms/newton_01.output @@ -1,12 +1,47 @@ DEAL:Newton::Check 0 98.0000 DEAL:Newton::Starting value 98.0000 +DEAL:Newton::Step 1 +DEAL:Newton::solution 10.0000 +DEAL:Newton::update 4.90000 +DEAL:Newton::residual 98.0000 +DEAL:Newton:: DEAL:Newton::Check 1 24.0100 +DEAL:Newton::Step 2 +DEAL:Newton::solution 5.10000 +DEAL:Newton::update 2.35392 +DEAL:Newton::residual 24.0100 +DEAL:Newton:: DEAL:Newton::Check 2 5.54095 +DEAL:Newton::Step 3 +DEAL:Newton::solution 2.74608 +DEAL:Newton::update 1.00888 +DEAL:Newton::residual 5.54095 +DEAL:Newton:: DEAL:Newton::Check 3 1.01785 +DEAL:Newton::Step 4 +DEAL:Newton::solution 1.73719 +DEAL:Newton::update 0.292957 +DEAL:Newton::residual 1.01785 +DEAL:Newton:: DEAL:Newton::Check 4 0.0858237 +DEAL:Newton::Step 5 +DEAL:Newton::solution 1.44424 +DEAL:Newton::update 0.0297124 +DEAL:Newton::residual 0.0858237 +DEAL:Newton:: DEAL:Newton::Check 5 0.000882829 +DEAL:Newton::Step 6 +DEAL:Newton::solution 1.41453 +DEAL:Newton::update 0.000312058 +DEAL:Newton::residual 0.000882829 +DEAL:Newton:: DEAL:Newton::Check 6 9.73804e-08 +DEAL:Newton::Step 7 +DEAL:Newton::solution 1.41421 +DEAL:Newton::update 3.44292e-08 +DEAL:Newton::residual 9.73804e-08 +DEAL:Newton:: DEAL:Newton::Check 7 0 DEAL:Newton::Convergence step 7 value 0 DEAL:: square root 1.41421 diff --git a/tests/algorithms/newton_01_compat.output b/tests/algorithms/newton_01_compat.output index 38c1c23723..48e97a1986 100644 --- a/tests/algorithms/newton_01_compat.output +++ b/tests/algorithms/newton_01_compat.output @@ -1,12 +1,47 @@ DEAL:Newton::Check 0 98.0000 DEAL:Newton::Starting value 98.0000 +DEAL:Newton::Step 1 +DEAL:Newton::solution 10.0000 +DEAL:Newton::update 4.90000 +DEAL:Newton::residual 98.0000 +DEAL:Newton:: DEAL:Newton::Check 1 24.0100 +DEAL:Newton::Step 2 +DEAL:Newton::solution 5.10000 +DEAL:Newton::update 2.35392 +DEAL:Newton::residual 24.0100 +DEAL:Newton:: DEAL:Newton::Check 2 5.54095 +DEAL:Newton::Step 3 +DEAL:Newton::solution 2.74608 +DEAL:Newton::update 1.00888 +DEAL:Newton::residual 5.54095 +DEAL:Newton:: DEAL:Newton::Check 3 1.01785 +DEAL:Newton::Step 4 +DEAL:Newton::solution 1.73719 +DEAL:Newton::update 0.292957 +DEAL:Newton::residual 1.01785 +DEAL:Newton:: DEAL:Newton::Check 4 0.0858237 +DEAL:Newton::Step 5 +DEAL:Newton::solution 1.44424 +DEAL:Newton::update 0.0297124 +DEAL:Newton::residual 0.0858237 +DEAL:Newton:: DEAL:Newton::Check 5 0.000882829 +DEAL:Newton::Step 6 +DEAL:Newton::solution 1.41453 +DEAL:Newton::update 0.000312058 +DEAL:Newton::residual 0.000882829 +DEAL:Newton:: DEAL:Newton::Check 6 9.73804e-08 +DEAL:Newton::Step 7 +DEAL:Newton::solution 1.41421 +DEAL:Newton::update 3.44292e-08 +DEAL:Newton::residual 9.73804e-08 +DEAL:Newton:: DEAL:Newton::Check 7 0 DEAL:Newton::Convergence step 7 value 0 DEAL:: square root 1.41421 diff --git a/tests/algorithms/theta_01.cc b/tests/algorithms/theta_01.cc new file mode 100644 index 0000000000..bb18059d1b --- /dev/null +++ b/tests/algorithms/theta_01.cc @@ -0,0 +1,146 @@ +// --------------------------------------------------------------------- +// $Id: theta_timestepping.cc 30050 2013-07-18 20:03:31Z maier $ +// +// Copyright (C) 2005 - 2013 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + + +// See documentation of ThetaTimestepping for documentation of this example + +#include +#include +#include + +#include +#include + +#include + +using namespace dealii; +using namespace Algorithms; + + +class Explicit + : public Operator > +{ +public: + Explicit(const FullMatrix &matrix); + void operator() (AnyData &out, const AnyData &in); +private: + SmartPointer, Explicit> matrix; + FullMatrix m; +}; + + +class Implicit + : public Operator > +{ +public: + Implicit(const FullMatrix &matrix); + void operator() (AnyData &out, const AnyData &in); +private: + SmartPointer, Implicit> matrix; + FullMatrix m; +}; + +// End of declarations + +int main() +{ + std::string logname = "output"; + std::ofstream logfile(logname.c_str()); + deallog.attach(logfile); + deallog.depth_console(0); + deallog.threshold_double(1.e-10); + + FullMatrix matrix(2); + matrix(0,0) = 0.; + matrix(1,1) = 0.; + matrix(0,1) = numbers::PI; + matrix(1,0) = -numbers::PI; + + OutputOperator > out; + out.initialize_stream(logfile); + + Explicit op_explicit(matrix); + Implicit op_implicit(matrix); + ThetaTimestepping > solver(op_explicit, op_implicit); + solver.timestep_control().start_step(0.1); + //solver.set_output(out); + + Vector value(2); + value(0) = 1.; + AnyData indata; + AnyData outdata; + Vector *p = &value; + outdata.add(p, "value"); + deallog << "Initial: " << value(0) << ' ' << value(1) << std::endl; + solver.notify(Events::initial); + solver(outdata, indata); + deallog << "Result: " << value(0) << ' ' << value(1) + << " Norm " << value.l2_norm() << std::endl; +} + + +Explicit::Explicit(const FullMatrix &M) + : + matrix(&M) +{ + m.reinit(M.m(), M.n()); +} + + +void +Explicit::operator() (AnyData &out, const AnyData &in) +{ + const double* step = in.read_ptr("Timestep"); + + if (this->notifications.test(Events::initial) || this->notifications.test(Events::new_timestep_size)) + { + m.equ(-*step, *matrix); + for (unsigned int i=0; inotifications.clear(); + unsigned int i = in.find("Previous iterate"); + m.vmult(*out.entry*>(0), *in.read_ptr >(i)); +} + + +Implicit::Implicit(const FullMatrix &M) + : + matrix(&M) +{ + m.reinit(M.m(), M.n()); +} + + +void +Implicit::operator() (AnyData &out, const AnyData &in) +{ + const double * step = in.read_ptr("Timestep"); + + if (this->notifications.test(Events::initial) || this->notifications.test(Events::new_timestep_size)) + { + m.equ(*step, *matrix); + for (unsigned int i=0; inotifications.clear(); + + unsigned int i = in.find("Previous time"); + m.vmult(*out.entry*>(0), *in.read_ptr >(i)); +} + + diff --git a/tests/algorithms/theta_01.output b/tests/algorithms/theta_01.output new file mode 100644 index 0000000000..7c4d9ddfc7 --- /dev/null +++ b/tests/algorithms/theta_01.output @@ -0,0 +1,13 @@ + +DEAL::Initial: 1.00000 0 +DEAL:Theta::Time step:0.100000 +DEAL:Theta::Time step:0.200000 +DEAL:Theta::Time step:0.300000 +DEAL:Theta::Time step:0.400000 +DEAL:Theta::Time step:0.500000 +DEAL:Theta::Time step:0.600000 +DEAL:Theta::Time step:0.700000 +DEAL:Theta::Time step:0.800000 +DEAL:Theta::Time step:0.900000 +DEAL:Theta::Time step:1.00000 +DEAL::Result: -0.999676 0.0254599 Norm 1.00000 diff --git a/tests/bits/parameter_handler_4.output b/tests/bits/parameter_handler_4.output index ad00d8923e..7262fc5e81 100644 --- a/tests/bits/parameter_handler_4.output +++ b/tests/bits/parameter_handler_4.output @@ -9,6 +9,7 @@ \begin{itemize} \item {\it Parameter name:} {\tt double} +\phantomsection\label{parameters:Testing/double} \index[prmindex]{double} @@ -24,6 +25,7 @@ {\it Possible values:} [Double -1.79769e+308...1.79769e+308 (inclusive)] \item {\it Parameter name:} {\tt int} +\phantomsection\label{parameters:Testing/int} \index[prmindex]{int} @@ -36,6 +38,7 @@ {\it Possible values:} [Integer range -2147483648...2147483647 (inclusive)] \item {\it Parameter name:} {\tt string list} +\phantomsection\label{parameters:Testing/string list} \index[prmindex]{string list} diff --git a/tests/bits/parameter_handler_4a.output b/tests/bits/parameter_handler_4a.output index a3e8ca1109..6dcf0739c9 100644 --- a/tests/bits/parameter_handler_4a.output +++ b/tests/bits/parameter_handler_4a.output @@ -9,6 +9,7 @@ \begin{itemize} \item {\it Parameter name:} {\tt double} +\phantomsection\label{parameters:Testing/double} \index[prmindex]{double} @@ -24,6 +25,7 @@ {\it Possible values:} [Double -1.79769e+308...1.79769e+308 (inclusive)] \item {\it Parameter name:} {\tt int} +\phantomsection\label{parameters:Testing/int} \index[prmindex]{int} @@ -36,6 +38,7 @@ {\it Possible values:} [Integer range -2147483648...2147483647 (inclusive)] \item {\it Parameter name:} {\tt string list} +\phantomsection\label{parameters:Testing/string list} \index[prmindex]{string list} @@ -59,6 +62,7 @@ \begin{itemize} \item {\it Parameter name:} {\tt double 2} +\phantomsection\label{parameters:Testing/Testing 2/double 2} \index[prmindex]{double 2} @@ -74,6 +78,7 @@ {\it Possible values:} [Double -1.79769e+308...1.79769e+308 (inclusive)] \item {\it Parameter name:} {\tt int 2} +\phantomsection\label{parameters:Testing/Testing 2/int 2} \index[prmindex]{int 2} @@ -86,6 +91,7 @@ {\it Possible values:} [Integer range -2147483648...2147483647 (inclusive)] \item {\it Parameter name:} {\tt string list 2} +\phantomsection\label{parameters:Testing/Testing 2/string list 2} \index[prmindex]{string list 2} diff --git a/tests/lac/gmres_eigenvalues.cc b/tests/lac/gmres_eigenvalues.cc index 94ff53077d..f638e6ae7d 100644 --- a/tests/lac/gmres_eigenvalues.cc +++ b/tests/lac/gmres_eigenvalues.cc @@ -28,7 +28,7 @@ template void test (unsigned int variant) { - const unsigned int n = variant < 3 ? 64 : 16; + const unsigned int n = variant % 2 == 0 ? 64 : 16; Vector rhs(n), sol(n); rhs = 1.; @@ -42,7 +42,7 @@ void test (unsigned int variant) matrix(i,i) = (i+1); else if (variant == 1) for (unsigned int i=0; i(0); diff --git a/tests/lac/gmres_eigenvalues.output b/tests/lac/gmres_eigenvalues.output index ff01b9bef3..4bcadb4367 100644 --- a/tests/lac/gmres_eigenvalues.output +++ b/tests/lac/gmres_eigenvalues.output @@ -1,17 +1,17 @@ -DEAL:double:0:GMRES::Starting value 8.00 +DEAL:double:0:GMRES::Starting value 8.0 DEAL:double:0:GMRES::Convergence step 56 value 0 -DEAL:double:0:GMRES::Eigenvalue estimate: (1.00,0.00) (2.00,0.00) (3.00,0.00) (4.00,0.00) (5.00,0.00) (6.00,0.00) (7.00,0.00) (8.00,0.00) (9.00,0.00) (10.0,0.00) (11.0,0.00) (12.0,0.00) (13.0,0.00) (14.0,0.00) (15.0,0.00) (16.0,0.00) (17.1,0.00) (18.2,0.00) (19.5,0.00) (20.7,0.00) (22.0,0.00) (23.3,0.00) (24.7,0.00) (26.0,0.00) (27.5,0.00) (28.9,0.00) (30.3,0.00) (31.8,0.00) (33.2,0.00) (34.7,0.00) (36.1,0.00) (37.6,0.00) (39.0,0.00) (40.4,0.00) (41.7,0.00) (43.0,0.00) (44.3,0.00) (45.5,0.00) (46.7,0.00) (47.9,0.00) (49.0,0.00) (50.0,0.00) (51.0,0.00) (52.0,0.00) (53.0,0.00) (54.0,0.00) (55.0,0.00) (56.0,0.00) (57.0,0.00) (58.0,0.00) (59.0,0.00) (60.0,0.00) (61.0,0.00) (62.0,0.00) (63.0,0.00) (64.0,0.00) -DEAL:double:0:cg::Starting value 8.00 +DEAL:double:0:GMRES::Eigenvalue estimate: (1.0,0.0) (2.0,0.0) (3.0,0.0) (4.0,0.0) (5.0,0.0) (6.0,0.0) (7.0,0.0) (8.0,0.0) (9.0,0.0) (10.,0.0) (11.,0.0) (12.,0.0) (13.,0.0) (14.,0.0) (15.,0.0) (16.,0.0) (17.,0.0) (18.,0.0) (19.,0.0) (21.,0.0) (22.,0.0) (23.,0.0) (25.,0.0) (26.,0.0) (27.,0.0) (29.,0.0) (30.,0.0) (32.,0.0) (33.,0.0) (35.,0.0) (36.,0.0) (38.,0.0) (39.,0.0) (40.,0.0) (42.,0.0) (43.,0.0) (44.,0.0) (46.,0.0) (47.,0.0) (48.,0.0) (49.,0.0) (50.,0.0) (51.,0.0) (52.,0.0) (53.,0.0) (54.,0.0) (55.,0.0) (56.,0.0) (57.,0.0) (58.,0.0) (59.,0.0) (60.,0.0) (61.,0.0) (62.,0.0) (63.,0.0) (64.,0.0) +DEAL:double:0:cg::Starting value 8.0 DEAL:double:0:cg::Convergence step 56 value 0 -DEAL:double:0:cg:: 1.00 2.00 3.00 4.00 5.00 6.00 7.00 8.00 9.00 10.0 11.0 12.0 13.0 14.0 15.0 16.1 17.3 18.4 19.7 21.0 22.3 23.7 25.1 26.6 28.0 29.5 31.0 32.5 34.0 35.5 37.0 38.4 39.9 41.3 42.7 44.0 45.3 46.6 47.7 48.9 50.0 51.0 52.0 53.0 54.0 55.0 56.0 57.0 58.0 59.0 60.0 61.0 62.0 63.0 64.0 -DEAL:double:1:GMRES::Starting value 8.00 -DEAL:double:1:GMRES::Convergence step 64 value 4.41e-10 -DEAL:double:1:GMRES::Eigenvalue estimate: (1.00,0.00) (16.0,0.00) (81.0,0.00) (256.,0.00) (625.,0.00) (1.30e+03,0.00) (2.40e+03,0.00) (4.10e+03,0.00) (6.56e+03,0.00) (1.00e+04,0.00) (1.46e+04,0.00) (2.07e+04,0.00) (2.86e+04,0.00) (3.84e+04,0.00) (5.06e+04,0.00) (6.55e+04,0.00) (8.35e+04,0.00) (1.05e+05,0.00) (1.30e+05,0.00) (1.60e+05,0.00) (1.94e+05,0.00) (2.34e+05,0.00) (2.80e+05,0.00) (3.32e+05,0.00) (3.91e+05,0.00) (4.57e+05,0.00) (5.31e+05,0.00) (6.15e+05,0.00) (7.07e+05,0.00) (8.10e+05,0.00) (9.24e+05,0.00) (1.05e+06,0.00) (1.19e+06,0.00) (1.34e+06,0.00) (1.50e+06,0.00) (1.68e+06,0.00) (1.87e+06,0.00) (2.09e+06,0.00) (2.31e+06,0.00) (2.56e+06,0.00) (2.83e+06,0.00) (3.11e+06,0.00) (3.42e+06,0.00) (3.75e+06,0.00) (4.10e+06,0.00) (4.48e+06,0.00) (4.88e+06,0.00) (5.31e+06,0.00) (5.76e+06,0.00) (6.25e+06,0.00) (6.77e+06,0.00) (7.31e+06,0.00) (7.89e+06,0.00) (8.50e+06,0.00) (9.15e+06,0.00) (9.83e+06,0.00) (1.06e+07,0.00) (1.13e+07,0.00) (1.21e+07,0.00) (1.30e+07,0.00) (1.38e+07,0.00) (1.48e+07,0.00) (1.58e+07,0.00) (1.68e+07,0.00) -DEAL:double:2:GMRES::Starting value 8.00 +DEAL:double:0:cg:: 1.0 2.0 3.0 4.0 5.0 6.0 7.0 8.0 9.0 10. 11. 12. 13. 14. 15. 16. 17. 18. 20. 21. 22. 24. 25. 27. 28. 30. 31. 33. 34. 35. 37. 38. 40. 41. 43. 44. 45. 47. 48. 49. 50. 51. 52. 53. 54. 55. 56. 57. 58. 59. 60. 61. 62. 63. 64. +DEAL:double:1:GMRES::Starting value 4.0 +DEAL:double:1:GMRES::Convergence step 16 value 0 +DEAL:double:1:GMRES::Eigenvalue estimate: (1.0,0.0) (16.,0.0) (81.,0.0) (2.6e+02,0.0) (6.3e+02,0.0) (1.3e+03,0.0) (2.4e+03,0.0) (4.1e+03,0.0) (6.6e+03,0.0) (1.0e+04,0.0) (1.5e+04,0.0) (2.1e+04,0.0) (2.9e+04,0.0) (3.8e+04,0.0) (5.1e+04,0.0) (6.6e+04,0.0) +DEAL:double:2:GMRES::Starting value 8.0 DEAL:double:2:GMRES::Convergence step 64 value 0 -DEAL:double:2:GMRES::Eigenvalue estimate: (-63.0,0.00) (-61.0,0.00) (-59.0,0.00) (-57.0,0.00) (-55.0,0.00) (-53.0,0.00) (-51.0,0.00) (-49.0,0.00) (-47.0,0.00) (-45.0,0.00) (-43.0,0.00) (-41.0,0.00) (-39.0,0.00) (-37.0,0.00) (-35.0,0.00) (-33.0,0.00) (-31.0,0.00) (-29.0,0.00) (-27.0,0.00) (-25.0,0.00) (-23.0,0.00) (-21.0,0.00) (-19.0,0.00) (-17.0,0.00) (-15.0,0.00) (-13.0,0.00) (-11.0,0.00) (-9.00,0.00) (-7.00,0.00) (-5.00,0.00) (-3.00,0.00) (-1.00,0.00) (2.00,0.00) (4.00,0.00) (6.00,0.00) (8.00,0.00) (10.0,0.00) (12.0,0.00) (14.0,0.00) (16.0,0.00) (18.0,0.00) (20.0,0.00) (22.0,0.00) (24.0,0.00) (26.0,0.00) (28.0,0.00) (30.0,0.00) (32.0,0.00) (34.0,0.00) (36.0,0.00) (38.0,0.00) (40.0,0.00) (42.0,0.00) (44.0,0.00) (46.0,0.00) (48.0,0.00) (50.0,0.00) (52.0,0.00) (54.0,0.00) (56.0,0.00) (58.0,0.00) (60.0,0.00) (62.0,0.00) (64.0,0.00) -DEAL:double:3:GMRES::Starting value 4.00 +DEAL:double:2:GMRES::Eigenvalue estimate: (-63.,0.0) (-61.,0.0) (-59.,0.0) (-57.,0.0) (-55.,0.0) (-53.,0.0) (-51.,0.0) (-49.,0.0) (-47.,0.0) (-45.,0.0) (-43.,0.0) (-41.,0.0) (-39.,0.0) (-37.,0.0) (-35.,0.0) (-33.,0.0) (-31.,0.0) (-29.,0.0) (-27.,0.0) (-25.,0.0) (-23.,0.0) (-21.,0.0) (-19.,0.0) (-17.,0.0) (-15.,0.0) (-13.,0.0) (-11.,0.0) (-9.0,0.0) (-7.0,0.0) (-5.0,0.0) (-3.0,0.0) (-1.0,0.0) (2.0,0.0) (4.0,0.0) (6.0,0.0) (8.0,0.0) (10.,0.0) (12.,0.0) (14.,0.0) (16.,0.0) (18.,0.0) (20.,0.0) (22.,0.0) (24.,0.0) (26.,0.0) (28.,0.0) (30.,0.0) (32.,0.0) (34.,0.0) (36.,0.0) (38.,0.0) (40.,0.0) (42.,0.0) (44.,0.0) (46.,0.0) (48.,0.0) (50.,0.0) (52.,0.0) (54.,0.0) (56.,0.0) (58.,0.0) (60.,0.0) (62.,0.0) (64.,0.0) +DEAL:double:3:GMRES::Starting value 4.0 DEAL:double:3:GMRES::Convergence step 16 value 0 -DEAL:double:3:GMRES::Eigenvalue estimate: (0.272,-1.14) (0.272,1.14) (1.52,-3.81) (1.52,3.81) (2.81,-3.71) (2.81,3.71) (4.50,-4.11) (4.50,4.11) (7.18,-5.38) (7.18,5.38) (10.3,-5.98) (10.3,5.98) (13.9,-5.92) (13.9,5.92) (18.5,-4.07) (18.5,4.07) -DEAL:double:3::Actual eigenvalues: (0.272,-1.14) (0.272,1.14) (1.52,-3.81) (1.52,3.81) (2.81,-3.71) (2.81,3.71) (4.50,-4.11) (4.50,4.11) (7.18,-5.38) (7.18,5.38) (10.3,-5.98) (10.3,5.98) (13.9,-5.92) (13.9,5.92) (18.5,-4.07) (18.5,4.07) +DEAL:double:3:GMRES::Eigenvalue estimate: (0.27,-1.1) (0.27,1.1) (1.5,-3.8) (1.5,3.8) (2.8,-3.7) (2.8,3.7) (4.5,-4.1) (4.5,4.1) (7.2,-5.4) (7.2,5.4) (10.,-6.0) (10.,6.0) (14.,-5.9) (14.,5.9) (19.,-4.1) (19.,4.1) +DEAL:double:3::Actual eigenvalues: (0.27,-1.1) (0.27,1.1) (1.5,-3.8) (1.5,3.8) (2.8,-3.7) (2.8,3.7) (4.5,-4.1) (4.5,4.1) (7.2,-5.4) (7.2,5.4) (10.,-6.0) (10.,6.0) (14.,-5.9) (14.,5.9) (19.,-4.1) (19.,4.1)