This is needed to guarantee that the levels of neighboring cells differ at most by 1 also the periodic boundaries. This is guaranteed by default within the domain but not at the boundaries (only if pbcs have been applied).
PM On Wednesday, 27 October 2021 at 11:23:15 UTC+2 [email protected] wrote: > A follow-up question will be: what if we only apply the periodicity after > copying the serial mesh to p:f:t? Why it is necessary to apply periodicity > before copying the mesh as well? > > On Wednesday, 20 October 2021 at 22:29:26 UTC+2 peterrum wrote: > >> The reason why you need to add periodicity twice (independent of the type >> of the base tria) is that during copying the mesh the periodicity is not >> applied. When I implemented p:f:T I did not find an easy way to do it. >> >> Hope this helps, >> PM >> >> On Wednesday, 20 October 2021 at 16:19:28 UTC+2 [email protected] >> wrote: >> >>> Hello, >>> >>> Many thanks for the reply. >>> >>> I went through these examples, it is quite helpful to understand the >>> concept of copying a serial or distributed triangulation to fully >>> distributed triangulation. I wanted some clarity in one of the examples >>> (*copy_serial_tria_04 >>> <https://github.com/dealii/dealii/blob/master/tests/fullydistributed_grids/copy_serial_tria_04.cc>*) >>> of >>> triangulation with periodic faces: >>> >>> >>> >>> >>> >>> >>> >>> >>> *// create serial triangulationTriangulation<dim> >>> basetria;GridGenerator::subdivided_hyper_cube(basetria, n_subdivisions);// >>> new: add periodicy on serial >>> meshadd_periodicy(basetria);basetria.refine_global(n_refinements);GridTools::partition_triangulation_zorder(Utilities::MPI::n_mpi_processes(comm), >>> >>> basetria);* >>> >>> >>> >>> >>> >>> *// create instance of pftparallel::fullydistributed::Triangulation<dim> >>> tria_pft(comm);// extract relevant information form serial >>> triangulationauto construction_data =* >>> >>> >>> >>> >>> >>> *TriangulationDescription::Utilities::create_description_from_triangulation(basetria, >>> >>> comm);// actually create >>> triangulationtria_pft.create_triangulation(construction_data);// new: add >>> periodicy on fullydistributed mesh (!!!)add_periodicy(tria_pft);* >>> >>> >>> As highlighted in red color, may I know why periodicity is added twice >>> to the triangulation (before and after copying )? And will same applies >>> for a p:d:t basetria as well because no such example was discussed for this >>> case? >>> >>> Thanks and Regards, >>> Aditya >>> On Tuesday, 19 October 2021 at 22:02:45 UTC+2 peterrum wrote: >>> >>>> Hi Aditya, >>>> >>>> what you found is the old development (and proof of concept) code. >>>> Since that time, we have integrated p:f:T into deal.II. All the tests >>>> there >>>> have been converted to proper deal.II tests in the folder >>>> https://github.com/dealii/dealii/tree/master/tests/fullydistributed_grids. >>>> I suggest you to take a look at the test copy_serial_tria_*.cc and >>>> copy_distributed_tria_*.cc. With increasing number of indices the tests >>>> become more complex (adding PBC, MG, ...). So start with 01 ;) >>>> >>>> Hope that helps! >>>> PM >>>> >>>> On Tuesday, 19 October 2021 at 21:50:18 UTC+2 [email protected] >>>> wrote: >>>> >>>>> Hi, >>>>> >>>>> I was trying to understand dealii-pft/step-5 >>>>> <https://github.com/peterrum/dealii-pft/blob/master/step-5/src/main.cpp> >>>>> and implement something similar in my code, but there are some things >>>>> which >>>>> are confusing for me! >>>>> >>>>> As from line 55, Author copies the data from serial triangulation to a >>>>> fully distributed triangulation and then reinitialize it: >>>>> >>>>> >>>>> >>>>> >>>>> >>>>> >>>>> >>>>> *// create instance of >>>>> pftparallel::fullydistributed::Triangulation<dim> tria_pft(comm);// >>>>> extract >>>>> relevant information form serial triangulationauto construction_data >>>>> =parallel::fullydistributed::Utilities::copy_from_triangulation(basetria, >>>>> tria_pft);// actually create >>>>> triangulationtria_pft.reinit(construction_data)*; >>>>> >>>>> As I tried to implement similar, firstly I cannot find any function >>>>> called *reinit *for triangulation file. Also, the way *construction_data >>>>> *has been created, I was unable to find function called >>>>> *copy_from_triangulation >>>>> *for doing the same! I am assuming that this is because I have not >>>>> included proper header file. >>>>> >>>>> I can see that author has included <*deal.II/distributed/tria_util.h*> >>>>> in the code, but it was not available atleast in dealii version 9.2.0 or >>>>> 9.3.0 . >>>>> >>>>> Can anyone suggest something for this? >>>>> >>>>> Looking forward to some valuable suggestions. >>>>> >>>>> >>>>> Thanks and Regards, >>>>> Aditya >>>>> >>>>> >>>>> >>>>> -- The deal.II project is located at http://www.dealii.org/ For mailing list/forum options, see https://groups.google.com/d/forum/dealii?hl=en --- You received this message because you are subscribed to the Google Groups "deal.II User Group" group. To unsubscribe from this group and stop receiving emails from it, send an email to [email protected]. To view this discussion on the web visit https://groups.google.com/d/msgid/dealii/ab5b6985-b1f1-42e1-b434-4716319c96aan%40googlegroups.com.
