EVALUATION OF GMSH MESHING ALGORITHMS IN PREPARATION OF HIGH-RESOLUTION WIND SPEED SIMULATIONS IN URBAN AREAS

This paper aims to evaluate the applicability of automatically generated meshes produced within the gmsh meshing framework in preparation for high resolution wind speed and atmospheric pressure simulations for detailed urban environments. The process of creating a skeleton geometry for the meshing process based on level of detail (LOD) 2 CityGML data is described. Gmsh itself offers several approaches for the 2D and 3D meshing process respectively. The different algorithms are shortly introduced and preliminary rated in regard to the mesh quality that is to be expected, as they differ inversely in terms of robustness and resulting mesh element quality. A test area was chosen to simulate the turbulent flow of wind around an urban environment to evaluate the different mesh incarnations by means of a relative comparison of the residuals resulting from the finite-element computational fluid dynamics (CFD) calculations. The applied mesh cases are assessed regarding their convergence evolution as well as final residual values, showing that gmsh 2D and 3D algorithm combinations utilizing the Frontal meshing approach are the preferable choice for the kind of underlying geometry as used in the on hand experiments. * Corresponding author


INTRODUCTION 1.1 Motivation
In (Langheinrich et al., 2017) a new method was proposed to generate high-resolution wind speed and pressure maps from remote sensing data.These were derived by using a computational fluid dynamics (CFD) framework to simulate the turbulent flow of wind over a 3D model of an urban area in North Rhine-Westphalia (NRW), Germany.While the experiments in (Langheinrich et al., 2017) were successfully conducted, the underlying building geometry was still very coarse as it was generated from raster images.
Introducing inaccuracies in the Navier-Stokes calculation of turbulent flow, this low-quality modelling of the Earth's surface poses as a problem that has to be addressed to generate reliable and accurate datasets.With NRW introducing the OpenGeoData platform in 2017, users have unconstrained access to spatial base data of the state, including level of detail (LOD) 2 building models provided in CityGML format.This enables the generation of a higher quality city model than the model used in the initial CFD experiments.
While an increased accuracy of the underlying 3D Model in general can be considered to have an influence on the improvement of the CFD results, it also creates new challenges concerning the overall computation pipeline, in particular the generation of a volumetric mesh.The latter divides the defined atmospheric cube around the city model into individual cells, which enable the finite-element based numerical calculation of the Navier-Stokes equations.
Generating a well-formed mesh poses as one of the most time consuming and error-prone steps in CFD processing (Ali et al., 2016).Despite a high accuracy of the building model a mesh poorly derived from this 3D data can lead to largely inaccurate calculation results.The performance hereby strongly depends on how adequate the mesh is chosen and generated in regard to the particular simulation.This paper aims to evaluate the performance of different algorithms implemented in the open source 3D-mesh generator gmsh (Geuzaine & Remacle, 2009).

Comparable literature
In (Schubert et al, 2008) an unstructured meshing method is developed, aligning building footprints to mesh holes (not taking building height into account) to model hydrodynamics and evaluate the feasibility of the mesh in regards to urban flooding.(Gargallo-Peiró et al., 2016) present an automatic method to generate an urban surface geometry on basis of cadaster information and simplifying the topology by generalizing the build up to city blocks.
Reviews on general mesh generation methods and best practices concerning structured and unstructured meshes are given in (Ali et al. 2016) and (Frey, 2005).(Reiter, 2008) gives an in-depth introduction on the configuration of the CFD framework and evaluation of results concerning urban environments, utilizing the FLUENT CFD simulation software.
The on hand work differs from previous approaches in that it takes the detailed 3D geometry of every single building object into account, which enables the identification of small scale turbulent wind effects.

Geometry generation
In order to enable the proper mesh generation with gmsh an underlying model, the skeleton geometry has to be defined.The gmsh software itself features capabilities for this task in two forms.A general gmsh inherent approach to construct the geometry by defining primitive features like vertices, lines, surfaces and volumes as well as the CAD-like openCASCADE (SAS, 2008) framework offering constructive geometry features.Both approaches proved not feasible for an automatic geometry generation utilizing Python.
Because of inferior results produced by the gmsh inherent construction tools, not enabling proper mesh generation a manual approach was chosen by making use of the OpenSource CAD solution FreeCAD.
The available Python script to extract building geometries from a CityGML file was used to generate a basic 3D city model in terms of vertices and straight lines.This model was then imported into FreeCAD.Here solid objects of each individual building were produced in a stepwise manner by defining bounding polygon and outer shell objects.To produce the volume defining the atmospheric region around the buildings for the test area, a Boolean difference operator was applied to the building objects and a cube of the size of the distance field to generate a single solid atmospheric volume omitting the three-dimensional regions occupied by the modelled buildings (Figure 1).

Figure 1. Boolean solid of test area
This volume geometry is then re-imported to gmsh to define physical boundaries needed for the following CFD calculation like inlet, outlet and symmetry surfaces, wall surfaces and the physical atmospheric volume.In this state the geometry is then ready to undergo the meshing process.
One important practical note is that it was necessary to reduce the geometry to the origin as the high UTM values of the coordinates especially featured in the y-dimension lead to numerical problems in the mesh conversion step between the gmsh mesh and its OpenFoam representation.

Meshing algorithms
The gmsh framework offers several meshing algorithms for the two-and three-dimensional case respectively: For two-dimensional applications the Mesh Adapt (Geuzaine & Remacle, 2009), Delaunay ( George & Frey, 2000) and Frontal (Rebay, 1993) are considered proven and robust implementations.For the three-dimensional case this applies to the 3D adapted Delaunay utilizing an initial union of all the volumes in the model (Si, 2004), its improved form of New Delaunay and a 3D Frontal approach using J. Schoeberl's Netgen algorithm (Schoeberl, 1997).
The other algorithms are considered experimental implementations and, as will be seen in the presentation of the results, are not feasible for proper mesh construction in the on hand use case.
In regards to the meshing performance the gmsh documentation (Geuzaine and Remacle, 2017) preliminary states the behaviour of the stable algorithms as described in Table 1.A clear tendency of the performance of the meshing process in the direction of this preliminary characterization will be revealed in course of the presentation of the results.
Table 1.Preliminary meshing algorithm characterization with 1 representing the worst and 3 the best performance

Data
The basic geometry as well as geolocalization of the buildings in the test area was extracted from a LOD2 CityGML dataset provided by the German federal state of North Rhine-Westphalia via the Open.NRW geodata project and platform.
The test area comprises an area of 370 by 270 by 100 meters in the city of Gelsenkirchen in NRW, Germany, including nine distinctive building structures and located 361960E, 5716029N (lower left corner) to 362330E, 5716299N (upper right corner) UTM Zone 32 (Figure 2).Regarding these particular extents of the test area, the actual dimensions of the region were chosen according to: The blocking ratio not exceeding 3%, where the blocking ratio is the ratio between building vertical surfaces exposed to the wind and the surface of the simulation inlet (Baetke et al., 1990 andReiter, 2008).


The minimum dimensions of the simulated field is 5 * H (with H as the height of the highest building) upstream of the buildings, 5 * H on each side of the buildings, 5 * H above the buildings and 15 * height downstream (Hall, 1997 andReiter, 2008).
The particular area was also chosen to feature distinctive object build-ups in terms of lower and higher shape complexity as well as featuring open and blocked areas while covering a relatively small area to enable an acceptable processing time regarding the overall experiments.

Meshing in gmsh
Before the actual meshing process was initiated the corner vertices of the building objects were defined as refinement points with a factor of 0.1.This ensures that the mesh is constructed with a higher resolution approaching these regions of the geometry.A refinement is necessary to enable a more detailed numerical calculation near the wall regions and the spaces between the buildings in order to reveal small scale effects concerning turbulent wind dynamics like wind tunnel effects and creation of sub-pressure regions in areas of turbulent eddies in future experiments.The refined mesh can be seen in Figure 3.

Figure 3. Refined mesh regions around and between building objects
The meshing algorithms were then combined in all possible combinations of 2D and 3D implementations.Table 2 shows the results of the meshing processes in terms of applicability to the CFD calculations.
Table 2. Meshing results in terms of usability.Green = no errors, Red = Error in OpenFoam mesh conversion, Blue = Gmsh meshing error Applicability of the mesh was judged before the actual CFD process by converting the gmsh mesh to an OpenFoam representation via the framework tool gmshToFoam and the subsequent use of the tool checkMesh to quickly review the particular mesh for discrepancies.
In Table 2 the algorithm combinations marked in green are considered applicable for CFD calculations.Combinations in red were meshed without errors in gmsh but could not be transformed error free to an OpenFoam representation of the mesh.The R-Tree algorithm was not able to produce meshes from the provided geometry at all.Following these insights only the mesh cases marked with green color in Table 2 were applied in the next steps of the experiments The numbers given in the cells of Table 2 serve as case identifiers in the further course of this work.

CFD calculations using openFoam
As the emphasis of this work is on the evaluation of the mesh quality resulting from the different gmsh meshing algorithms, a less complex configuration for the CFD was chosen to keep simulation difficulties at an acceptable level.
Therefore a simplified version of the Navier-Stokes equations of turbulent flow is used by the following postulations:


The air, affected by the turbulent flow, can be treated as an incompressible fluid in the on hand simulation as the simulated wind speed is at 30 m/s which adds up to a low Mach-Number of 0.087 (Tritton, 2012). A constant value is considered for the boundary conditions of wind speed and air pressure resulting in a steady-state realization of the equations.
Further a standard RANS (Reynolds-averaged Navier-Stokes) turbulence model is chosen in form of the Lam and Bremhorst low-Reynolds number k-epsilon turbulence model (Lam & Bremhorst, 1981).
The initial boundary conditions are set to 30 m/s for wind speed and the pressure to 1013.25 hPa.The simulation is iterated over a period of 400 seconds with an iteration step size of 1 second.

Evaluation approach
The evaluation of the particular mesh quality in terms of numerical CFD calculations poses as a not too easy task.In order to judge the individual feasibility the mesh quality has to be measured by the results of the actual CFD process.To do so by means of absolute values, the results should be compared to reference data in form of in-situ or wind channel measurements, both of which are not available for the area of interest.
Therefore a relative approach was chosen to evaluate the performance of each individual mesh.The OpenFoam framework iteratively outputs the residuals of the target values (wind speeds in x-, y-and z-dimension as well as pressure) and the parameters of the k-ε turbulence model of the Navier-Stokes equations of turbulent flow.For CFD, RMS residual levels of 10 -4 are considered to be loosely converged, levels of 10 -5 are considered to be well converged and levels of 10 -6 are considered to be tightly converged (Kuron, 2018).
By these criteria the results of the computational fluid dynamics calculations will be compared to derive the individual relative mesh quality.

RESULTS
Figures 5-10 show the residual convergence plots of the CFD simulation over the period of 400 iteration steps for the applied meshes (for case number reference see Table 2).It can clearly be seen that mesh cases 1 and 6 feature the highest final residual values as well as the highest continuous value instability seen by the strong and frequent fluctuations of the graph.Cases 12 and 13 on the other hand throughout exhibit a smooth convergence curve.By contemplating Table 2 and 3 this conforms to the preliminary characterization of the algorithms showing that the 2D Mesh Adapt algorithms generally generates less robust meshes while solutions utilizing the Frontal approach in two and three dimensions offer the highest robustness.
To further evaluate the overall feasibility of the individual meshing algorithm combinations Table 3 lists the final residuals for every variable.Judging from the overall number of good residual values per case, it can be seen that the meshes of highest relative quality are generated by combinations using the 3D frontal approach (case 11, 12 and 13) which also explains why case 11 which is meshed in 2D with the inferior Mesh Adapt algorithm shows a higher robustness than the two meshes from case 1 and 6 as described above.The comparison leads to the conclusion that for the on hand experiments, meshing case 12 can be considered the best quality mesh in terms of CFD applicability as it accumulates the highest count of top 3 low residual values.This corresponds to the algorithm combination of 2D Delaunay and 3D Frontal meshing.The resulting mesh is shown in Figure 4.
Figure 4. Crinkled cross-section of the case 12 mesh overlayed with the converged wind field at a height 2m above ground

CONCLUSION
While the presented analyses shows that for the conducted experiment the 2D Delaunay -3D Frontal approach proves to be the best performing meshing algorithm combination in gmsh for the non-manual mesh preparation of an urban geometry for CFD applications, further experiments are planned in preparation for the application described in the beginning of this paper.
Analysis has shown that for the general case unstructured meshes as produced by all of the gmsh algorithms are inferior in their quality to structured meshes (Ali et al., 2016).Therefore the on hand results have to be set into relation with the mesh quality of structured approaches for an urban environment geometry model in future experiments.
Further, it is expected that a more complex CFD simulation configuration for the Navier-Stokes turbulent flow of air is necessary to model critical turbulent behaviour during a wind event in a city.In future it has to be tested how and if meshes resulting from the presented methods hold for these more demanding simulation setups in terms of feasibility.

Figure 2 .
Figure 2. Satellite image of the test region (Source: Google Earth)

Figure 5 .
Figure 5. Residuals of wind speed in x-dimension in m/s

Table 3 .
Final residuals of CFD process per mesh case.Green = best result in column, Yellow = 2 nd best result and Orange = 3 rd best result