[dune-pdelab] [dune-pdelab-commit] dune-pdelab r471 - trunk/dune/pdelab/gridfunctionspace

Bernd Flemisch bernd at iws.uni-stuttgart.de
Mon May 17 14:27:22 CEST 2010


Hey Christian,

thank you for the explanation. I try build a minimal example.

Kind regards
Bernd

On 05/17/2010 02:13 PM, Christian Engwer wrote:
> Up to recently the you had two options regarding the blockin of you
> dofs. Assume two cells A,B with two physical values x,y and two DOFs
> 1,2.
>
> you can have Elementwise blocking of physical value wise
> blocking. This results two different numberings:
>
> Elementwise : [Ax1, Ay1, Ax2, Ay2,   Bx1, By1, By2, By2]
> Phys. wise  : [Ax1, Ax2, Bx1, Bx2,   Ay1, Ay2, By1, By2]
>
> For elementwise blocking it was necessary that all physical values
> have the same number of dofs. This rendered certain setups
> impossible. Think of Stokes equation. 
>
> For a Stoeks DG discretitzation you have a function space like this:
>   (P_k)^dim x P_{k-1}
>
> In order to benefit from you DG structure you want to do element wise
> blocking and you want a structure in the following way:
>  [Av1_1, .. Av1_n, Av2_1, .. Av2_n, Av3_1, .. Av3_n, Ap_1, .. Ap_m, Bv1_1, ...]
>
> The changes in the code are intended to allow exactly this
> behaviour. I tried to keep them backwards compatible.
>
> You please describe which exact setup you have and how the indices
> changed. Please give a boiled down test case ;-)
>
> Christian
>
> On Mon, May 17, 2010 at 11:17:59AM +0200, Bernd Flemisch wrote:
>   
>> Dear Dune-PDELab, hey Christian,
>>
>> this commit unfortunately breaks my code in the sense that I obtain
>> nonsense solutions. I guess that some assumption I put in on the
>> indexing is not valid any more. Could you maybe explain your commit?
>>
>> Thank you. Kind regards
>> Bernd
>>
>> On 05/07/2010 08:05 PM, christi at conan.iwr.uni-heidelberg.de wrote:
>>     
>>> Author: christi
>>> Date: Fri May  7 20:05:07 2010
>>> New Revision: 471
>>> URL: http://svn.dune-project.org/websvn/listing.php?repname=dune-pdelab&path=/&rev=471&sc=1
>>>
>>> Log:
>>> allow to change the blocking structure
>>>
>>> Modified:
>>>    trunk/dune/pdelab/gridfunctionspace/gridfunctionspace.hh
>>>
>>> Modified: trunk/dune/pdelab/gridfunctionspace/gridfunctionspace.hh
>>> ==============================================================================
>>> --- trunk/dune/pdelab/gridfunctionspace/gridfunctionspace.hh	Fri May  7 11:34:26 2010	(r470)
>>> +++ trunk/dune/pdelab/gridfunctionspace/gridfunctionspace.hh	Fri May  7 20:05:07 2010	(r471)
>>> @@ -4,22 +4,22 @@
>>>  #define DUNE_PDELAB_GRIDFUNCTIONSPACE_HH
>>>  
>>>  #include <cstddef>
>>> -#include<map>
>>> +#include <map>
>>>  #include <ostream>
>>> -#include<set>
>>> -#include<vector>
>>> +#include <set>
>>> +#include <vector>
>>>  
>>> -#include<dune/common/exceptions.hh>
>>> -#include<dune/common/geometrytype.hh>
>>> -#include<dune/common/stdstreams.hh>
>>> -#include<dune/grid/common/genericreferenceelements.hh>
>>> +#include <dune/common/exceptions.hh>
>>> +#include <dune/common/geometrytype.hh>
>>> +#include <dune/common/stdstreams.hh>
>>> +#include <dune/grid/common/genericreferenceelements.hh>
>>>  
>>>  #include <dune/localfunctions/common/localkey.hh>
>>>  
>>> -#include"../common/countingptr.hh"
>>> -#include"../common/multitypetree.hh"
>>> -#include"../common/cpstoragepolicy.hh"
>>> -#include"../common/geometrywrapper.hh"
>>> +#include "../common/countingptr.hh"
>>> +#include "../common/multitypetree.hh"
>>> +#include "../common/cpstoragepolicy.hh"
>>> +#include "../common/geometrywrapper.hh"
>>>  
>>>  #include"localfunctionspace.hh"
>>>  
>>> @@ -175,12 +175,8 @@
>>>  	/** \brief Tag indicating an arbitrary number of unkowns per entity.
>>>       *
>>>       * class used to pass compile-time parameter to the GridFunctionSpace.
>>> -     *
>>>       */
>>> -	struct GridFunctionGeneralMapper
>>> -	{
>>> -	  enum {dummy=0} ;
>>> -	};
>>> +	struct GridFunctionGeneralMapper {};
>>>  
>>>  
>>>      // Empty constraints assembler class
>>> @@ -494,17 +490,12 @@
>>>        std::set<unsigned int> codimUsed;
>>>  	};
>>>  
>>> -
>>> -
>>>  	/** \brief Tag indicating a fixed number of unkowns per entity (known at compile time).
>>>       *
>>>       * class used to pass compile-time parameter to the GridFunctionSpace.
>>>       *
>>>       */
>>> -	struct GridFunctionRestrictedMapper
>>> -	{
>>> -	  enum {dummy=1} ;
>>> -	};
>>> +	struct GridFunctionRestrictedMapper {};
>>>  
>>>      //! \}
>>>  
>>> @@ -782,7 +773,6 @@
>>>      template<typename IIS>
>>>  	struct GridFunctionStaticSize
>>>  	{
>>> -	  enum {dummy=2} ;
>>>        typedef IIS IntersectionIndexSet;
>>>  	};
>>>  
>>> @@ -1114,35 +1104,57 @@
>>>      //! \brief Indicates lexicographics ordering of the unknowns for composite
>>>      //! grid function spaces.
>>>      //!
>>> -	//! this class may be used to pass compile-time
>>> -	//! parameters to the implementation of 
>>> +    //! this class may be used to pass compile-time
>>> +    //! parameters to the implementation of 
>>>      //! \link PowerGridFunctionSpace PowerGridFunctionSpace \endlink or
>>>      //! \link CompositeGridFunctionSpace CompositeGridFunctionSpace \endlink
>>> -	struct GridFunctionSpaceLexicographicMapper
>>> -	{
>>> -	  enum {dummy=0} ;
>>> -	};
>>> +    struct GridFunctionSpaceLexicographicMapper {};
>>>  
>>>      //! \brief Indicates using block-wise ordering of the unknowns for composite
>>>      //! grid function spaces.
>>>      //!
>>> -	//! this class may be used to pass compile-time
>>> -	//! parameters to the implementation of 
>>> +    //! The exact blocking structure can be passed as template parameters
>>> +    //!
>>> +    //! this class may be used to pass compile-time
>>> +    //! parameters to the implementation of 
>>>      //! \link PowerGridFunctionSpace PowerGridFunctionSpace \endlink or
>>>      //! \link CompositeGridFunctionSpace CompositeGridFunctionSpace \endlink
>>> -	struct GridFunctionSpaceBlockwiseMapper
>>> -	{
>>> -	  enum {dummy=0} ;
>>> -	};
>>> +    template<int s0 = 1, int s1 = 1, int s2 = 1, int s3 = 1, int s4 = 1, int s5 = 1, int s6 = 1, int s7 = 1, int s8 = 1, int s9 = 1>
>>> +    struct GridFunctionSpaceComponentBlockwiseMapper
>>> +    {
>>> +      static const int size[];
>>> +      static const int offset[];
>>> +    };
>>> +    template<int s0, int s1, int s2, int s3, int s4, int s5, int s6, int s7, int s8, int s9>
>>> +    const int GridFunctionSpaceComponentBlockwiseMapper<s0,s1,s2,s3,s4,s5,s6,s7,s8,s9>::
>>> +      size[] = { s0, s1, s2, s3, s4, s5, s6, s7, s8, s9 };
>>> +    template<int s0, int s1, int s2, int s3, int s4, int s5, int s6, int s7, int s8, int s9>
>>> +    const int GridFunctionSpaceComponentBlockwiseMapper<s0,s1,s2,s3,s4,s5,s6,s7,s8,s9>::
>>> +      offset[] = { 0, s0, s0+s1, s0+s1+s2, s0+s1+s2+s3, s0+s1+s2+s3+s4,
>>> +                   s0+s1+s2+s3+s4+s5, s0+s1+s2+s3+s4+s5+s6, s0+s1+s2+s3+s4+s5+s6+s7,
>>> +                   s0+s1+s2+s3+s4+s5+s6+s7+s8, s0+s1+s2+s3+s4+s5+s6+s7+s8+s9 };
>>> +    
>>> +    //! \brief Indicates using block-wise ordering of the unknowns for composite
>>> +    //! grid function spaces.
>>> +    //!
>>> +    //! this class may be used to pass compile-time
>>> +    //! parameters to the implementation of 
>>> +    //! \link PowerGridFunctionSpace PowerGridFunctionSpace \endlink or
>>> +    //! \link CompositeGridFunctionSpace CompositeGridFunctionSpace \endlink
>>> +    struct GridFunctionSpaceBlockwiseMapper : GridFunctionSpaceComponentBlockwiseMapper<> {};
>>>  
>>>      //! \}
>>>  
>>>  	template<typename T, int k, typename P>
>>>  	class PowerGridFunctionSpace;
>>>  
>>> -    // product of identical grid function spaces
>>> -    // base class that holds implementation of the methods
>>> -    // this is the default version with lexicographic ordering
>>> +    //! product of identical grid function spaces
>>> +    //! base class that holds implementation of the methods
>>> +    //! this is the default version with lexicographic ordering
>>> +    //!
>>> +    //! \tparam T PLEASE DOCUMENT
>>> +    //! \tparam k PLEASE DOCUMENT
>>> +    //! \tparam P PLEASE DOCUMENT
>>>  	template<typename T, int k, typename P>
>>>  	class PowerGridFunctionSpaceBase 
>>>        : public PowerNode<T,k,CountingPointerStoragePolicy>,
>>> @@ -1297,7 +1309,7 @@
>>>      private:
>>>        void setup ()
>>>        {
>>> -        Dune::dinfo << "power grid function space(lexicographic version):"
>>> +        Dune::dinfo << "PowerGridFunctionSpace(lexicographic version):"
>>>                      << std::endl;
>>>          Dune::dinfo << "( ";
>>>          offset[0] = 0;
>>> @@ -1325,18 +1337,17 @@
>>>      // product of identical grid function spaces
>>>      // base class that holds implementation of the methods
>>>      // specialization for blockwise ordering
>>> -	template<typename T, int k>
>>> -	class PowerGridFunctionSpaceBase<T,k,GridFunctionSpaceBlockwiseMapper> 
>>> +	template<typename T, int k, int s>
>>> +	class PowerGridFunctionSpaceBase<T,k,GridFunctionSpaceComponentBlockwiseMapper<s> >
>>>        : public PowerNode<T,k,CountingPointerStoragePolicy>,
>>>          public Countable
>>>  	{
>>> -      friend class PowerGridFunctionSpace<T,k,GridFunctionSpaceBlockwiseMapper>;
>>> -
>>> +      friend class PowerGridFunctionSpace<T,k,GridFunctionSpaceComponentBlockwiseMapper<s> >;
>>>  	public:
>>>        //! export traits class
>>>        typedef PowerCompositeGridFunctionSpaceTraits<typename T::Traits::GridViewType, 
>>>                                                      typename T::Traits::BackendType,
>>> -                                                    GridFunctionSpaceBlockwiseMapper, k>
>>> +                                                    GridFunctionSpaceComponentBlockwiseMapper<s>, k>
>>>        Traits;
>>>  
>>>        //! extract type of container storing Es
>>> @@ -1408,7 +1419,7 @@
>>>        template<int i>
>>>  	  typename Traits::SizeType subMap (typename Traits::SizeType j) const
>>>  	  {
>>> -		return j*k+i;
>>> +        return (j%s)+(j/s)*k*s+i*s;
>>>  	  }
>>>  
>>>        //------------------------------
>>> @@ -1478,21 +1489,26 @@
>>>      private:
>>>        void setup ()
>>>        {
>>> -        Dune::dinfo << "power grid function space(blockwise version):"
>>> +        Dune::dinfo << "PowerGridFunctionSpace(blockwise version):"
>>>                      << std::endl;
>>>          Dune::dinfo << "( ";
>>>          offset[0] = 0;
>>>          maxlocalsize = 0;
>>>          for (int i=0; i<k; i++)
>>> -          {
>>> -            childSize[i] = this->getChild(i).globalSize();
>>> -            Dune::dinfo << childSize[i] << " ";
>>> -            offset[i+1] = offset[i]+childSize[i];
>>> -            maxlocalsize += this->getChild(i).maxLocalSize();
>>> -          }
>>> +        {
>>> +          childSize[i] = this->getChild(i).globalSize();
>>> +          offset[i+1] = offset[i]+childSize[i];
>>> +          Dune::dinfo << childSize[i] << "[" << offset[i] << "] ";
>>> +          maxlocalsize += this->getChild(i).maxLocalSize();
>>> +        }
>>>          Dune::dinfo << ") total size = " << offset[k]
>>>                      << " max local size = " << maxlocalsize
>>>                      << std::endl;
>>> +        /* check the local block size */
>>> +        if (this->getChild(0).maxLocalSize()%s != 0)
>>> +          DUNE_THROW(Exception,
>>> +            "number of DOFs (" << this->getChild(0).maxLocalSize() << ") per component "
>>> +            "must be a multiple of the BlockSize (" << s << ")");
>>>          for (int i=1; i<k; i++)
>>>            if (childSize[i]!=childSize[0])
>>>              DUNE_THROW(Exception, "components must be of equal size");
>>> @@ -2181,7 +2197,7 @@
>>>      private:
>>>        void setup ()
>>>        {
>>> -        Dune::dinfo << "composite grid function space(lexicographic version):"
>>> +        Dune::dinfo << "CompositeGridFunctionSpace(lexicographic version):"
>>>                      << std::endl;
>>>  
>>>          CompositeGridFunctionSpaceBaseVisitChildMetaProgram<CompositeGridFunctionSpaceBase,BaseT::CHILDREN,0>::
>>> @@ -2216,14 +2232,15 @@
>>>      // P is the ordering parameter
>>>      // Ti are all grid function spaces
>>>  	template<typename T0, typename T1, typename T2, typename T3,
>>> -			 typename T4, typename T5, typename T6,
>>> -			 typename T7, typename T8>
>>> -	class CompositeGridFunctionSpaceBase<GridFunctionSpaceBlockwiseMapper,T0,T1,T2,T3,T4,T5,T6,T7,T8>
>>> +			 typename T4, typename T5, typename T6, typename T7, typename T8,
>>> +             int s0, int s1, int s2, int s3, int s4, int s5, int s6, int s7, int s8, int s9>
>>> +	class CompositeGridFunctionSpaceBase<GridFunctionSpaceComponentBlockwiseMapper<s0,s1,s2,s3,s4,s5,s6,s7,s8,s9>,
>>> +                                         T0,T1,T2,T3,T4,T5,T6,T7,T8>
>>>  	  : public CompositeNode<CountingPointerStoragePolicy,T0,T1,T2,T3,T4,T5,T6,T7,T8>,
>>>  		public Countable
>>>  	{
>>> -      friend class CompositeGridFunctionSpace<GridFunctionSpaceBlockwiseMapper,
>>> -                                              T0,T1,T2,T3,T4,T5,T6,T7,T8>; // for setup
>>> +      typedef GridFunctionSpaceComponentBlockwiseMapper<s0,s1,s2,s3,s4,s5,s6,s7,s8,s9> BlockwiseMapper;
>>> +      friend class CompositeGridFunctionSpace<BlockwiseMapper, T0,T1,T2,T3,T4,T5,T6,T7,T8>; // for setup
>>>  
>>>        typedef CompositeNode<CountingPointerStoragePolicy,T0,T1,T2,T3,T4,T5,T6,T7,T8> BaseT;
>>>  
>>> @@ -2231,7 +2248,7 @@
>>>        //! export traits class
>>>        typedef PowerCompositeGridFunctionSpaceTraits<typename T0::Traits::GridViewType, 
>>>                                                      typename T0::Traits::BackendType,
>>> -                                                    GridFunctionSpaceBlockwiseMapper,
>>> +                                                    BlockwiseMapper,
>>>                                                      NonEmptyChilds<T0,T1,T2,T3,T4,T5,
>>>                                                                     T6,T7,T8>::value>
>>>        Traits;
>>> @@ -2301,7 +2318,15 @@
>>>        template<int i>
>>>  	  typename Traits::SizeType subMap (typename Traits::SizeType j) const
>>>  	  {
>>> -		return j*BaseT::CHILDREN+i;
>>> +        // make the block sizes and offsets available in an array
>>> +        static const int blockSize[] = { s0, s1, s2, s3, s4, s5, s6, s7, s8, s9 };
>>> +        static const int blockOffset[] = { 0, s0, s0+s1, s0+s1+s2, s0+s1+s2+s3, s0+s1+s2+s3+s4,
>>> +                                           s0+s1+s2+s3+s4+s5, s0+s1+s2+s3+s4+s5+s6, s0+s1+s2+s3+s4+s5+s6+s7,
>>> +                                           s0+s1+s2+s3+s4+s5+s6+s7+s8, s0+s1+s2+s3+s4+s5+s6+s7+s8+s9 };
>>> +        return (j%BlockwiseMapper::size[i])
>>> +          +(j/BlockwiseMapper::size[i])*BlockwiseMapper::offset[BaseT::CHILDREN]
>>> +          +BlockwiseMapper::offset[i];
>>> +        return (j%blockSize[i])+(j/blockSize[i])*blockOffset[BaseT::CHILDREN]+blockOffset[i];
>>>  	  }
>>>  
>>>        //------------------------------
>>> @@ -2357,15 +2382,24 @@
>>>      private:
>>>        void setup ()
>>>        {
>>> -        Dune::dinfo << "composite grid function space(blockwise version):"
>>> +        Dune::dinfo << "CompositeGridFunctionSpace(blockwise version):"
>>>                      << std::endl;
>>>  
>>>          CompositeGridFunctionSpaceBaseVisitChildMetaProgram<CompositeGridFunctionSpaceBase,BaseT::CHILDREN,0>::
>>>            setup(*this,childGlobalSize,childLocalSize);
>>>  
>>> +        // make the block sizes available in an array
>>> +        static const int blockSize[] = { s0, s1, s2, s3, s4, s5, s6, s7, s8, s9 };
>>> +        // check for compatible sizes
>>>          for (int i=1; i<BaseT::CHILDREN; i++)
>>> -          if (childGlobalSize[i]!=childGlobalSize[0])
>>> +        {
>>> +          if (childLocalSize[i]%blockSize[i]!=0)
>>> +            DUNE_THROW(Exception,
>>> +              "number of DOFs (" << childLocalSize[i] << ") per component "
>>> +              "must be a multiple of the BlockSize (" << blockSize[i] << ")");
>>> +          if (childGlobalSize[i]/blockSize[i]!=childGlobalSize[0]/blockSize[0])
>>>              DUNE_THROW(Exception, "components must be of equal size");
>>> +        }
>>>  
>>>          Dune::dinfo << "( ";
>>>          offset[0] = 0;
>>> @@ -2389,7 +2423,6 @@
>>>        mutable std::vector<typename Traits::SizeType> childglobal;
>>>  	};
>>>  
>>> -
>>>      //! \addtogroup GridFunctionSpace \{
>>>  
>>>      //! \brief grid function space composed of other grid function spaces
>>> @@ -3199,7 +3232,7 @@
>>>        GridFunctionSubSpace (const GFS& gfs) 
>>>          : GridFunctionSubSpaceIntermediateBase<GFS,k,GFS::isLeaf>(gfs)
>>>        {
>>> -        Dune::dinfo << "grid function subspace:" << std::endl;
>>> +        Dune::dinfo << "GridFunctionSubSpace:" << std::endl;
>>>          Dune::dinfo << "root space size = " << gfs.globalSize()
>>>                      << " max local size = " << this->maxLocalSize()
>>>                      << std::endl;
>>>
>>> _______________________________________________
>>> dune-pdelab-commit mailing list
>>> dune-pdelab-commit at dune-project.org
>>> http://lists.dune-project.org/mailman/listinfo/dune-pdelab-commit
>>>
>>>   
>>>       
>>     
>   

-- 
_____________________________________________________________________

Bernd Flemisch                               phone: +49 711 685 69162
IWS, Universität Stuttgart                   fax:   +49 711 685 60430
Pfaffenwaldring 61                  email: bernd at iws.uni-stuttgart.de
D-70569 Stuttgart                  url: www.hydrosys.uni-stuttgart.de
_____________________________________________________________________





More information about the dune-pdelab mailing list