When we were approached by the folks from FRIENDSHIP Systems to work with them on integrating Pointwise into their optimization framework it was easy to agree because the project touched on many facets of CFD that interest us: collaboration, automatic meshing, OpenFOAM, and optimization.

In this first phase of what should be a multi-phase effort, FRIENDSHIP Systems’ Carsten Fuetterer was able to demonstrate successfully how to automate mesh generation with Pointwise, integrate Pointwise into FRIENDSHIP-Framework (FFW) and use the two (along with other tools) to optimize in 2D a vertical axis wind turbine. The work described below and most of the written words are Carsten’s.

The Optimization Process

Geometry is created in FFW based on design variables and the optimization algorithm. Pointwise uses that geometry to create a high quality mesh which is used as input for an OpenFOAM CFD calculation. Gnuplot is used to create plots from the OpenFOAM results, which are then imported into FFW to evaluate the objective function for the optimization which in this case is the coefficient of performance, COP.

Each component of this process is described in the following sections.

Geometry Creation in FFW

As a first step, a fully parametric 3D model of the VAWT was created in FFW, where the airfoil shape and the shape of the 3D model can be varied using design variables such as pitch angle, camber, camber position, chord length, thickness, radius, and twist.

The twist design variable can change the VAWT from a straight blade type to a helix type. Also, all design variables can be varied in all directions. For example, it’s possible to change the chord length from the midspan to the beginning and to the end of the blades. The pitch angle is defined as zero when the leading edge is parallel to the dotted circle as shown below. The export format for the geometry is IGES or STL.

The Objective Function

The performance of wind turbines is often evaluated with the coefficient of performance, COP. This coefficient describes the ratio of the power that is extracted by the turbine to the power available in the wind.

In the equation shown above, τ is the average torque delivered by the wind turbine, ω is the rotational speed, ρ the density of the air, h is the height of the turbine (in the 2D case it is the height of one cell), d is the diameter of the turbine and Vw the design wind velocity.

The theoretical maximum for COP can be expressed by the Betz limit which is around 59%. Because the torque changes during a single rotation, the average torque of the last rotation is calculated in FFW and used as input for the COP evaluation.

Mesh Generation with Pointwise

The mesh of the 2D geometry is created with Pointwise. This program allows the generation of high quality meshes and is able to run in batch mode with the help of Glyph scripts. For beginners, the creation of the scripts is easy because it can be recorded during the set up. For more experienced users the scripts can be written directly, to get more compact and robust scripts.

Best practices suggest that a hybrid mesh would be the preferred choice. This consists of a structured mesh around the airfoils to resolve the boundary layer, which is needed for the selected turbulence model, and an unstructured mesh in the far-field.

The boundary layer is resolved with 25 cells where the first cells are 0.0002 m high to achieve y+ values of around 1-5. The growth factor is defined to be 1.1. The spacings on the connectors of the pressure and suction side are 0.003. Due to high gradients, the spacing is set to 0.0012 at the leading edge, and 0.0015 at the trailing edge. This technique ensures that fine resolution is achieved where it’s needed, while decreasing the overall cell count. Furthermore, it results in a smooth transition between the structured and unstructured cells. The trailing edge is resolved with 6 cells. The resulting mesh around the airfoil is shown in the images below.

The rest of the domain is meshed with unstructured cells, with a smooth transition from the fine cells at the airfoils to the more coarse cells in the outer domain. The mesh of the whole domain is shown below.

The mesh at the interface between the stationary and rotating region is shown below. Since a sliding mesh approach is used for the CFD calculation, the nodes at each side of the interface have to be separate. Otherwise it would lead to mesh deformation during the rotation.

The near-wall mesh consists of 59,518 cells and the outer mesh consists of 31,050 cells, which leads to an overall mesh size of 90,568 cells. The mesh generation is very fast taking less the one minute.

Another big advantage of the mesh generation with Pointwise is that the mesh can be directly exported to an OpenFOAM file format. All boundary conditions, types and names can be set in Pointwise. Some small modifications of the OpenFOAM text files are necessary, but this can easily be done automatically.

CFD with OpenFOAM

The CFD calculation was done with the open source toolbox OpenFOAM 2.2.0. The computations are transient with the sliding mesh approach in the rotating region. Due to turbulence and separation effects at the airfoils, the k-omega SST turbulence model was selected. This model combines a good description of the fluid effects near the airfoils as well in the far field. The solver was pimpleDyMFoam with a fixed time step of 0.0001 s which leads to CFL numbers of 0.04 -7. The interface is from the cyclicAMI type.

The forces and torques are recorded during the calculation, so that the torque can be extracted after the simulation. Gnuplot is used to plot the torque and residuals. The simulation of two VAWT revolutions took around 2 hours and 15 minutes. Typical flowfield images are shown for pressure followed by velocity.

Optimization with the FRIENDSHIP-Framework

FFW is able to use the scripts of Pointwise and OpenFOAM in order to run an automatic variation and optimization study. Therefore two so-called “Software Connectors” are used (not to be confused with Pointwise’s connectors), with one for each program. The first software connector exports the IGES geometry file and runs the Pointwise script. The result files are the exported Polymesh files for OpenFOAM. These files are then linked to the second Software Connector for OpenFOAM.

The Software Connector for OpenFOAM is shown below. As input files, it contains the linked file from the Pointwise Software Connector, the text files used for OpenFOAM, and the script called “Allrun”. This shell script triggers the commands for the different OpenFOAM tools. Furthermore the script runs Gnuplot to create plots when the simulation is done. The result files are the two plots from Gnuplot, the case.foam and the forces.dat files. In the Results Files part of the Software Connector, the last values of the forces.dat file can be extracted and used as an objective. Due to the fact that the torque has to be averaged over the last rotation, the actual objective values are extracted from the table forces.dat from the Result Files part. A parameter can calculate the number of time steps per rotation, which can be used to average the last rows for one revolution. This method is only possible when the time step is fixed in OpenFOAM.

From the averaged torque, the coefficient of performance can be calculated, which is used as a single objective function for the optimization algorithm.

For the optimization, many strategies are possible. One strategy often used with FFW is a Design of Experiments (Sobol) which scans the design space for global minima. This DoE is often followed by a tangent search (T-Search), which is a gradient based algorithm for local minima.

In this study, a genetic algorithm NSGA2 is used. Generally, it is similar to the combined approach of DoE and T-Search as mentioned above. It creates populations which are based on the best candidate designs from the previous generation. But also there is the mutation part, where design variables are combined randomly. This mutation makes sure that the algorithm ultimately finds the global minimum. A disadvantage is that the time to converge can take longer, due to the mutations. And a compromise has to be found between scanning the design space and the time to converge.

Optimization Results

The optimization was done with a wind velocity of 10 m/s, a tip speed ratio of 0.5 and a radius of 0.6 m. The following design variables were chosen, Pitch angle, Camber, Camber position, Thickness, Chord length.

As a first step a parameter screening of the pitch angle of the blades was executed, to reduce the size of the design space. The results are shown below.

The results of the optimization run with the NSGA2 algorithm are shown below. As can be seen, the algorithm converges quite well by about 90 runs. If more designs were run, additional improvements might be realized, since a small increase of the COP can be seen and the two best designs are within the last 10 runs. But as a first guess the improvement with more runs would probably be below 1% or 2%.

The best design has a COP of 26%. From the initial 25 designs it can be seen that the algorithm scans the design space and how rapidly it can find better designs. After the initial 25 designs, the rate of improvement is more modest and increases by only 3% to the best design after 81 runs. The few designs (6) which have a COP less than 22% are due to the mutation factors and are necessary to scan the design space, to avoid local minima and ultimately result in faster convergence. The fact that only six designs are below 22% after 25 runs, could lead to the assumption, that the mutation factor is too small to scan the design space effectively. Vice versa the crossover factor could be too high, so that the algorithm just crosses the best designs thus far.

Finally, the normalized variables are shown over the runs in comparison to the normalized COP. This graph shows detailed information of how the optimization algorithm worked and if the bounds for the design variables are set correctly. For camber, camber position, pitch and chord length, the convergence is quite fast. This can be directly correlated to the fast convergence of the COP. After 40 runs the thickness seemed to be converged as well, but after that it increases suddenly to the maximum value. It can be seen, that this increase leads to an increase in COP. This is a good example of how the algorithm leaves a local minimum due to mutation factors. The fact, that the maximum thickness is at a value of one leads to the assumption, that the upper bound for the profile thickness is too low, so that a further improvement could be reached when the thickness is increased further. This upper bound is also a limiting factor for the chord length. But both, the thickness and the chord length have direct impact on the weight of the wind turbine, which will results in much more resistance during its rotation. This effect was not studied.

Compared to larger VAWT, which achieve COP values of 35-40%, this setup with low tip speed ratio and small radius didn’t perform as well. This leads to the assumption, that to improve the COP, the tip speed ratio (TSR) has to be increased. But with the specified TSR there are already approximately 1.3 full rotations per second, which is already quite high for real world applications. Another way to increase the TSR could be to increase the radius of the wind turbine.

Conclusions

The study created a prototype of an optimization process for vertical axis wind turbines centered on FRIENDSHIP-Framework and using the commercial mesh generator Pointwise along with open source software OpenFOAM and Gnuplot. The next step will be to extend the analysis and optimization to 3D.

We thank Carsten for all his work on this project.