Ortega et al 2011a - Revision history

Cinmemj: Cinmemj moved page Draft Samper 265118390 to Ortega et al 2011a

2018-05-29T12:51:48Z

Cinmemj moved page Draft Samper 265118390 to Ortega et al 2011a

Cinmemj at 12:48, 29 May 2018

2018-05-29T12:48:36Z

Cinmemj at 12:42, 29 May 2018

2018-05-29T12:42:03Z

Cinmemj at 12:16, 29 May 2018

2018-05-29T12:16:57Z

Cinmemj at 12:15, 29 May 2018

2018-05-29T12:15:46Z

Cinmemj at 12:14, 29 May 2018

2018-05-29T12:14:18Z

Cinmemj at 12:11, 29 May 2018

2018-05-29T12:11:35Z

Show changes

Cinmemj at 11:39, 29 May 2018

2018-05-29T11:39:41Z

Cinmemj at 11:39, 29 May 2018

2018-05-29T11:39:20Z

Cinmemj at 11:39, 29 May 2018

2018-05-29T11:39:02Z

@@ Line 391: / Line 391: @@
 Once the new positions of the discrete points are obtained the boundary normal vectors and other discretization related data are updated accordingly. In addition, although the cloud connectivity remains unchanged after domain deformation, the relative position of its points does not. Consequently, the local approximation must be re-computed in the affected clouds.
-=6. ADAPTIVITY=
+==6. ADAPTIVITY==
 Unsteady compressible flow problems typically involve evolving discontinuities resulting from flow unsteadiness or body motions. Accurate solutions require these localized flow features to be properly resolved and this requires a domain resolution which could be unaffordable in many engineering analyses. An efficient way to deal with these problems with accuracy while keeping the computational cost acceptable is by using ''h''-adaptivity. Moreover, ''h''-adaptivity constitutes a good complement to domain deformation techniques because it can regenerate in a local manner highly distorted zones in the analysis domain improving the quality of the discretization. The reduced discretization restrictions found in meshless methods makes them an ideal context to implement ''h''-adaptive strategies. In this work the main lines of the methodology developed by the authors in [9, 30] is followed, although some modifications are introduced to increase the effectiveness and robustness of the technique. Next, the main aspects of the adaptive strategy are described.
-==6.1 Refinement/coarsening indicator==
+===6.1 Refinement/coarsening indicator===
 Local clouds in the analysis domain where either refinement or coarsening is required are identified as follows. First, a sensor measuring in an approximate manner the curvature of the density solution field is computed for all points in the domain by
@@ Line 405: / Line 405: @@
 {| style="text-align: left; margin:auto;width: 100%;"
 |-
-| style="text-align: center;" | <math>r_{i} \; =\; \sum _{j=1}^{n_{nei} }\left|{l}_{ij} \cdot \left({\nabla }\rho _{j} \; -\; {\nabla }\rho _{i} \right)\right|  </math>
+| style="text-align: center;" | <math>r_{i} \; =\; \sum _{j=1}^{n_{nei} }\left|{\bf l}_{ij} \cdot \left({\nabla }\rho _{j} \; -\; {\nabla }\rho _{i} \right)\right|  </math>
 |}
 | style="width: 5px;text-align: right;white-space: nowrap;" | (23)
 |}
-being <math>{n}_{nei}</math> the number of points in the layer of nearest neighbours of <math>{\bf x}_{i}</math> and <math display="inline">{l}_{ij} \; =\; {x}_{j} \; -\; {x}_{i} </math>. Next, a smoothed non-dimensional refinement indicator is computed at each star point as
+being <math>{n}_{nei}</math> the number of points in the layer of nearest neighbours of <math>{\bf x}_{i}</math> and <math display="inline">{\bf l}_{ij} \; =\; {\bf x}_{j} \; -\; {\bf x}_{i} </math>. Next, a smoothed non-dimensional refinement indicator is computed at each star point as
 <span id='ZEqnNum698202'></span>
@@ Line 426: / Line 426: @@
 Finally, the logarithm of the smoothed indicator <span id='cite-ZEqnNum698202'></span>[[#ZEqnNum698202|(24)]] is taken in order to compress the scale of the distribution and its mean value <math display="inline">\varphi _{m} </math> and standard deviation <math display="inline">s_{\varphi } </math> are computed. These values are used to determine local clouds in which refinement or coarsening is required: the star point <math display="inline">{\bf x}_{i}</math> is tagged for refinement when <math display="inline">\varphi _{i} \; >\; \varphi _{m} \; +\; N_{ref} s_{\varphi } </math> and, conversely, the point <math display="inline">{\bf x}_{i}</math> is marked to be removed if <math display="inline">\varphi _{i} \; <\; \varphi _{m} \; -\; N_{rem} s_{\varphi } </math>. The threshold parameters <math display="inline">N_{ref} </math> and <math display="inline">N_{rem} </math> must be set according to the problem under study, typical values employed in the examples shown in this work are <math display="inline">N_{ref} =1</math> and <math display="inline">N_{rem} =0</math>. It was observed that this indicator performs better than that employed previously in [9, 30] in capturing strong and smooth features of the flow; in addition, it presents a reduced dependence of the user defined threshold parameters.
-==6.2 Removal and insertion of points==
+===6.2 Removal and insertion of points===
 The adaptive stage begins with the removal of points tagged for deletion. The coarsening of surface grids is performed by collapsing edges with marked nodes and interior volume points are simply deleted from the vertices list. After the removal of points the data structure is updated. It should be noticed that the removal of points is restricted only to the existing points that have been inserted in prior refinement levels.
@@ Line 432: / Line 432: @@
 The surface grids are refined after the coarsening stage. To this end, surface elements having tagged all its nodes are selected and a new point is inserted at its centroid (and the element is subdivided) if the distance from the latter to any other point in the discretization is greater than a minimum element size (set as the minimum nodal distance of its nodes). The minimum nodal distance, which determines the desired level of refinement, is computed for the original points as a percentage of the distance between the point and the nearest neighbour in its cloud; for new points added during refinement stages this distance is obtained from interpolation of the original distribution. Once new boundary points are inserted, its positions (originally coincident with the centroid of the underlying element) are slightly improved by interpolation [31] to avoid excessive faceting of the surfaces in successive element subdivisions. Finally, edge swapping is applied on the affected boundary elements to improve its connectivity and facilitate the insertion and removal of points in future adaptive passes.
-<span id='OLE_LINK1'></span><span id='OLE_LINK2'></span>In the interior volume, when a cloud of points is selected to be refined, the Voronoi vertices surrounding the star point <math display="inline">{\bf x}_{i}</math> are computed by means of its Delaunay grid of nearest neighbours. Next, a new point is added at each Voronoi vertex location if the distance from the latter to any existing point is greater than the minimum distance for the cloud and if, in addition, boundary constraints are satisfied (basically, the ray from the star point to any new point added to the cloud cannot pierce any boundary surface). In order to perform the many spatial search operations required during the adaptive stage efficiently, a dynamically updated octree structure is employed.
+In the interior volume, when a cloud of points is selected to be refined, the Voronoi vertices surrounding the star point <math display="inline">{\bf x}_{i}</math> are computed by means of its Delaunay grid of nearest neighbours. Next, a new point is added at each Voronoi vertex location if the distance from the latter to any existing point is greater than the minimum distance for the cloud and if, in addition, boundary constraints are satisfied (basically, the ray from the star point to any new point added to the cloud cannot pierce any boundary surface). In order to perform the many spatial search operations required during the adaptive stage efficiently, a dynamically updated octree structure is employed.
-==6.3 Update==
+===6.3 Update===
-<span id='OLE_LINK4'></span>Once coarsening and refinement steps are finished, a few steps of a Laplacian smoothing are carried out on the affected area. After that, the clouds of points and the local approximation concerning the new points are constructed. In addition, the data concerning existing clouds of points affected by deletion, insertion of new points or smoothing are re-constructed. Finally, the flow and other problem variables at the new points are obtained as an average of the variables at their previously existing first nearest neighbours.
+Once coarsening and refinement steps are finished, a few steps of a Laplacian smoothing are carried out on the affected area. After that, the clouds of points and the local approximation concerning the new points are constructed. In addition, the data concerning existing clouds of points affected by deletion, insertion of new points or smoothing are re-constructed. Finally, the flow and other problem variables at the new points are obtained as an average of the variables at their previously existing first nearest neighbours.
 =7. NUMERICAL APPLICATIONS=

@@ Line 360: / Line 360: @@
 Although system [[#eq-20|(20)]] has excellent properties and can be solved by a few Gauss-Seidel iterations with little computational cost, numerical experiments have demonstrated that in some cases this step could be skipped assuming <math display="inline">{\bf U}^{h} \, \approx \, \hat{{\bf U}}</math> without causing any negative impact on the accuracy of the numerical solution. This behaviour could be explained to a large extent due to the fact that <math display="inline">(\cdot )^{h}</math> and (<math display="inline">\hat{\cdot }</math>) parameters become very closer in well-behaved finite point approximations (the shape functions tend to interpolate point data).
-=5. POINTS MOVEMENT STRATEGY=
+==5. POINTS MOVEMENT STRATEGY==
 <span id='OLE_LINK12'></span><span id='OLE_LINK13'></span>The simulation of unsteady problems involving deforming or moving bodies requires the domain discretization to conform continuously to the instantaneous body shape. According to the analysis geometry and the motion involved, the domain deformation technique could vary from simple rigid rotations or translations to more complex operations if the body modifies its shape or relative motions between several bodies must be accounted for. An overview of the strategies usually proposed to deal with these issues is presented in [20]. In this work a spring network approach is adopted [28, 29]. Therefore, the displacement of any interior point <math display="inline">{\bf x}_i</math> in response to instantaneous displacements of body points are obtained by enforcing the static equilibrium of the forces exerted by all the points <math display="inline">{\bf x}_j</math> connected through the <math display="inline">{\bf x}_i</math>'s layer of nearest neighbours (outer boundary points are considered to be fixed). This leads to a system of equations written in terms of displacements which is solved by Jacobi iterations,
@@ Line 369: / Line 369: @@
 {| style="text-align: center; margin:auto;"
 |-
-| <span id='ZEqnNum360001'></span>  <math>{\delta }_i^{k+1}\mbox{ }=\mbox{ }\frac{{\sum }_jk_{ij}\mbox{ }{\delta }_j^k}{{\sum }_jk_{ij}}</math>
+| <span id='ZEqnNum360001'></span>  <math>{\boldsymbol\delta }_i^{k+1}\mbox{ }=\mbox{ }\frac{{\sum }_jk_{ij}\mbox{ }{\boldsymbol\delta }_j^k}{{\sum }_jk_{ij}}</math>
 |}
 | style="width: 5px;text-align: right;white-space: nowrap;" | (21)
@@ Line 381: / Line 381: @@
 {| style="text-align: center; margin:auto;"
 |-
-| <math>x_i^{new}=x_i^{old}+{\boldsymbol\delta }_i^{}</math>
+| <math>{\bf x}_i^{new}={\bf x}_i^{old}+{\boldsymbol\delta }_i</math>
 |}
 | style="width: 5px;text-align: right;white-space: nowrap;" | (22)
@@ Line 387: / Line 387: @@
-and its velocities are simply estimated by  <math display="inline">{\bf w}_i={\delta }_i/\Delta t</math> being  <math display="inline">\Delta t</math> the physical time increment employed in the simulation. Notice that an accurate solution of the system <span id='cite-ZEqnNum360001'></span>[[#ZEqnNum360001|(21)]] is not required; in this work 10-30 iterations showed enough to propagate the body displacements inside the domain achieving a smooth distortion of the discretization. For very large body displacements the link stiffness can be raised by increasing parameter ''p.''
+and its velocities are simply estimated by  <math display="inline">{\bf w}_i={\boldsymbol\delta }_i/\Delta t</math> being  <math display="inline">\Delta t</math> the physical time increment employed in the simulation. Notice that an accurate solution of the system <span id='cite-ZEqnNum360001'></span>[[#ZEqnNum360001|(21)]] is not required; in this work 10-30 iterations showed enough to propagate the body displacements inside the domain achieving a smooth distortion of the discretization. For very large body displacements the link stiffness can be raised by increasing parameter ''p.''
 Once the new positions of the discrete points are obtained the boundary normal vectors and other discretization related data are updated accordingly. In addition, although the cloud connectivity remains unchanged after domain deformation, the relative position of its points does not. Consequently, the local approximation must be re-computed in the affected clouds.

@@ Line 345: / Line 345: @@
 where <math display="inline">{{\mathcal C}}</math> denotes the Courant number, ''c''<math display="inline">{}_{i}</math> is the speed of sound and the rest of the variables have been previously defined. In addition to local time steeping, implicit residual smoothing [25] is applied to accelerate the convergence of the system [[#eq-17|(17)]] to the steady state. Notice that the explicit treatment of the term <math display="inline">3{\bf U}_{{i}}^{{*}} /2\Delta t</math> could lead to numerical instabilities if the physical time step is small [26]. Fortunately, this problem can be easily overcome by treating implicitly that term in Eq. [[#eq-18|(18)]], ''cf.'' [27] for further details.
-The difference between <math display="inline">(·)^{h}</math> and (<math display="inline">\hat{\cdot }</math>) parameters has been explained in Section 2. Taking into account that <math display="inline">{\bf R}({\bf U}_{i})</math> is a function of <math>{\bf U}_{j}^{h}\nabla {\bf x}_{j}\in \Omega_{i}</math>, the following linear system has to be solved to recover the <math display="inline">(·)^{h}</math> parameters from the approximate ones
+The difference between <math display="inline">(\cdot )^{h}</math> and (<math display="inline">\hat{\cdot }</math>) parameters has been explained in Section 2. Taking into account that <math display="inline">{\bf R}({\bf U}_{i})</math> is a function of <math>{\bf U}_{j}^{h}\nabla {\bf x}_{j}\in \Omega_{i}</math>, the following linear system has to be solved to recover the <math display="inline">(\cdot )^{h}</math> parameters from the approximate ones
 <span id="eq-20"></span>
@@ Line 358: / Line 358: @@
 |}
-Although system [[#eq-20|(20)]] has excellent properties and can be solved by a few Gauss-Seidel iterations with little computational cost, numerical experiments have demonstrated that in some cases this step could be skipped assuming <math display="inline">{\bf U}^{h} \, \approx \, \hat{{\bf U}}</math> without causing any negative impact on the accuracy of the numerical solution. This behaviour could be explained to a large extent due to the fact that <math display="inline">(·)^{h}</math> and (<math display="inline">\hat{\cdot }</math>) parameters become very closer in well-behaved finite point approximations (the shape functions tend to interpolate point data).
+Although system [[#eq-20|(20)]] has excellent properties and can be solved by a few Gauss-Seidel iterations with little computational cost, numerical experiments have demonstrated that in some cases this step could be skipped assuming <math display="inline">{\bf U}^{h} \, \approx \, \hat{{\bf U}}</math> without causing any negative impact on the accuracy of the numerical solution. This behaviour could be explained to a large extent due to the fact that <math display="inline">(\cdot )^{h}</math> and (<math display="inline">\hat{\cdot }</math>) parameters become very closer in well-behaved finite point approximations (the shape functions tend to interpolate point data).
 =5. POINTS MOVEMENT STRATEGY=

@@ Line 291: / Line 291: @@
 |}
-where <math display="inline">-{\bf R}({\bf U}_{i})</math> is the right-hand side of the system [[#eq-13|(13)]] and the subscripts <math display="inline">(·)^{h}</math> and <math display="inline">(\hat{\cdot })</math> have been omitted for the sake of clarity. Next, for a given increment in physical time <math display="inline">\Delta t</math>, a modified (unsteady) residual is defined from Eq. [[#eq-15|(15)]] as
+where <math display="inline">-{\bf R}({\bf U}_{i})</math> is the right-hand side of the system [[#eq-13|(13)]] and the subscripts <math display="inline">(\cdot )^{h}</math> and <math display="inline">(\hat{\cdot })</math> have been omitted for the sake of clarity. Next, for a given increment in physical time <math display="inline">\Delta t</math>, a modified (unsteady) residual is defined from Eq. [[#eq-15|(15)]] as
 <span id="eq-16"></span>

@@ Line 291: / Line 291: @@
 |}
-where <math display="inline">-{\bf R}({\bf U}_{i})</math> is the right-hand side of the system [[#eq-13|(13)]] and the subscripts <math display="inline">{(·)}^{h}</math> and (<math display="inline">\hat{\cdot }</math>) have been omitted for the sake of clarity. Next, for a given increment in physical time <math display="inline">\Delta t</math>, a modified (unsteady) residual is defined from Eq. [[#eq-15|(15)]] as
+where <math display="inline">-{\bf R}({\bf U}_{i})</math> is the right-hand side of the system [[#eq-13|(13)]] and the subscripts <math display="inline">(·)^{h}</math> and <math display="inline">(\hat{\cdot })</math> have been omitted for the sake of clarity. Next, for a given increment in physical time <math display="inline">\Delta t</math>, a modified (unsteady) residual is defined from Eq. [[#eq-15|(15)]] as
 <span id="eq-16"></span>
@@ Line 345: / Line 345: @@
 where <math display="inline">{{\mathcal C}}</math> denotes the Courant number, ''c''<math display="inline">{}_{i}</math> is the speed of sound and the rest of the variables have been previously defined. In addition to local time steeping, implicit residual smoothing [25] is applied to accelerate the convergence of the system [[#eq-17|(17)]] to the steady state. Notice that the explicit treatment of the term <math display="inline">3{\bf U}_{{i}}^{{*}} /2\Delta t</math> could lead to numerical instabilities if the physical time step is small [26]. Fortunately, this problem can be easily overcome by treating implicitly that term in Eq. [[#eq-18|(18)]], ''cf.'' [27] for further details.
-The difference between <math display="inline">{(·)}^{h}</math> and (<math display="inline">\hat{\cdot }</math>) parameters has been explained in Section 2. Taking into account that <math display="inline">{\bf R}({\bf U}_{i})</math> is a function of <math>{\bf U}_{j}^{h}\nabla {\bf x}_{j}\in \Omega_{i}</math>, the following linear system has to be solved to recover the <math display="inline">{(·)}^{h}</math> parameters from the approximate ones
+The difference between <math display="inline">(·)^{h}</math> and (<math display="inline">\hat{\cdot }</math>) parameters has been explained in Section 2. Taking into account that <math display="inline">{\bf R}({\bf U}_{i})</math> is a function of <math>{\bf U}_{j}^{h}\nabla {\bf x}_{j}\in \Omega_{i}</math>, the following linear system has to be solved to recover the <math display="inline">(·)^{h}</math> parameters from the approximate ones
 <span id="eq-20"></span>
@@ Line 358: / Line 358: @@
 |}
-Although system [[#eq-20|(20)]] has excellent properties and can be solved by a few Gauss-Seidel iterations with little computational cost, numerical experiments have demonstrated that in some cases this step could be skipped assuming <math display="inline">{\bf U}^{h} \, \approx \, \hat{{\bf U}}</math> without causing any negative impact on the accuracy of the numerical solution. This behaviour could be explained to a large extent due to the fact that <math display="inline">{(·)}^{h}</math> and (<math display="inline">\hat{\cdot }</math>) parameters become very closer in well-behaved finite point approximations (the shape functions tend to interpolate point data).
+Although system [[#eq-20|(20)]] has excellent properties and can be solved by a few Gauss-Seidel iterations with little computational cost, numerical experiments have demonstrated that in some cases this step could be skipped assuming <math display="inline">{\bf U}^{h} \, \approx \, \hat{{\bf U}}</math> without causing any negative impact on the accuracy of the numerical solution. This behaviour could be explained to a large extent due to the fact that <math display="inline">(·)^{h}</math> and (<math display="inline">\hat{\cdot }</math>) parameters become very closer in well-behaved finite point approximations (the shape functions tend to interpolate point data).
 =5. POINTS MOVEMENT STRATEGY=

← Older revision	Revision as of 12:51, 29 May 2018
(No difference)

@@ Line 169: / Line 169: @@
 Even though not any point discretization is valid for meshless computations (especially when the problem exhibits complex geometrical/flow features), the fact that topology is not a matter of concern allows achieving proper point discretizations faster and cheaper. Hence, it is expected that any point generation technique employed makes the most of this advantage. A highly effective technique is that developed by Calvo [17], in which the point generation starts from a user-supplied boundary grid and inserts new points in the centre of empty spheres filling &#x03a9; through an optimization driven point insertion procedure. This incremental quality technique based on unconstrained Delaunay tetrahedralization allows achieving a quality point discretization with approximated cost O(n) [18, 19]. In addition to the low computational cost, this technique has other advantages that reduce model developing time significantly. For instance, the input boundary mesh may contain any type of elements and be non-conformal. Furthermore, the point spacing inside the domain does not need to be specified, it is automatically assigned by computing linear variations between grid boundary sizes. These characteristics make the Calvo’s point generation technique highly competitive with respect to other approaches used in mesh-based discretization methods.
-==2.3 Construction of local clouds==
+===2.3 Construction of local clouds===
 Cloud construction techniques in LSQ-based meshless approximations must be aimed at achieving point distributions which facilitate (from a mathematical point of view) the solution of the minimization problem. Moreover, these procedures should be able to deal with all the geometrical features found in practical application problems with robustness and a low computational cost. In this regard, we follow the technique proposed by Löhner ''et al.'' [6], which has proven effective. The procedure is summarized in the next.

@@ Line 165: / Line 165: @@
 with ''d''<math display="inline">{}_{j}</math> = '''x''' ''<math>{}_{j}</math>-'' '''x''' ''<math>{}_{i}</math>'', <math display="inline">\alpha </math> = <math display="inline">\beta </math>/''w'' and <math display="inline">\beta </math> = <math display="inline">\gamma </math> ''d<math>{}_{max}</math>'' (<math display="inline">\gamma >1</math>). The parameters ''w, k'' and <math display="inline">\gamma </math> govern the functional shape of the weighting function and could have important effects on the resulting numerical approximation. Some guidelines for its setting are given in [15].
-==2.2 Discretization of the analysis domain: the point generator==
+===2.2 Discretization of the analysis domain: the point generator===
 Even though not any point discretization is valid for meshless computations (especially when the problem exhibits complex geometrical/flow features), the fact that topology is not a matter of concern allows achieving proper point discretizations faster and cheaper. Hence, it is expected that any point generation technique employed makes the most of this advantage. A highly effective technique is that developed by Calvo [17], in which the point generation starts from a user-supplied boundary grid and inserts new points in the centre of empty spheres filling &#x03a9; through an optimization driven point insertion procedure. This incremental quality technique based on unconstrained Delaunay tetrahedralization allows achieving a quality point discretization with approximated cost O(n) [18, 19]. In addition to the low computational cost, this technique has other advantages that reduce model developing time significantly. For instance, the input boundary mesh may contain any type of elements and be non-conformal. Furthermore, the point spacing inside the domain does not need to be specified, it is automatically assigned by computing linear variations between grid boundary sizes. These characteristics make the Calvo’s point generation technique highly competitive with respect to other approaches used in mesh-based discretization methods.

@@ Line 148: / Line 148: @@
 :* The approximation coefficients &#x03b1; can be obtained by inverting matrix '''A''' in Eq. <span id='cite-ZEqnNum188430'></span>[[#ZEqnNum188430|(5)]] but ill-conditioning can lead to numerical inaccuracies depending on the spatial distribution of points in the local cloud. While 2D approximations are normally insensitive to this problem and can be solved satisfactorily by matrix inversion, numerical misbehaviours related to cloud geometrical features are frequently seen in 3D. For such cases, the authors presented in [15] an alternative QR-based approach which has proven to improve accuracy and robustness.
-==2.1 The weighting function==
+===2.1 The weighting function===
 Different possibilities exist for defining the weighting function in Eq. <span id='cite-ZEqnNum595040'></span>[[#ZEqnNum595040|(3)]]. The typical choice in the FPM is a Gaussian function located at the star point of the cloud. Namely,