CFD Online Discussion Forums

CFD Online Discussion Forums (
-   OpenFOAM Running, Solving & CFD (
-   -   Lund Recycled Method for LES (flat plate) (

Joachim September 23, 2013 13:27

Lund Recycled Method for LES (flat plate)
5 Attachment(s)
Hey all!

I recently implemented Lund's boundary condition in OpenFOAM.
The BC seems to work just fine, but for some reason, my velocity profile is not that great.

Initially, I used exactly the parameters described in the original paper (same mesh, (y+=1.2, ...) second order accuracy), but I got pretty bad results (see Fig. 1).

Then, I decided to refine the mesh (basically refine the mesh in each direction by a factor 2 -> y+=0.5, ...), and obtained the results shown in Fig. 2.

For these simulations, I use the simple Smagorinsky model, with and without Van Driest damping functions to evaluate its impact on the solution (pretty significant!). I also tested the local dynamic Smagorinsky model, but the results were not that great (skin friction underestimated).

As seen in Fig. 3 and 4, it seems that the boundary condition has been correctly implemented. We can see the self-similar solution between inlet and recycled plane for the inner and outer layers. However, the simulation does not converge towards the correct profile (in Fig. 1 and 2)!

The fine simulation is not that bad, but several people (Jewkes, Lund, etc) obtained very good results on the exact same case on a way coarser grid (the first one)!

Has anyone already faced similar issues with LES in OpenFOAM? Where could this come from?

Thank you very much for your time!


cfdonline2mohsen September 23, 2013 14:23

1 Attachment(s)
Dear Joachim

Hello and so nice talking to you again.
your 2nd Fig. from left seems perfect to me. your results are exactly coincident with the Lund's results up to Y+=50 and after that there is a little discrepancy but they have the same trend.

There is something that seems very interesting for me. when using The simple Smagorinsky model with Van Driest damping functions i.e,:
Attachment 25567
your results are completely coincident up to Y+=50 but after Y+>60 The simple Smagorinsky model with cubeRootDelta are completely coincident with Lund results.
So if we combine these two deltas, your desired Results will be obtained !!!
have a look at the source of the vanDriestDelta.C:

The purpose of the Van Driest damping function is to correct the nusgs near the wall but why it will affect the results for 60<Y+<1000??

So may be this is the problem not your mesh!

have you tried the modified version of dynamic Smagorinsky model by Alberto Passalacqua? what about dynOneEqEddy which is dynamic and also solve an extra transport equation for ksgs?
By the way, your results are for which version of Lund's method? displacement thickness or momentum thickness?

Joachim September 23, 2013 16:26

2 Attachment(s)
Hey Mohsen!

Yep, something is definitely strange with this damping function. :)
It actually makes sense that I must not use it in the viscous sublayer, because my mesh is so fine (y+ = 0.5 with hyperbolic tangent stretching to increase the grid density close to the wall). However, it might be that my mesh is too coarse in the log/outer layers (see attached figures).

I have used Alberto's dynamic model, but it did not give very good results on the base grid. I did not try it on my fine mesh though...

I can't use the dynamic one equation model at this point, because I would need to recycle the turbulent energy. I don't know which self-similar profile I should assume for it. It might be interesting to investigate this, but I really don't have time right now. :)
I will use the one equation model for the actual simulation.

In these simulations, I imposed the boundary layer thickness, as was done in the original paper. I also developed a version of the BC based on the golden section method to impose the momentum thickness, but the results were really similar.

I think I will just extract the plane from my base case and see if it gives good results for my "main" simulation. :D

cfdonline2mohsen September 23, 2013 17:07

I have a question about your sample plane:
How did you define 2 cyclic planes inside the domain (interior patches) with very little distance (I think 0.001)? did you use CreateBaffles utility?


I think I will just extract the plane from my base case and see if it gives good results for my "main" simulation
is your recycling section seperate from your main domain? How did you run 2 simultaneous cases at the same time and use your recycling outlet for your main inlet at each timestep?
in the Eugene de Villers's PhD dissertation, the recycling section is attached to the main domain and there is a single domain. which one is better?

Joachim September 23, 2013 17:17

I actually used ICEM CFD to generate my mesh. For some reasons, blockMesh was unhappy with my little gap.

Indeed, it is possible to attach the two domains. One good thing with this approach is that you don't have to worry about the outlet contaminating your inflow simulation. In addition, you don't have the hassle of having to extract your plane and re-inject it in the next simulation.

I did not use this method for two reasons:

1. My final objective is to verify my implementation of a new HRLES model. As mentioned previously, I cannot recycle the turbulent kinetic energy, so I cannot use this model for my inflow simulation.

2. I will probably run quite a few simulations later on, when improving my turbulence model. I don't want to run the inflow simulation over and over again. One time is enough! :D

cfdonline2mohsen September 23, 2013 17:45

Thank you so much for your helpful comments.

Indeed, it is possible to attach the two domains. One good thing with this approach is that you don't have to worry about the outlet contaminating your inflow simulation.
Yeah, I agree with you.

I don't want to run the inflow simulation over and over again. One time is enough!
What do you mean? I think inlet must be time dependent in LES and one of the objectives of this method is to create time-dependent turbulent structures at each time step. So you are going to use the outlet of this method (a fixed value) for you inlet of the main calculation at all time steps??

Joachim September 23, 2013 17:47

Nop, I am extracting the velocity, turbulent kinetic energy and pressure from a given plane (in the middle roughly) at each time step. I am saving that somewhere, and then I use it in my main simulation using the timeVaryingIDontRememberTheName boundary condition. :)

cfdonline2mohsen September 23, 2013 18:27

at each time step!
No Way!:eek::eek:

I am saving that somewhere
Huge amount of data!
I have worked with timeVaryingMappedFixedValue before.
It is really troublesome to correct the data at each time step for the required format of timeVaryingMappedFixedValue.

pardon me for asking a lot of questions but this topic is very interesting for me.
And my final question for today (I promise:p):
so you mean that first you will run the driver section for generating inflow and save the results from a given plane at each time step.
after finishing your initial runs for the driver section you will start your main domain run by using these saved data at each time step. am I right?

Joachim September 23, 2013 18:30

Yep, just like it is done in the pitzDaily tutorial. The amount of data is not that huge actually. You do not need to save the whole flow field, only one plane. :)

cfd.with.openfoam November 27, 2013 15:19


LES Inflow - In this area I am trying to wrap my head around the terms recycling and rescaling. Following is my current understanding -

1. Recycling - It is basically helpful when we have a fully developed flow e.g in a channel or a pipe. Is that correct? This is where directMapped comes into picture.

2. Rescaling - It seems to be helpful in a developing BL e.g BL over a flat plate. Here we can maintain a specific BL thickness at the inlet by rescaling the flow field data from a downstream plane by e.g using the law of the wall in the inner region and the defect law in the outer region (from TS Lund 1998).

My questions - Do we have a utility in OpenFOAM that can do this rescaling and maintain a specific BL thickness at the Inlet? You mentioned that you have developed a BC that is based on the paper from TS Lund. Would you be able to share that with me?

In this rescaling technique it seems that I would still need a mean profile for different variables at the inlet of my domain. After that I would need to perturb the domain e.g using the perturbU utility. After following these steps I should expect a self sustaining turbulent BL. Is that correct?

There are some more questions but the message is already very long.

Thanks very much.

Joachim November 27, 2013 16:36


As you said, when you recycle the flow, you do not modify it. It is like a "copy-paste" operation. When you re-scale it, you basically "copy-modify-paste". In the case of Lund's method, you assume some kind of self-similar turbulent profile (law of the wall close to the surface and law of the wake in the outer region).

The boundary condition hasn't been implemented in the official version of OpenFOAM. When I did it a few months ago, I developed three versions of the code:

1. a first one that imposes the boundary layer thickness at the inlet (the original model of Lund)
2. a second one that imposes the boundary layer momentum thickness at the inlet (suggested by Lund. This has to be done iteratively. Basically, you impose different boundary layer thicknesses and everytime, you check the resulting momentum thickness. I used a Golden-section approach to converge quickly).
3. a final algorithm that imposes the boundary layer displacement thickness. This was suggested by Jewkes in 2011 or 2012 (I don't remember exactly).

To be honest, I obtained very similar results with the three approaches. It basically depends on your application. I like to impose the momentum thickness, which was a bit more elegant than the other ones I thought.

You are right. You will also need to perturb the flow if you want to obtain a turbulent profile. You don't need anything fancy, just use Spalding's law or something similar (if you are very picky, I guess you could couple it with Cole's law in the wake, using the weighting function described by Lund).

Send me your private address in an email and I'll send you the codes if you want. It was really frustrating because I spent a lot of time on this, but it still did not give results as good as the ones obtained by Lund or Jewkes. I never found out why! (maybe you will ;) ). Still, it works, and will give you a pretty decent profile for your main simulation.

Good luck!


cfd.with.openfoam November 27, 2013 17:47

Dear Joachim,

Thank you for your prompt response. Can you please send the code/s to this address -

Are you currently carrying out your main simulations, using the inflow data generated by precursor simulations, employing Lund's rescaling method, that you developed? If yes, then do you use an inbuilt BC for mapping the saved data to your main simulation inlet? I have seen that there are BC's named timeVaryingProfile (cant remember the exact name). Can you please shed some light on the complete process that you follow?

Awaiting your response
Thanks very much

cfd.with.openfoam November 27, 2013 19:07

Dear Joachim,

I just looked at another post of yours (no. 7 in the link below) where you thank Perry Johnson.


I also came across a couple of papers from him (links below) where he talks about the rescaling method. My question is - were you able to discuss the Lund BC problem with him as it seems that he might be able to take a quick look at the code and point out the problem.


The issue mentioned in link no. 1 is also quite interesting as to how one can couple swak4Foam (or groovyBC) with directMapped (or turbulentInlet or Lund BC). How do you manage to handle such an issue? One thing that I can think of is just to initialise the whole domain with a profile that we would like to have at the inlet, perturb it and then Lund BC will take care of the rest? Any comments?

Thanks very much.

Thomas_Lloyd December 9, 2013 06:14

Lund recycle technique
Hi Joachim,
I am interested in using the Lund recycle technique for examining tidal channel turbulence profiles, and comparing the results to a synthetic inflow turbulence generator that I have been using. I just finished my PhD which was investigating tidal turbines using LES.
Maybe you could share the code with me and I could help check it? If you are happy to, please send it to

jnilsson January 9, 2014 06:23

Dear Joachim,

I'm interested in your OpenFOAM implementation of the Lund BC. I started to implement my own version of it at what appears to be, roughly speaking, the same time as you. I have implemented a rather simplified version of the BC, but had then to turn my attention to other things. My version seems to help my main simulation, but it is not as physically accurate as a full/proper implementation that you have made. I was now considering going back to continue to work on my implementation of the BC when I saw your post. I'm a PhD student and use the BC for my work/thesis. If you would like to share your implementation with me I would be very happy. Hopefully I can contribute as well. I have sent you my e-mail address over pm.

Thank you in advance

lixx January 9, 2014 09:57

Dear Joachim,

I am also working to implement Lund's method and its variant, the method by Kataoka. Could you please share me your implementation so that I can see if I can learn from it and, if possible, help contribute to it? My email address is

Many thanks.


cfd.with.openfoam January 18, 2014 16:42

Hello Everyone,

If anybody needs Joachim's code then please send me an email at

I wanted to upload the code here but I don't have his permission yet. Hopefully he will be OK with me distributing the code via email :)


flames January 26, 2014 20:12

Dear All,

I am also interested in the Lund's method and am reading his JCP paper in 1998. As you know, in that appendix, he explained the random fluctuation inflow generation method. I did not understand how the amplitude tensor a_{ij} in Eq. (27) come from. Does anybody know how they derived? GIve me some hints?


Note (Jan. 27, 2014):

The above are the elements of the lower triangular matrix of Reynolds stress tensor after Cholesky decomposition. very standard. But another question arises: why a_ij here can be used to represent the amplitude of the random signals and then the constructed velocity has the desired second order one point statistics? Anyone knows? any comments are very welcome!

The link for Cholesky decomposition:

ChrisWe April 20, 2015 02:22

Hello FOAMers,

I'm also interested in Joachim's implementation of Lund's BC in OpenFOAM. Is it still possible to get the code from anybody who used it? I want to do some DNS with it and it would help me so much to get a running code for this. My email is:

If it is not possible anymore to get this code, please let me know...

Thank you in advance,

PS: Thank you so much, this helped me tremendously!

ChrisWe April 24, 2015 05:19

Dear all,

when I tried to run Joachims BC (scaledMappedVelocityFixedValueDelta), I got some ERRORS, that are identical to my problem described here

To run the scaledMappedVelocityFixedValueDelta BC, I copied the:
- Make/files
- Make/options
- scaledMappedVelocityFixedValueFvPatchField.C
- scaledMappedVelocityFixedValueFvPatchField.H

into the following folder:

Compiling the new BC with wclean and wmake libso did not show any errors and there is a created in:

Changing the boundary field for a simple test case in the 0/U file to:

        type                scaledMappedVelocityFixedValueDelta;
        value              uniform (10 0 0);
        interpolationScheme cell;
        setAverage          true;
        average            (10 0 0);

and inserting in the controlDict (Make/files: ...LIB = $(FOAM_USER_LIBBIN)/libScaledMappedVelocityBCDelta)

libs ("")
I ran my simulation and got a FATAL IO ERROR

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Create time
Create mesh for time = 0
Reading field p
Reading field U

Unknown patchField type scaledMappedVelocityFixedValueDelta for patch type genericPatch

Valid patchField types are :

When I tried to do the same with the "mapped" BC (copying to my personal folder and try to run it without changes), I got the identical error. When I tried the same with the "parabolicInlet" BC (no generic BC), all worked fine.

Is there anybody, who knows, why the scaledMappedVelocityFixedValueDelta BC could not be detected by OpenFOAM2.3.1? Maybe someone has a running testcase that he/she can share with me?

Which OpenFOAM version did you use to run Lunds BC?

Thank you in advance,

All times are GMT -4. The time now is 21:31.