[Dune-devel] Destiny of MultiIndex
Oliver Sander
sander at igpm.rwth-aachen.de
Tue Mar 17 16:45:53 CET 2015
Hi Martin,
interesting how different people view this as different things.
> personally, I never understood the class MultiIndex in the StructuredGridFactory. To me, a multiindex is an array of integers with "vector space operations" (have a look at the MultiIndex in
> dune-spgrid, if you like). Unlike in StructuredGridFactory, it is not an iterator of any kind, i.e., ++multiIndex does not make sense to me.
I have always looked at multi-indices as dyadic representations of numbers,
where the base of each digit can be chosen separately. Then, ++multiIndex
makes perfect sense -- you simply increment the numbers.
>
> What you are looking for (as far as I underdand it) is a special kind of a set of multiindices: Those within within a "cube". This could be realized as an iterator range of multiindices. But for
> multiindices, other ranges are of interest too (e.g., all multiindices within a "simplex").
I've never thought about the simplex problem; you are right that this can be of interest, too.
Without having giving it much thought I tend to think that your approach (i.e., having algorithmically
defined sets of multi-indices, and iterators over these sets) is more general, and hence preferable.
Is your implementation easy to look at and understand?
Cheers,
Oliver
>
> Apart from such implementation details, I'd be very happy to have such a concept in dune-common, as it would reduce the SPGrid code by a few hundred (basically redundant) lines of code.
>
> Best,
>
> Martin
>
>
> On 03/16/2015 10:26 PM, Oliver Sander wrote:
>> Am 16.03.2015 um 19:33 schrieb Dominic Kempf:
>>> Sounds like a good idea to me. Maybe we could/should extend the interface
>>> then? postfix operator++, the operator--s, logical and stream operators
>>> come to my mind.
>>
>> Definitely. Also, one usually needs a way to iterate over the entire range
>> of the multi-index. Previously I have used a method 'cycle', which says how
>> many times I can call operator++ before I getting an overflow, but there may
>> be better solutions.
>>
>>>
>>> On Mon, Mar 16, 2015 at 7:17 PM, Oliver Sander <sander at igpm.rwth-aachen.de>
>>> wrote:
>>>
>>>> Hi all,
>>>>
>>>>>
>>>>> Move MultiIndex class from StructuredGridFactory to a separate header
>>>>>
>>>>> The TensorGridFactory needs the same class, so I move it to a new
>>>> header.
>>>>
>>>> a smart move. I have needed (and reimplemented) such a class again and
>>>> again.
>>>> Maybe we should move multiindex.hh into dune-common right away?
>>>> --
>>>> Oliver
>>>>
>>>>
>>>>>
>>>>> dune/grid/utility/CMakeLists.txt | 1 +
>>>>> dune/grid/utility/Makefile.am | 1 +
>>>>> dune/grid/utility/multiindex.hh | 57
>>>> ++++++++++++++++++++++++++++++
>>>>> dune/grid/utility/structuredgridfactory.hh | 53
>>>> +++------------------------
>>>>> dune/grid/utility/tensorgridfactory.hh | 2 ++
>>>>> 5 files changed, 65 insertions(+), 49 deletions(-)
>>>>> create mode 100644 dune/grid/utility/multiindex.hh
>>>>>
>>>>>
>>>>>
>>>>> diff --git a/dune/grid/utility/CMakeLists.txt
>>>> b/dune/grid/utility/CMakeLists.txt
>>>>> index 1c42604..15d2311 100644
>>>>> --- a/dune/grid/utility/CMakeLists.txt
>>>>> +++ b/dune/grid/utility/CMakeLists.txt
>>>>> @@ -8,6 +8,7 @@ set(HEADERS
>>>>> gridtype.hh
>>>>> hierarchicsearch.hh
>>>>> hostgridaccess.hh
>>>>> + multiindex.hh
>>>>> parmetisgridpartitioner.hh
>>>>> persistentcontainer.hh
>>>>> persistentcontainerinterface.hh
>>>>> diff --git a/dune/grid/utility/Makefile.am
>>>> b/dune/grid/utility/Makefile.am
>>>>> index e14c152..d81c3e2 100644
>>>>> --- a/dune/grid/utility/Makefile.am
>>>>> +++ b/dune/grid/utility/Makefile.am
>>>>> @@ -10,6 +10,7 @@ gridutility_HEADERS = \
>>>>> gridtype.hh \
>>>>> hierarchicsearch.hh \
>>>>> hostgridaccess.hh \
>>>>> + multiindex.hh \
>>>>> parmetisgridpartitioner.hh \
>>>>> persistentcontainer.hh \
>>>>> persistentcontainerinterface.hh \
>>>>> diff --git a/dune/grid/utility/multiindex.hh
>>>> b/dune/grid/utility/multiindex.hh
>>>>> new file mode 100644
>>>>> index 0000000..b6b281f
>>>>> --- /dev/null
>>>>> +++ b/dune/grid/utility/multiindex.hh
>>>>> @@ -0,0 +1,57 @@
>>>>> +#ifndef DUNE_GRID_UTILITY_MULTIINDEX_HH
>>>>> +#define DUNE_GRID_UTILITY_MULTIINDEX_HH
>>>>> +
>>>>> +/** \file
>>>>> + * \brief Implements a multiindex with arbitrary dimension and fixed
>>>> index ranges
>>>>> + * This is used by various factory classes.
>>>>> + */
>>>>> +
>>>>> +#include<array>
>>>>> +
>>>>> +namespace Dune
>>>>> +{
>>>>> + namespace FactoryUtilities
>>>>> + {
>>>>> + template<std::size_t dim>
>>>>> + class MultiIndex : public std::array<unsigned int,dim>
>>>>> + {
>>>>> + // The range of each component
>>>>> + std::array<unsigned int,dim> limits_;
>>>>> +
>>>>> + public:
>>>>> + /** \brief Constructor with a given range for each digit */
>>>>> + MultiIndex(const std::array<unsigned int,dim>& limits) :
>>>> limits_(limits)
>>>>> + {
>>>>> + std::fill(this->begin(), this->end(), 0);
>>>>> + }
>>>>> +
>>>>> + /** \brief Increment the MultiIndex */
>>>>> + MultiIndex<dim>& operator++()
>>>>> + {
>>>>> + for (int i=0; i<dim; i++)
>>>>> + {
>>>>> + // Augment digit
>>>>> + (*this)[i]++;
>>>>> +
>>>>> + // If there is no carry-over we can stop here
>>>>> + if ((*this)[i]<limits_[i])
>>>>> + break;
>>>>> +
>>>>> + (*this)[i] = 0;
>>>>> + }
>>>>> + return *this;
>>>>> + }
>>>>> +
>>>>> + /** \brief Compute how many times you can call operator++ before
>>>> getting to (0,...,0) again */
>>>>> + size_t cycle() const
>>>>> + {
>>>>> + size_t result = 1;
>>>>> + for (int i=0; i<dim; i++)
>>>>> + result *= limits_[i];
>>>>> + return result;
>>>>> + }
>>>>> + };
>>>>> + }
>>>>> +}
>>>>> +
>>>>> +#endif
>>>>> diff --git a/dune/grid/utility/structuredgridfactory.hh
>>>> b/dune/grid/utility/structuredgridfactory.hh
>>>>> index 8832149..8e2700a 100644
>>>>> --- a/dune/grid/utility/structuredgridfactory.hh
>>>>> +++ b/dune/grid/utility/structuredgridfactory.hh
>>>>> @@ -19,6 +19,7 @@
>>>>> #include <dune/common/shared_ptr.hh>
>>>>>
>>>>> #include <dune/grid/common/gridfactory.hh>
>>>>> +#include <dune/grid/utility/multiindex.hh>
>>>>> #include <dune/grid/yaspgrid.hh>
>>>>>
>>>>> namespace Dune {
>>>>> @@ -34,59 +35,13 @@ namespace Dune {
>>>>>
>>>>> static const int dimworld = GridType::dimensionworld;
>>>>>
>>>>> - /** \brief dim-dimensional multi-index. The range for each
>>>> component can be set individually
>>>>> - */
>>>>> - class MultiIndex
>>>>> - : public array<unsigned int,dim>
>>>>> - {
>>>>> -
>>>>> - // The range of each component
>>>>> - array<unsigned int,dim> limits_;
>>>>> -
>>>>> - public:
>>>>> - /** \brief Constructor with a given range for each digit */
>>>>> - MultiIndex(const array<unsigned int,dim>& limits)
>>>>> - : limits_(limits)
>>>>> - {
>>>>> - std::fill(this->begin(), this->end(), 0);
>>>>> - }
>>>>> -
>>>>> - /** \brief Increment the MultiIndex */
>>>>> - MultiIndex& operator++() {
>>>>> -
>>>>> - for (int i=0; i<dim; i++) {
>>>>> -
>>>>> - // Augment digit
>>>>> - (*this)[i]++;
>>>>> -
>>>>> - // If there is no carry-over we can stop here
>>>>> - if ((*this)[i]<limits_[i])
>>>>> - break;
>>>>> -
>>>>> - (*this)[i] = 0;
>>>>> -
>>>>> - }
>>>>> - return *this;
>>>>> - }
>>>>> -
>>>>> - /** \brief Compute how many times you can call operator++ before
>>>> getting to (0,...,0) again */
>>>>> - size_t cycle() const {
>>>>> - size_t result = 1;
>>>>> - for (int i=0; i<dim; i++)
>>>>> - result *= limits_[i];
>>>>> - return result;
>>>>> - }
>>>>> -
>>>>> - };
>>>>> -
>>>>> /** \brief Insert a structured set of vertices into the factory */
>>>>> static void insertVertices(GridFactory<GridType>& factory,
>>>>> const FieldVector<ctype,dimworld>&
>>>> lowerLeft,
>>>>> const FieldVector<ctype,dimworld>&
>>>> upperRight,
>>>>> const array<unsigned int,dim>& vertices)
>>>>> {
>>>>> -
>>>>> - MultiIndex index(vertices);
>>>>> + FactoryUtilities::MultiIndex<dim> index(vertices);
>>>>>
>>>>> // Compute the total number of vertices to be created
>>>>> int numVertices = index.cycle();
>>>>> @@ -166,7 +121,7 @@ namespace Dune {
>>>>> cornersTemplate[i] += unitOffsets[j];
>>>>>
>>>>> // Insert elements
>>>>> - MultiIndex index(elements);
>>>>> + FactoryUtilities::MultiIndex<dim> index(elements);
>>>>>
>>>>> // Compute the total number of elementss to be created
>>>>> int numElements = index.cycle();
>>>>> @@ -234,7 +189,7 @@ namespace Dune {
>>>>>
>>>>> // Loop over all "cubes", and split up each cube into dim!
>>>>> // (factorial) simplices
>>>>> - MultiIndex elementsIndex(elements);
>>>>> + FactoryUtilities::MultiIndex<dim> elementsIndex(elements);
>>>>> size_t cycle = elementsIndex.cycle();
>>>>>
>>>>> for (size_t i=0; i<cycle; ++elementsIndex, i++) {
>>>>> diff --git a/dune/grid/utility/tensorgridfactory.hh
>>>> b/dune/grid/utility/tensorgridfactory.hh
>>>>> index 971bc4e..8aa47cf 100644
>>>>> --- a/dune/grid/utility/tensorgridfactory.hh
>>>>> +++ b/dune/grid/utility/tensorgridfactory.hh
>>>>> @@ -19,6 +19,8 @@
>>>>> #include<memory>
>>>>> #include<vector>
>>>>>
>>>>> +#include<dune/grid/utility/multiindex.hh>
>>>>> +
>>>>> namespace Dune
>>>>> {
>>>>> // forward declaration of TensorGridFactoryCreator, which is the real
>>>> factory
>>>>>
>>>>> _______________________________________________
>>>>> Dune-Commit mailing list
>>>>> Dune-Commit at dune-project.org
>>>>> http://lists.dune-project.org/mailman/listinfo/dune-commit
>>>>>
>>>>
>>>>
>>>>
>>>> _______________________________________________
>>>> Dune-devel mailing list
>>>> Dune-devel at dune-project.org
>>>> http://lists.dune-project.org/mailman/listinfo/dune-devel
>>>>
>>>>
>>>
>>
>>
>>
>>
>> _______________________________________________
>> Dune-devel mailing list
>> Dune-devel at dune-project.org
>> http://lists.dune-project.org/mailman/listinfo/dune-devel
>>
>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: signature.asc
Type: application/pgp-signature
Size: 473 bytes
Desc: OpenPGP digital signature
URL: <https://lists.dune-project.org/pipermail/dune-devel/attachments/20150317/ac6a6ac2/attachment.sig>
More information about the Dune-devel
mailing list