| Thread Previous • Date Previous • Date Next • Thread Next |
On 09/01/09 09:38 AM, Johan Jansson wrote:
Hi!This looks interesting, could you please post references to the papers you mentioned about the algorithm?
Basic concept is from Rivara, but I'd like to avoid recursive algorithms, which has usually efficiency problems.
I've found in some papers [1], where they use similar algorithm but with iteration procedure. Propably something like this is also implemented in FEtk.org (GAMer package). I write "probably" because the code is still unpublished, and I can't find any documentation about it.
[1] N. A. Baker, D. Sept, M. J. Holst, J. A. McCammon: The adaptive multilevel finite element solution of the Poisson–Boltzmann equation on massively parallel computers,IBM J. RES. & DEV. VOL. 45 NO. 3/4 MAY/JULY 2001
Algorithm with I've implemented. It wasn't directly taken from any paper:
while( there are any marked cells ):
forech cell which is marked:
find the longest edge
if the edge hasn't been bisected:
bisect longest edge e
remember bisected edge and new vertex
else :
use vertex which was remembered
bisect cell
foreach bisected edge:
mark all cells which incidence the edge
The weak point of the algorithm is the storing of bisected edges, and
serching in them. I use double index STL mapping for that purpose -
which works pretty well, but that's the slowest part of the algorithm.
I have some ideas, which can make it faster, but it hasn't been tested yet. Stay tuned ...
The mesh quality is fine, because _only the longest edge_ is bisected.
I've been working on an implementation of the Rivara recursive bisection algorithm, which looks like this:function bisect(K): Split longest edge e While K_i(e) is non-conforming bisect(K_i)where K_i(e) are the cells incident on edge e and a non-conforming cell K is defined as having a neighboring cell that has a hanging vertex incident on an edge of K.The implementation defines a dynamic mesh data structure, which I think is required for this algorithm. The implementation can be found here:http://www.nada.kth.se/cgi-bin/jjan/hgwebdir.cgi/unicorn-kth/file/6e94418ae942/lib/unicorn/MeshAdaptivity.h http://www.nada.kth.se/cgi-bin/jjan/hgwebdir.cgi/unicorn-kth/file/6e94418ae942/lib/MeshAdaptivity.cppI will try to put these in DOLFIN as soon as possible, or if anybody wants to help me.This is also a robust algorithm with bounds on mesh quality in 2D and with a variant (Arnold, no longer longest edge) in 3D. In our tests, the plain Rivara algorithm seems sufficient in 3D as well. The question is though, what is the relation between the Rivara algorithm to the algorithm you've implemented? Is your algorithm just an iterative version of the recursion, or something else? If it's the same algorithm, perhaps there can be an advantage of keeping the fixed mesh.
Yes. I think that my algorithm is mainly an iterative version of recursion. The advantage is the simplicity, robustness, and mesh quality is much better than the old simplest bisection used in Dolfin.
We've also started working together with Gaetan Compere in the Gmsh group with integrating their work on local mesh operations (bisection, collapse, edge/face swapping) with FEniCS. See below for a reference to their methods/implementation. This work is very interesting, since then bad cells can be improved, and you also get coarsening. Depending on robustness, this could make algorithms like Rivara obsolete, since you then don't need to guarantee quality for refinement.
That sounds interesting. It would be great to have possibility to manipulate the mesh on the local level.
regrds. BArtek
For the time being though, it would be useful to have these robust algorithms in DOLFIN.These are the papers I mentioned above:Gaetan Compere, Emilie Marchandise, and Jean-Francois Remacle, Transient adaptivity applied to two-phase incompressible flows, J. Comput. Phys., 227 (2008), pp. 1923-1942.M-C. Rivara, Local modification of meshes for adaptive and/or multigrid finite-element methods, Journal of Computational and Applied Mathematics, 36 (1992), pp. 7889.D. N. Arnold, A. Mukherjee, and L. Pouly, Locally adapted tetrahedral meshes using bisection, SIAM Journal on Scientific Computing, 22 (2000), pp. 431448.Best, Johan Bartosz Sawicki wrote:Please find attached patch with new iterative refinement algorithm.The iteration process is robust, usually less that 10 iterations is needed to converge. I've tested it with triangular and tetrahedral meshes up to 100k cells.This is my first bigger commit into the project. So please read it carefully. I tried to write fast and efficient code, but maybe someone has any ideas how to improve it. Waiting for any opinions...regrds. BArtek On 07/01/09 04:15 PM, Bartosz Sawicki wrote:On 31/12/08 12:51 AM, Anders Logg wrote:On Tue, Dec 30, 2008 at 04:54:05PM -0700, Bartosz Sawicki wrote:I'm starting to work on implementation of more advanced algorithm for mesh refinement. At this moment the simplest version of the Longest-Edge-Bisection is used inside the Mesh::refine().I plan to develop some algorithm based on Edge-Bisection with two set of markers and iterative swapping. There are some papers about this method. It is thought to be robust and produce good quality mesh.This sounds good. Our current refinement algorithm sometimes produces bad quality mesh (very thin triangles) if called many times without intermediate calls to mesh.smooth().Today I was stopped by the first problem related with edge global index. We know that every bisection of triangle produces 3 new edges, and removes 1 old edge. The question is, if there is method to find what will be the global index of the new edges. If I understand correctly edge index can created (mesh.init(1)) only after mesh object is created and fulfilled by Vertexes and Cells. Maybe I'm wrong and there is simple way to get edge index even without calling mesh.init(1)?Unfortunately no. It would be very good if there were a way to knowglobal edge numbers. It would also help with the parallel implementationcurrently in progress. But there might be some clever numbering strategy that could work. If anyone knows of such a numbering strategy (uniquely assigning global edge numbers based on global vertex numbers) we could implement that as a post-processing step in TopologyComputation. One first suggestion would be to numberthe first edge of the first vertex 0, then the second edge connected to the first vertex etc.I have finished first working version of new mesh refinement algorithm. The code is quite dirty now, but I will try to make it more readable for everyone. I'll send you a patch.Edge numbering problem was solved by using std:map<Array<int>, int, cmp>. I've used two vertex number as a no-changing edge identifier. This works quite fast, but I'm open for new ideas. Maybe someone has better idea how to build MeshFunction connected with edges, and transform it into different mesh?Similar problem would be crucial for implementation of MeshFunctions transformations connected with mesh refinement. As boundary indicators, and material properties should be preserved when we change the mesh.I plan to implement this in the next step. regrd. BArtek------------------------------------------------------------------------_______________________________________________ DOLFIN-dev mailing list DOLFIN-dev@xxxxxxxxxx http://www.fenics.org/mailman/listinfo/dolfin-dev_______________________________________________ DOLFIN-dev mailing list DOLFIN-dev@xxxxxxxxxx http://www.fenics.org/mailman/listinfo/dolfin-dev------------------------------------------------------------------------ _______________________________________________ DOLFIN-dev mailing list DOLFIN-dev@xxxxxxxxxx http://www.fenics.org/mailman/listinfo/dolfin-dev
| Thread Previous • Date Previous • Date Next • Thread Next |
