Roy, thanks for the very helpful demonstration of the coarsen and refine flagging in libmesh and clarifying my doubts on the view of a mesh.
I have another question that stems from your explanation. If I start from a mesh and refine it once and coarsen it once, do you still store the old refined dofs in the mesh and just make them inactive. When exactly would you release this allocation ? Also, this does call the EquationSystem.reinit() in order to associate all system variables to the new dofs. Is this correct ? I can see this getting very messy quickly since there is way too much bookkeeping that is needed to get the V-cycle right. > Ideally the external user shouldn't need to do anything that he isn't > already other than set some solver options or choose a different > solver object. I do not use Libmesh or Petsc as a black-box but have inherited from the structures exposed by these libraries and have my own wrappers on top to give me flexibility. So I do not mind getting my hands dirty a little, if need be, in order to get a geometric multi-grid framework working. The question though is how much of this can be done without changing too many interfaces in existing libmesh objects ?! Again, I appreciate all the help Roy. > -----Original Message----- > From: Roy Stogner [mailto:[EMAIL PROTECTED] > Sent: Wednesday, November 12, 2008 7:02 PM > To: Vijay S. Mahadevan > Cc: [email protected] > Subject: Re: [Libmesh-users] Multigrid techniques with libmesh > > > On Wed, 12 Nov 2008, Vijay S. Mahadevan wrote: > > > I am curious on the concept of the 'view' of a mesh. > > > Basically, we store hierarchic meshes like this: > > *-----------*-----------* > *-----*=====*=====*=====* > *==*==* > > Where the mesh is 1D, there are 5 active elements and 3 ancestor > elements. > > We can also store meshes like this: > > *-----------*===========* > *=====*=====*.....*.....* > *..*..* > > Where there are 3 active elements, 1 ancestor element, and 4 > "subactive" elements". Looping over active elements just gets you the > ones in the active view, but the subactive elements and their > associated Dof objects still exist. We use this for coarsening > non-Lagrange elements properly, but it could make geometric multigrid > possible too. > > > For example when I have an unstructured grid, and when I coarsen > > uniformly twice and refine uniformly twice, would I get the exact > > same mesh?! > > No. For instance, the second mesh above is a uniform coarsening (so > far as such a thing exists) of the first mesh, but uniformly refining > it will produce two new elements that weren't in the first mesh. > > We could, however, create a "rerefine" method that jumps you straight > back to the finest mesh by making any active element with children > into an ancestor and making any subactive element without children > active. > > > If the dof_map and mesh can provide views that are consistent with > > the refinement or coarsening based on an initial mesh, then I do think > > this should be possible without too much programming effort for an > > external user. > > Ideally the external user shouldn't need to do anything that he isn't > already other than set some solver options or choose a different > solver object. > > But that implies the existence of an internal developer with the time > and expertise to do the library work, and such a person may not exist > right now. > > > Of course, you would also have to propagate the mesh_data (Or again, > > am I the only one using this ?) which stores material > > information/element or create a view of it also. > > You may be the only one using this. ;-) But yes, we'd want to handle > its restriction as well. > --- > Roy ------------------------------------------------------------------------- This SF.Net email is sponsored by the Moblin Your Move Developer's challenge Build the coolest Linux based applications with Moblin SDK & win great prizes Grand prize is a trip for two to an Open Source event anywhere in the world http://moblin-contest.org/redirect.php?banner_id=100&url=/ _______________________________________________ Libmesh-users mailing list [email protected] https://lists.sourceforge.net/lists/listinfo/libmesh-users
