Finite Element Method - 2D Mesh Generator - Metfem2D

Generator Optimization via the Metropolis method

Internal mesh elements

Presented method is based on the following algorithm[1, 2]:

1. Define the element size h and consequently the element area A.
2. Initialize the configuration of triangles and then select the internal nodes Pint = { pi: pi ∈ P ∧ pi ∉ Γ } i. e. these nodes does not belong to the domain boundary Γ.
3. For each node pi in Pint find its subdomain Ωi defined as a set of triangles Ti to which the node pi belongs.
4. Perform the Metropolis approach to every internal node pi within its subdomain Ωi. The Metropolis algorithm is adopted in order to adjust an area of each triangle in the node's subdomain to prescribed value A by shifting the position of the node pi (Fig. 5 demonstrates robustness of the Metropolis approach; compare the node distribution in a) and in b)). That adjustment is governed by the following rules:
1. Find an area of each triangle Ak (where k = 1, 2, …, K) in Ωi together with the vectors rji = pi - pj for each pj ∈ Ωi connected to node pi
2. Calculate the length of each triangle edge |rji| and its deviation δ |rji| from the designed element size h i. e. δ |rji| = |rji| - h
3. Calculate the new position of the node pnewi as

 pnewi = pi - ∑ j Fj δ |rji| versor(rji)
where versor(rji) (i. e. an unit vector) has the standard meaning as rji/|rji| and Fj are weights corresponding to magnitude of j-th force applying to node pi. Finally, they were set to the constant value F.
4. Find an area of each triangle Anewk ∈ Ωi after shifting pi → pnewi
5. Apply an energetic measure E to a sub – mesh Ωi. That quantity could be understand in terms of a square deviation of a mesh element area from the prescribed element area A. Therefore, in the presented paper the δ E is defined as a sum of a discrepancy between each triangle area Ak and A after moving node pi and prior it, respectively
 δE = ∑ k ( (Anewk - A)2 - (Ak - A)2 ) k = 1, 2, …, K.

If the obtained value of an energetic change is lower than zero the change is accepted. Otherwise, the Metropolis rule is applied i. e. the following condition is checked

e-δE⁄T > r
where r is an uniformly distributed random number on the unit interval (0, 1) and T denotes temperature.
6. The above – presented algorithm is repeated unless an assumed tolerance will be achieved.
In order to reach a better convergence of the presented method several other improvements could be adopted. For instance, the change in the length of the triangle edge could be an additional measure of mesh approximation goodness. That condition will ensure a lack of elongated mesh elements i. e. elements with very high ratio of its edge lengths (to see such skinny elements look at Fig. 5a)).

Figure shows an application of the Metropolis algorithm. Picture a) presents initial positions of new nodes just after generating them whereas picture b) shows their positions after node shifts according to the procedure described in Sec. 4.1 with the following two values of the force strength Fj: 0.006 and 0.1 applied to each internal node j = 1, …, 7 and temperature set as T = 0.01. The table presents the total number of Metropolis steps that was required to obtain the final result shown in b).

References

1. ^ O. C. Zienkiewicz, R. L. Taylor and J. Z. Zhu, The Finite Element Method: Its Basis and Fundamentals, Sixth edition., Elsevier 2005
2. ^ N. Metropolis, S. Ulam, The Monte Carlo Method, J. Amer. Stat. Assoc., 44, No. 247., pp. 335-341, 1949