From: bangerth Date: Tue, 26 Jun 2007 05:04:09 +0000 (+0000) Subject: Introduce Utilities::fixed_power. Use it in data_out_base. Move another exception... X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=4cf01537f74c7df6e55b594386d0aed9bb4d3cb2;p=dealii-svn.git Introduce Utilities::fixed_power. Use it in data_out_base. Move another exception away from that file. git-svn-id: https://svn.dealii.org/trunk@14806 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/base/include/base/utilities.h b/deal.II/base/include/base/utilities.h index 422d07e2b6..af5ac3472a 100644 --- a/deal.II/base/include/base/utilities.h +++ b/deal.II/base/include/base/utilities.h @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 2005, 2006 by the deal.II authors +// Copyright (C) 2005, 2006, 2007 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -14,6 +14,7 @@ #define __deal2__utilities_h #include +#include #include #include @@ -141,6 +142,22 @@ namespace Utilities const double sigma); + /** + * Calculate a fixed power, provided as a + * template argument, of a number. + * + * This function provides an efficient way + * to calculate things like + * t^N where N is + * a known number at compile time. + * + * Use this function as in + * fixed_power@ (n). + */ + template + T + fixed_power (const T t); + /** * A namespace for utility functions that * probe system properties. @@ -208,6 +225,37 @@ namespace Utilities } } + +// --------------------- inline functions + +namespace Utilities +{ + template + inline + T fixed_power (const T n) + { + Assert (N>0, ExcNotImplemented()); + switch (N) + { + case 1: + return n; + case 2: + return n*n; + case 3: + return n*n*n; + case 4: + return n*n*n*n; + default: + T result = n; + for (unsigned int d=1;d +#include #include #include #include @@ -45,52 +46,54 @@ DEAL_II_NAMESPACE_OPEN -//TODO: Is it reasonable to have undocumented exceptions? -DeclException2 (ExcUnexpectedInput, - std::string, std::string, - << "Unexpected input: expected line\n <" - << arg1 - << ">\nbut got\n <" - << arg2 << ">"); +// we need the following exception from a global function, so can't declare it +// in the usual way inside a class +namespace +{ + DeclException2 (ExcUnexpectedInput, + std::string, std::string, + << "Unexpected input: expected line\n <" + << arg1 + << ">\nbut got\n <" + << arg2 << ">"); +} -DeclException4 (ExcIncompatibleDimensions, - int, int, int, int, - << "Either the dimensions <" << arg1 << "> and <" - << arg2 << "> or the space dimensions <" - << arg3 << "> and <" << arg4 - << "> do not match!"); //----------------------------------------------------------------------// //Auxiliary data //----------------------------------------------------------------------// -static const char* gmv_cell_type[4] = -{ - "", "line 2", "quad 4", "hex 8" -}; -static const char* ucd_cell_type[4] = +namespace { - "", "line", "quad", "hex" -}; + + static const char* gmv_cell_type[4] = + { + "", "line 2", "quad 4", "hex 8" + }; -static const char* tecplot_cell_type[4] = -{ - "", "", "quadrilateral", "brick" -}; + static const char* ucd_cell_type[4] = + { + "", "line", "quad", "hex" + }; + + static const char* tecplot_cell_type[4] = + { + "", "", "quadrilateral", "brick" + }; #ifdef DEAL_II_HAVE_TECPLOT -static unsigned int tecplot_binary_cell_type[4] = -{ - 0, 0, 1, 3 -}; + static unsigned int tecplot_binary_cell_type[4] = + { + 0, 0, 1, 3 + }; #endif -static unsigned int vtk_cell_type[4] = -{ - 0, 3, 9, 12 -}; + static unsigned int vtk_cell_type[4] = + { + 0, 3, 9, 12 + }; //----------------------------------------------------------------------// //Auxiliary functions @@ -99,104 +102,80 @@ static unsigned int vtk_cell_type[4] = //linearly at the point (xstep, ystep, zstep)*1./n_subdivisions. //If the popints are saved in the patch->data member, return the //saved point instead -template -inline -void -compute_node( - Point& node, - const DataOutBase::Patch* patch, - const unsigned int xstep, - const unsigned int ystep, - const unsigned int zstep, - const unsigned int n_subdivisions) -{ - if (patch->points_are_available) - { - unsigned int point_no=0; - // note: switch without break ! - switch (dim) - { - case 3: - Assert (zstepdata(patch->data.size(0)-spacedim+d,point_no); - } - else - { - // perform a dim-linear interpolation - const double stepsize=1./n_subdivisions, - xfrac=xstep*stepsize; - - node = (patch->vertices[1] * xfrac) + (patch->vertices[0] * (1-xfrac)); - if (dim>1) - { - const double yfrac=ystep*stepsize; - node*= 1-yfrac; - node += ((patch->vertices[3] * xfrac) + (patch->vertices[2] * (1-xfrac))) * yfrac; - if (dim>2) - { - const double zfrac=zstep*stepsize; - node *= (1-zfrac); - node += (((patch->vertices[5] * xfrac) + (patch->vertices[4] * (1-xfrac))) - * (1-yfrac) + - ((patch->vertices[7] * xfrac) + (patch->vertices[6] * (1-xfrac))) - * yfrac) * zfrac; - } - } - } -} + template + inline + void + compute_node( + Point& node, + const DataOutBase::Patch* patch, + const unsigned int xstep, + const unsigned int ystep, + const unsigned int zstep, + const unsigned int n_subdivisions) + { + if (patch->points_are_available) + { + unsigned int point_no=0; + // note: switch without break ! + switch (dim) + { + case 3: + Assert (zstepdata(patch->data.size(0)-spacedim+d,point_no); + } + else + { + // perform a dim-linear interpolation + const double stepsize=1./n_subdivisions, + xfrac=xstep*stepsize; + node = (patch->vertices[1] * xfrac) + (patch->vertices[0] * (1-xfrac)); + if (dim>1) + { + const double yfrac=ystep*stepsize; + node*= 1-yfrac; + node += ((patch->vertices[3] * xfrac) + (patch->vertices[2] * (1-xfrac))) * yfrac; + if (dim>2) + { + const double zfrac=zstep*stepsize; + node *= (1-zfrac); + node += (((patch->vertices[5] * xfrac) + (patch->vertices[4] * (1-xfrac))) + * (1-yfrac) + + ((patch->vertices[7] * xfrac) + (patch->vertices[6] * (1-xfrac))) + * yfrac) * zfrac; + } + } + } + } -//TODO: Should this go to the utilities namespace? - -// Compute n^dim, where dim is a template parameter -template -inline -unsigned int template_power (const unsigned int n) -{ - Assert (dim>0, ExcNotImplemented()); - switch (dim) - { - case 1: - return n; - case 2: - return n*n; - case 3: - return n*n*n; - case 4: - return n*n*n*n; - default: - unsigned int result = n; - for (unsigned int d=1;d -static -void -compute_sizes(const std::vector >& patches, - unsigned int& n_nodes, - unsigned int& n_cells) -{ - n_nodes = 0; - n_cells = 0; - for (typename std::vector >::const_iterator patch=patches.begin(); - patch!=patches.end(); ++patch) - { - n_nodes += template_power(patch->n_subdivisions+1); - n_cells += template_power(patch->n_subdivisions); - } + template + static + void + compute_sizes(const std::vector >& patches, + unsigned int& n_nodes, + unsigned int& n_cells) + { + n_nodes = 0; + n_cells = 0; + for (typename std::vector >::const_iterator patch=patches.begin(); + patch!=patches.end(); ++patch) + { + n_nodes += Utilities::fixed_power(patch->n_subdivisions+1); + n_cells += Utilities::fixed_power(patch->n_subdivisions); + } + } + } @@ -1305,7 +1284,7 @@ void DataOutBase::write_cells( } // finally update the number // of the first vertex of this patch - first_vertex_of_patch += template_power(n_subdivisions+1); + first_vertex_of_patch += Utilities::fixed_power(n_subdivisions+1); } } @@ -1331,7 +1310,7 @@ DataOutBase::write_data ( Assert ((patch->data.n_rows() == n_data_sets && !patch->points_are_available) || (patch->data.n_rows() == n_data_sets+spacedim && patch->points_are_available), ExcDimensionMismatch (patch->points_are_available ? (patch->data.n_rows() + spacedim) : patch->data.n_rows(), n_data_sets)); - Assert (patch->data.n_cols() == template_power(n), + Assert (patch->data.n_cols() == Utilities::fixed_power(n), ExcInvalidDatasetSize (patch->data.n_cols(), n)); std::vector floats(n_data_sets); @@ -1339,7 +1318,7 @@ DataOutBase::write_data ( // Data is already in // lexicographic ordering - for (unsigned int i=0; i(n); ++i) + for (unsigned int i=0; i(n); ++i) if (double_precision) { for (unsigned int data_set=0; data_set > &patches, const unsigned int n1 = (dim>0) ? n : 1; const unsigned int n2 = (dim>1) ? n : 1; const unsigned int n3 = (dim>2) ? n : 1; - unsigned int cells_per_patch = template_power(n); + unsigned int cells_per_patch = Utilities::fixed_power(n); unsigned int dx = 1; unsigned int dy = n; unsigned int dz = n*n; @@ -1775,7 +1754,7 @@ void DataOutBase::write_gnuplot (const std::vector > &patche Assert ((patch->data.n_rows() == n_data_sets && !patch->points_are_available) || (patch->data.n_rows() == n_data_sets+spacedim && patch->points_are_available), ExcDimensionMismatch (patch->points_are_available ? (patch->data.n_rows() + spacedim) : patch->data.n_rows(), n_data_sets)); - Assert (patch->data.n_cols() == template_power(n), + Assert (patch->data.n_cols() == Utilities::fixed_power(n), ExcInvalidDatasetSize (patch->data.n_cols(), n_subdivisions+1)); Point this_point; @@ -1993,7 +1972,7 @@ void DataOutBase::write_povray (const std::vector > &patches Assert ((patch->data.n_rows() == n_data_sets && !patch->points_are_available) || (patch->data.n_rows() == n_data_sets+spacedim && patch->points_are_available), ExcDimensionMismatch (patch->points_are_available ? (patch->data.n_rows() + spacedim) : patch->data.n_rows(), n_data_sets)); - Assert (patch->data.n_cols() == template_power(n_subdivisions+1), + Assert (patch->data.n_cols() == Utilities::fixed_power(n_subdivisions+1), ExcInvalidDatasetSize (patch->data.n_cols(), n_subdivisions+1)); for (unsigned int i=0; i > &patches Assert ((patch->data.n_rows() == n_data_sets && !patch->points_are_available) || (patch->data.n_rows() == n_data_sets+spacedim && patch->points_are_available), ExcDimensionMismatch (patch->points_are_available ? (patch->data.n_rows() + spacedim) : patch->data.n_rows(), n_data_sets)); - Assert (patch->data.n_cols() == template_power(n), + Assert (patch->data.n_cols() == Utilities::fixed_power(n), ExcInvalidDatasetSize (patch->data.n_cols(), n_subdivisions+1)); @@ -3240,7 +3219,7 @@ void DataOutBase::write_tecplot_binary (const std::vector > // finally update the number // of the first vertex of this patch - first_vertex_of_patch += template_power(n); + first_vertex_of_patch += Utilities::fixed_power(n); } @@ -3522,7 +3501,7 @@ DataOutBase::write_gmv_reorder_data_vectors (const std::vectordata.n_rows() == n_data_sets && !patch->points_are_available) || (patch->data.n_rows() == n_data_sets+spacedim && patch->points_are_available), ExcDimensionMismatch (patch->points_are_available ? (patch->data.n_rows() + spacedim) : patch->data.n_rows(), n_data_sets)); - Assert (patch->data.n_cols() == template_power(n_subdivisions+1), + Assert (patch->data.n_cols() == Utilities::fixed_power(n_subdivisions+1), ExcInvalidDatasetSize (patch->data.n_cols(), n_subdivisions+1)); for (unsigned int i=0;idata.n_cols();++i, ++next_value) diff --git a/deal.II/doc/news/changes.html b/deal.II/doc/news/changes.html index cc516410f3..48d5cf736e 100644 --- a/deal.II/doc/news/changes.html +++ b/deal.II/doc/news/changes.html @@ -437,6 +437,15 @@ inconvenience this causes.

base

    +
  1. New: The function Utilities::fixed_power<n>(q) + calculates q to the power of n where + n is a constant known at compile time. It allows to + calculate powers efficiently at compile time, most often things like a + number to the power dim. +
    + (WB 2007/06/23) +

    +
  2. New: The deal.II intermediate format has been changed. Since files in this format are only meant to be processed by exactly the same version of deal.II, this should