As we said past week, I kept on investigating about the influencing parameters in the SEM simulation. I also checked my post-processing MATLAB script, and I inserted a post-processing in the fortran code as well. Here you are my results.

Steps in the results Prova dimmi dove la metti

We noticed last week the presence of some suspicious steps in the results. As we can easily see in this picture in fact there are several steps even after many seconds dozens (that means some thousands iterations). This phenomenon is much stronger as the eddies density decreases.

A furter investigation was required then, particularly taking a look at the instant velocity, to check if there is any big jump between two "very" different values. This analysis was completed and gave as a result the following graph. This graphs enphatize the area where one of these steps occurs (precisely is the jump you can see in the dens-019.90 line around the 40th seconds), and it easy to see the relation between the velocity and the step.


Question 1: what can cause such a phenonemon?

Such a behaviour may happen when two or more eddies, with similar spatial intensity, passes very close to the point considerated. In such a case, their inputs sum themselves and we see such a jump in the second-order averaged speed.

Question 2: after zooming the instant velocity shows a non uniform fluctuation. Is it right?

This second question is not directly related to our previous one, but came out from the instant velocity graph. We see infact a strange behaviour. The instant velocity shows a straight behaviour with fluctuation applied over it. These fluctuation seems "not to be enough". So we need to quantificate this enought. Past week I tried to find a relation between the eddies number (and so their density) and the convergence of the results. Mainly (as you can see later) I demonstrated that to get a full converged simulation you need to set-up a very big time-step number, that implies a very long CPU time requirement. So this method failed. This week I show a different method, based on the PDF of the instant velocity. I made a MATLAB script (available as a attached file) that allow you to create the PDF of the instant velocity. For this case the output is the following. PDFcomparing.jpg

The PDF shows exactly this non physical behaviour. We can easily see that there is a peak right around the 10 value (blue line), which can be explained with the hight number of point that have such a velocity. As we thought, such a density is not enough for our eddy box.

I plotted in red the PDF of a case with 1000 eddies. In this case we have a totally different behaviour, without any peak that shows the best velocity distribution.

This, then, seems to be another productive way to know if the eddies density is hight enough.

Convergence of the SEM

Another problem stressed past week was the non-convergence of some results, or, in some cases, the convergence to differents value. This week I ran some very long simulations (100.000 time steps --> 500 s --> 8 min) which showed that our simulation had not already converged. In fact, after some dozens-thousands time steps we really reach our expected results. The main problem is that it took from 24 to 36 hours to complete such simulations!! Anyway, here you are the results.


These three pictures shows theconvergence of a 100,000 time stes simulation. The simulation was set in this way:

N = 1000

x = 2 * %$ \sigma $%

dt = 0.005 s

We can see that all the parameters converge. Obviously the higher is the averaging order the longer it takes.

< Convergence of %$ <F_u> $% terms

^ Convergence of the %$<S_u>$% terms

\ Convergence of the %$<u^2>$% terms

Box Volume

After all the discussion, we see that one of the most important parameters so far is the eddies density. To increase such a parameter we can both try to increase the eddy number (not totally convenient because of increasing CPU time requirements), or reduce the box volume, which seems to be the easiest way to reach a good results. Obiously there are some limitations:

  • y and z dimensions are limited by the grid itself (that in this method is given by as a kind of external condition)
  • x dimension is limithed by the eddies %$ \sigma $%

In fact, in our situation (a 2D grid), I think that when the posterior part of the eddy has passed the grid it, obviously, cannot have any influeces on our results. So I investigated over this parameter, to see if there are substantial results difference between choosing a 4 %$ \ sigma $% x-lenght or a 2 %$ \sigma $% x-lenght.


This first graphs (average Su for x = 2 %$\sigma$% in red and x = 4 %$\sigma$% in blue) shost that there is no substantial difference beetween the behaviour of the two simulation. As I were expecting, chose a lenght is long as twice %$ \sigma $% does not change our results (but we have to be carefully that every vortex has already passed the grid befor it is re-generated). I got similar results for all the other not-plotted parameters


This second simulation show one main difference. As we go back to the theory we see that:

%$ F_{u_i} = \frac{<u_{i}^{'4}>}{{<u_{i}^{'2}>}^{2}} = \frac {1}{N^2 R_{ii}^2}<{(\sum_{k=1}^N X_{i}^{(k)})}^4> = 3 + \frac{1}{N} (4 F_{f}^3 F_\epsilon \frac {V_B}{\sigma^3} - 3) $%

From this equation we can explain the difference between the two simulation. In fact the asymptotic value is related to the volume of the box, that in our cases varies. We have a steady-result closer to 3, that is the asymptotic value reached when N --> %$\infty$%

PDFcomparisonX.jpg This last graph is the PDF of the two signals. Again we can see a positive influences of the volume decreasing. In fact the shape indicates a higher concentration towards the mean value for the x = 4 %$\sigma$% simulation.

All these results are very important because allow us to decrease the volume (and so increase the density), with the same CPU time. A further reduction in the x dimension, unluckily, involve a wrong implementation of the SEM method, which may appear as increasing mean-velocity and so on.


My personal conclusion at the moment are:

  1. The SEM method is a "slowly" converging method. This mainly because of the hight eddies density required which imply a hight number of eddies. As a consequence the CPU time required by the method is very long. To recude it we can:
    • reducing the box volume (particularly in the x direction), so to decrease the eddies number required
    • chosing the lowest possible number of eddies (using the PDF function as shown)
    • accepting a not-full-converged SEM solution
  2. As showed by Nicolas Jarrain, the final velocity field of the SEM method may be close to the LES one, but there is still a gap. He suggested a Divergence Free scheme for the SEM method, which could bring to some better results (I am still implementing it). The real question, in my opinion, is: does it really worth a while to get a fully converged solution in the SEM method, if it is not really what we are looking for?

Extra contents

The following paragraphs are added after the 16/11/'09 meeting. There is a part to validate the post processing script and a part related to the comparison between theoretical and numerical results.

Post-processing validation

Discussion 1 - What do the calculated quantities represent?

As we discussed during the meeting, to exclude any error in the preprocessing i run sum text with it. Here you are some results.(I have also attached the post-processing matlab code, if you want to check it!)


This first image shows three different normalized velocity:

  • %$ <u^{2}> $% in the blue line
  • %$ <u`^{2}> $% in the red line
  • %$ <u> $% in the green line

As you can see from my post-processing file the ways I calculate them are the followings:

%$ <u(i)> = [u(i) + <u(i-1)> * (i-1)] / i $%

%$ <u^2(i)> = [u^2(i) + <u^2(i-1)> * (i-1)] / i $%

%$ <u`^2(i)> = <u^2(i)> - <u(i)>^2 $%

Discussion 2 - Output of a already known signal


As aimput signal I used the following one:

%$ Signal = 9 + 2 * rand () $%

The generated signal is a equal probability [9,11] signal, with a average of 10. I passed this signal to my post-processor and, as expected, the output value of the first two order averagin are exactly as expected.

I tried the to calculated the 'hypothetical' Reynolds stresses:

%$ <u'u'> = <u^2> - <u>^2 $%

and again the numerical results are very confident with the theoretical ones.

Theoretical vs numerical results

The last discussion is about the comparison between theoretical and numerical results. In fact, mainly for the highest averaged order, we wanted to stress if there is a confident or not confident behaviour.

In the following graph (I plotted it as bigger as possible to be more understandable), you can see that for this simulation thenumerical results agrees with a good approximation to the theoretical ones.


I plot here instead the same variables plotted for a simulation with twice the eddies number. From the theory we know that we should have a asympthotic behaviour towards 3 as the eddies number increases. This behaviour is confirmed, as we can easily see that in this second case the asymptoticvalue is closer to 3 than the previous one. Further, this is another proof that the numerical results converge to the theoretical ons (even if they require thousands time-steps)


Extra content II

As discussed past week, I added the comparison between my and A. Revell results.

uu_rugg_vs_alistair.jpg Fu_rugg_vs_alistair.jpg

Basically, there are no many differences between the lower order (the two line are exactly one over the other), instead for the highests order the only differences are at the initial part (may be caused by numerical approximation since one method is Matlab-based, and the other one is Fortran-based). But even now, after about 1500 iteration the difference between these two metods is negligible.

Current Tags:
create new tag
, view all tags
Topic attachments
I Attachment Action Size Date Who Comment
elsem matlab_postprocessing.m manage 15.2 K 2009-11-17 - 12:21 UnknownUser Used matlab postprocessing file
elsem pdfplot.m manage 1.4 K 2009-11-16 - 13:50 UnknownUser Matlab script that allow you to print the PDF of a incoming signal
Topic revision: r4 - 2009-12-01 - 09:18:15 - RuggeroPoletto
Main Web
26 Aug 2019


Manchester CfdTm

Ongoing Projects


Previous Projects


Useful Links:

User Directory
Photo Wall
Upcoming Events
Add Event

Computational Fluid Dynamics and Turbulence Mechanics
@ the University of Manchester
Copyright © by the contributing authors. Unless noted otherwise, all material on this web site is the property of the contributing authors.