[dune-pdelab] periodic boundary conditions

[Oleh Krehel]

Thanks, Christian.
It works now. But what exactly is the overlap?
(Documentation says something like "overlap is overlap")

2011/11/2 Christian Engwer <christian.engwer at uni-muenster.de>:
> On Wed, Nov 02, 2011 at 01:12:17PM +0100, Oleh Krehel wrote:
>> Hi,
>>
>> I tried this code:
>>
>>   const int d=2;
>>   Dune::FieldVector<double,d> L(1);
>>   Dune::FieldVector<int,d> N(2);
>>   Dune::FieldVector<bool,d> periodic(true);
>>   int overlap=0;
>        ^^^^^^^^^
>
>>   Dune::YaspGrid<d> grid(L,N,periodic,overlap);
>>
>> Then when I traverse on the leafview the instersections of each element,
>> the sets of is->neighbor() and is->boundary() don't overlap, they're
>> just like the non-periodic ones.
>
> You asked for a grid without overlap
>
>> Can someone please explain, maybe there's something I've missed?
>>
>> regards,
>> Oleh
>>
>>
>>
>>
>>
>> 2011/10/31 Christian Engwer <christian.engwer at uni-muenster.de>:
>> > Hi!
>> >
>> > first you have to settle the question of the grid. I don't have a
>> > complete overview of the support of periodic boundaries of the
>> > different grid managers.
>> >
>> > I know:
>> > - YaspGrid handles what you want, but without inclusions
>> > - UG does have some kind of periodic support, but I don't know what
>> >  exactly is supported and you can't use easily fromm DUNE [patches
>> >  wellcome ;-)]
>> > - Alugrid from 1.23 on support periodic boundaries in 3D
>> >
>> > Depending on the grid you are using the periodicity is implemented
>> > similar to parallel coupling, which uses domain decomposition with
>> > overlap or ghosts. In this case you should use the periodic
>> > constraints. If the grid implements periodicity as periodicity in the
>> > topology of the grid, you don't have to change anything (not sure if
>> > any grid supports this).
>> >
>> > The last option is implement periodictity in PDELab without support of
>> > the grid. You should the write your own constraints-assembler
>> > local-operator, which adds dirichlet constraints on one side. There is
>> > no real tutorial on this. You will have to dive into it...
>> >
>> > Christian
>> >
>> >> I need a fully periodic cube: left-right, top-bottom in 2D,
>> >> and additionally front-back in 3D.
>> >> The cube has an arbitrarily shaped inclusion(like a circle or a sphere) inside,
>> >> with Dirichlet b.c.
>> >> I was looking at section 8.4 of the tutorial explaining constraints
>> >> for Dirichlet b.c.,
>> >> but the transformation matrix there is local, when I need the global one for the
>> >> periodic b.c.
>> >> Can you point me to the source file that does what you described above?
>> >>
>> >> regards,
>> >> Oleh
>> >>
>> >>
>> >>
>> >> 2011/10/31 Christian Engwer <christian.engwer at uni-muenster.de>:
>> >> > Dear Oleh,
>> >> >
>> >> >> I need to implement periodic boundary conditions for a transport problem.
>> >> >> Google shows that there's no such functionality in DUNE at the moment.
>> >> >> Perhaps someone can give me some pointers or suggestions on how
>> >> >> to implement periodic b.c. in a straightforward way.
>> >> >
>> >> > it depends on which BC you need.
>> >> >
>> >> > Straight forward: You have a mesh like this:
>> >> >
>> >> > o-x-x-x-o
>> >> > | | | | |
>> >> > o-x-x-x-o
>> >> > | | | | |
>> >> > o-x-x-x-o
>> >> >
>> >> > where the 'o' nodes are associated with each other. This kind of mesh
>> >> > is available with yaspgrid. It is implemented as a periodic mesh with
>> >> > overlap on the left and right side.
>> >> >
>> >> > If you need something more general it will get tricky and you should
>> >> > explain in detail what you want to achieve.
>> >> >
>> >> > Cheers
>> >> > Christian
>> >> >
>> >> >
>> >>
>> >
>
>> #define DUNE_DEVEL_MODE
>> //#define NDEBUG
>> #include "../config.h"
>> #include <iostream>
>> #include <dune/common/mpihelper.hh>
>> #include <dune/common/exceptions.hh>
>> #include <boost/exception/all.hpp>
>> #include <dune/grid/yaspgrid.hh>
>>
>> template <class GV>
>> void traverseGrid(const GV& gv){
>>   typedef typename GV::template Codim<0>::Iterator element_iterator;
>>   typedef typename GV::IntersectionIterator intersection_iterator;
>>   element_iterator it=gv.template begin<0>(),
>>                   end=gv.template end<0>();
>>   for(;it!=end;++it){
>>     Dune::GeometryType gt=it->type();
>>     std::cout<<"\nvisiting leaf "<<gt
>>              <<" with first vertex at "<<it->geometry().corner(0);
>>     intersection_iterator is=gv.ibegin(*it),
>>                        isend=gv.iend(*it);
>>     for(;is!=isend;++is){
>>       if(is->neighbor())
>>         std::cout<<"\nneighbor";
>>       if(is->boundary())
>>         std::cout<<"\nboundary";
>>     }
>>   }
>> }
>>
>> void dowork(int level){
>>   const int d=2;
>>   Dune::FieldVector<double,d> L(1);
>>   Dune::FieldVector<int,d> N(2);
>>   Dune::FieldVector<bool,d> periodic(true);
>>   int overlap=0;
>>   Dune::YaspGrid<d> grid(L,N,periodic,overlap);
>>   //grid.globalRefine(level);
>>   typedef Dune::YaspGrid<d>::LeafGridView GV;
>>   const GV& gv=grid.leafView();
>>   traverseGrid(gv);
>> }
>>
>> int main(int argc, char** argv)
>>     try
>>     {
>>       Dune::MPIHelper::instance(argc,argv);
>>       int level=2;
>>       if(argc==2)
>>         sscanf(argv[1],"%d",&level);
>>       dowork(level);
>>       std::cout<<"\n;\n";
>>       return 0;
>>     }
>>     catch (Dune::Exception& e){
>>       std::cerr << "Dune reported error: " << e << std::endl;
>>     }
>>     catch (boost::exception& e){
>>       std::cerr << boost::diagnostic_information(e);
>>     }
>>     catch (...){
>>       std::cerr << "Unknown exception thrown!" << std::endl;
>>     }
>
>