Documentation for Externals

[bstoeve]

Hi,

we are new to the neko code. It is clear, it is still heavily under development. Nevertheless, we are trying to test it for our applications. Now, we discovered several points, which made it difficult to enter the usage of the code with the exiting guides.
We would be willing to contribute on this point, however, we are a bit shy on just adding stuff, without contacting the main contributors. Who can we contact or could contact us, to get us a bit more into the project?

With kind regards

[timofeymukha]

Hi! Glad to hear you are trying out the code. We have a public chat channel on this Zulip:
https://nekodev.zulipchat.com/
Perhaps, it is a good place to discuss things, or we could just continue here in the issue. Very interesting to hear what kind of documentation is missing from the perspective of an external user.

[timofeymukha]

Dear @bstoeve , I've converted your issue to a "Discussion". We can continue here or on zulip.

[bstoeve]

Thanx a lot! Maybe for now, I would keep it here, since it can be read by externals, making it easier for others to approach the code.

[bstoeve]

Now, to our background: We have been using OpenFoam for a long time. I personally, also have some experience with nektar, but this is some time ago. We haven't used nek5000.
The installation works fine, if you follow the instruction. So this was fine.
The first issue, we ran into, was, that we had no idea, that the files used are nek5000 files in terms of the format. Also, to read the files you need a new version of paraview, not the standard ubuntu 22.04 version. After we figured this out (this took some time) we were able to run the example files.

[bstoeve]

The next thing was: What kind of mesh format does this code need? Well, we guessed nek5000. But of course we never generated a mesh for that before. But we managed.

[bstoeve]

Currently, we are on setting up a case. The general description of the case file contains everything in the case file. Now, if you want specific boundary and initial conditions, obviously you have to set the option "user". Then this is defined in a fortran code file of the case. This is all you get to know. There is no documentation on this file though.
Since the is no documantation about the meshing proceedure (you have to know this from nek5000 and its re2 files), this is all a bit confusing.

[timfelle]

The meshing is indeed generally poorly described in the documentations.
What i have done is design a mesh with some external meshing software (Coreform cubit in my case) and save the mesh as an Exodus-II file.
That file can be converted by a tool in the Nek5000 tool to the re2 file (https://github.com/Nek5000/Nek5000/tree/master/tools/exo2nek).
From the re2 file, Neko supplies a tool to run the conversion: https://github.com/ExtremeFLOW/neko/tree/develop/contrib/rea2nbin

I know it is a bit of a detour for creating meshes, but it have worked fine for me. I know that others use GMesh to create a mesh and then do a similar chain of conversions as i described above, which is automated by the authors of neko through their gmsh2nek tool: https://github.com/ExtremeFLOW/neko/tree/develop/contrib/gmsh2nek

[padakpra]

Hello

I am working with @bstoeve in trying to simulate flow around an airfoil using NEKO. I tried to mesh an airfoil in a structured mesh as shown in the files below using gmsh. I tried following the instructions as mentioned in the earlier post.

Initially tried converting a 2d mesh into gmsh2nek, it identified line elements and boundary groups. When tried executing the further conversion into .nmsh using rea2nbin, it popped up the following errors, but the file was generated in .nmsh format. When I tried to execute the file using a simple case file, it crashed.

Conversion of .msh to .re2

Conversion of .re2 into .nmsh

Please let me know where the error is - mesh file or the case file. I presume it is the mesh which is causing trouble and I am having trouble visualizing the converted binary format mesh .nmsh.

[timofeymukha]

Hi! I generally had more luck with using gmsh2nek directly on a 3d mesh. Neko is anyway a 3D solver, although you can try to mimic 2D by having 1 element in the 3rd direction and applying periodic conditions. I believe you hit a pain point here, these mesh conversions are not incredibly reliable in my experience.

Your input regarding the lack of documentation on I/O formats are on point, and this will be addressed. @vbaconnet is currently working on the documentation of the user file routines.

[padakpra]

Thank you for the reply, so I tried to execute the basic mesh conversion using a rectangular domain generated in gmsh. The conversion from .msh to .re2 was fine as shown below. When I tried doing the next part i.e. from .re2 to .nmsh I am getting the following error. Not sure where the error is. Is this the right forum to post these or do we have any other forum?

Conversion from .msh to .re2

Conversion from .re2 to .nmsh

[padakpra]

Ignore this, I have managed to convert the file.
Will try to execute the same with an airfoil inside. Thanks

[vbaconnet]

A quick fix for this could be to change the BoundaryIDs of your gmsh mesh to be strictly lower than 20. If you are able to, try to change the numbering to go from 1 to 6 instead of 43 to 48.

I think you are getting this error because the maximum number of labeled zones in Neko is fixed to 20. When the rea2nbin converter reads the boundary conditions it sees IDs 43,44,45,46,47,48 and throws the error here.

EDIT: didn't see your message above, but keeping this here in case someone runs into the same issue

[padakpra]

Thank you for the quick response on the queries. I tried to create few meshes with a rectangular domain with quad elements and then with an airfoil as shown below.

I discretized the domain in the following way to get quads and a structured mesh throughout the domain. When I run the case file the script stops at Fluid. Is it because of my discretization or what is the source of error?

The error:

Discretization of domain and volumes:

Thank you

[timofeymukha]

timofey "dot" mukha "at" kaust "dot" edu "dot" sa

[timfelle]

tife "at" dtu "dot" dk

[padakpra]

Thank you for the meet. @timfelle, you recommended running the case by lowering the polynomial order and having a large time step etc. can you mention it here and elaborate it

[padakpra]

Also, can you direct me to the link to the documentation under development?

[timfelle]

You are very welcome, I think it was @timofeymukha who mentioned the lower order polynomials.

The latest documentation can always be read here:
https://github.com/ExtremeFLOW/neko/tree/develop/doc/pages

The user file documentation we mentioned have been pulled in now, so it is in that link.

[padakpra]

Hello

Thank you for the support and insight, I was able to run the simulation as shown below. Re=10000; v = 1m/s; t = 10s with an ahw9300 airfoil. Though the mesh wasn't optimal the simulation ran pretty well.

Following are the points I have to mention:

1. The simulation ran when I specifically numbered the boundary ID's from 1 and not from any random number.
2. I am trying to use open-source tools to create high quality meshes around airfoil (construct 2D, gmsh, PyAero) but as the mesh generation follows a specific route manipulating meshes with constraints such as boundary ID's is a big issue.
3. Can I know if NEKO supports simualtion run continuation from the already run simulation? (let's say I want to extend the current simulation from t=10 to t =20)
4. OpenFOAM has lot of mesh manipulation tools inside and I have seen a repository foamtoNeko but it is mentioned as unusable. Can I know if any work is being done on this tool?
5. Also, will there be additions to the available default boundary conditions, i.e. slip wall?

Thank you
Pranay

[timofeymukha]

Hi, sorry I just saw this. Great to see that you got something running! Regarding your points.

3. Yes, you use checkpoints for that. To restart from a checkpoint, I think one uses the "restart_file" keyword. Take a look at the documentation.
4. Haven't put time in that one in the last year, so unfortunately a bit frozen. Not sure what the issue there was in the end, it almost worked...
5. Yes :-). Slip walls should hopefully be there in a matter of a month or two.

[padakpra]

Hello @timofeymukha and @timfelle

Good evening

Is there any way to visualise the re2 boundaries after we have converted using the gmsh2nek conversion. I want to see the assignment of boundaries to corresponding geometry faces. I tried to visualise the simulated results in paraview, but I cannot see the boundary details. I just see everything as a single mesh.

There is a hyperbolic mesh generator for airfoils - Construct2D and I am trying to use the mesh generated with this in NEKO.
https://sourceforge.net/projects/construct2d/

Thank you
Pranay

[padakpra]

Thank you @timfelle and @timofeymukha for your info. I guess I will set the output_boundary option and see if I can work with it

[padakpra]

@timofeymukha, I have defined the option in the casefile and got bdry0.f0000 and bdry0.nek5000 as output and I can visualise the nek5000 file, but I cannot see the boundaries. Also, I don't have an idea to visualize the f0000 file.

[timofeymukha]

When you open the .nek5000 e.g. in paraview it actually gets the data from the 0.f0... files. So the whole field has just one value or how does it look?

[padakpra]

it shows the entire mesh as a single block as the mesh was composed in that way but no info on the boundaries. Also, can you redirect if there are some options in paraview for boundary visualisation.

[timofeymukha]

Yeah the idea is that you look at the values of the "pressure" array. And then the boundary nodes will have a different value from the internal ones.

[padakpra]

Hello @timofeymukha and @timfelle

I tried simulating the airfoil using c-grid topology using couple of mesh development approaches below. In both the cases the solver diverged with iterations in all the directions reaching the maximum

1. A hyperbolic grid was generated using construct2D and then script was developed to split the farfield boundary for inlet and top and bottom walls. The boundary id's were redefined so as to begin from 1 and corresponding elements to be updated with the reassigned id's.
   below is the mesh and corresponding error file, case file.

Error:

casefile:

2. To see whether my conversion had any error, I simulated a c-grid using gmsh and it had the same error as above.

Error:

casefile:

Not sure what is causing the error, but another thing I have observed is that all the examples have rectangular domain and is it the case for the curved inlet? Please let me know what the source of error is.

PS: I have added snaps instead of files as it would be easy to follow.

Thank you
Pranay

[timofeymukha]

Some suggestions would be: set the time order to 2 for now, use periodic boundaries instead of symmetry on the sides. Set procjection_space_size to 0 for the pressure_sover.

There is a solution file dumped before the time loop. Did you look at that in paraview to see if things look reasonable, i.e. that your solution is 1 everywhere? You could also try to set the inflow value to something different to see if the value is correctly picked up at the inlet boundary -- basically jsut to visualize the inlet and confirm that it is working.

[timofeymukha]

There looks to be a strong jump in element size when crossing from the half-cicle part of the mesh to the rectangular one. It would be a good idea to never have jumps in size that are larger than 1.5x or so.

[padakpra]

Thank you for the suggestions @timofeymukha. I tried to rectify the mesh and run a case which initially ran for considerable timesteps and then the pressure diverged. I tried changing the projection space from 0 to 20 with intermediate values in between but it did not work. The mesh has periodic boundary conditions on the sidewall. I verified the boundary file which indicated different pressure values.

Can you also touch on the projection space (both for velocity and pressure), what it signifies and how should the value be selected. Also, I don't see the solution file before the time loop. If needed, I can attach the mesh file and the case file for your reference. Let me know if you see any problem with the above.

Regards
Pranay

[timofeymukha]

Hi! Well, the fact that you are getting the solver to strart crunching is very good. Since you know about when the simulation will crash, you can tweak the output frequency such that it will give you a file right when things are about to go wrong. Then you can look at that in Paraview and hopefully the problematic region will be exposed.

Regarding projections, it is an acceleration method for find the solution: you search for it in the space of previous solutions. That is at least my understanding. However, this process seems to not be very stable. So I would say for problematic cases it is better to turn it off.

Many of us will be at a conference starting tomorrow. But I think next week if you are still having troubles it would be good if you can share your case so we can have a quick look.

[padakpra]

Can you guide me in how to visualise the f0000 files?

[padakpra]

Hello @timofeymukha and @timfelle

Good afternoon.
Finally we have managed to run a simulation with coarse grid as shown below, where the simulation runs for 1.2 sec and the solver diverges as we see a problem near the inlet and outlet interface (Top and bottom) also there is a wave like behaviour in the inlet sector. Perhaps I will try to create another sector with inlet and define the top and bottom sectors as symmetry. Our understanding is the code couldnot comprehend the velocity vector which is parallel to the cells at the top and might be causing this. Can you share your view on why this might be happening.

Also I have used variable time step in the script, which should be implemented but I donot see the timestep being increased when CFL number drops. Is it beacuse I have an older version 0.7.99?

Thank you

[padakpra]

Hello

Can i know whether the bug is rectified and the possible release date of version_0.8.0?

Another update : I managed to run a rectangular domain with periodic BC on the z walls and there seems to be an instability rising as shown below around 4.5 sec in a 20sec simulation. One element in the z direction.

I managed to get another run going with z direction wall BC to sym and currently the simualtion is running around 11 sec 🤞

[timofeymukha]

Could you share the case and mesh for the periodic bc one? I would like to take a look at what is going on.

[padakpra]

Here you go. Let me know if you need anything else.

Thanks

2d_cylinder.zip

[padakpra]

Hello @timofeymukha, the 20sec run with symmetric BC completed successfully. I am trying to analyse the following points

1. Is it possible to extract normalised wall distance from subroutines near the walls for analysis of y+. I beleieve NEK5000 has subroutines inside which extracts this parameter. I tried to see if paraview has a filter inside to analyse simulated fields but couldn't figure it out.

2. I am trying to plot Cp, Cl, Cd coefficient distributions and I have seen that there is a subroutine calc_lift_from_field inside contrib. Can I use this to figure it out or is there any other way?

Thank you

[timofeymukha]

Hi!

1. Nothing built into Neko, you would have to do it yourself at a post-processing stage e.g. with Python or so by looking at the velocity field or directly at its gradient.
2. Yes, i think you should be able to compute some of this stuff with the utility. @MartinKarp would be the question to ask details, but he is a busy now since he is defending is PhD in like 10 days :-)

[padakpra]

Hello @timofeymukha,

This is regarding the nekton map file utilised for .rea2nbin file where I managed to generate .ma2 file (using genmap from NEK5000 for elements < 150000) and .co2 file (using gencon from NEK5000 for elements > 150000) but when I run the rea2nbin it looks for a nekton map file but it cannot identify these files. Can you suggest if this map file is necessary (which I believe is as it has the connectivity data inside it) and how to make sure the rea2nbin reads it.

tried renaming the file to .map but there seems to be an error

Thank you

[timofeymukha]

I am pretty sure the map file warning can be ignored. @MartinKarp do you know?

[MartinKarp]

Yeah, the map file can be ignored. I would recommend not using it actually. The global numbering is generated in neko during the conversion for re2 to nmsh. Then this is used to generate the connectivity when the simulation is run.

[padakpra]

Understood.. thanks a lot :)
As the convertor was showing the message, thought that this was needed while execution.

[padakpra]

Hello @timofeymukha

Good morning.
We are planning to run a DNS on an airfoil utilising 4xNvidia GPUs and would like to clarify these points.

1. We managed to run a test case to utilise all the 4 GPU's but I am not sure it has utilsied all the GPU's. The simulation consists of ~14000 elements. I have attached the scripts for the execution in the cluster below.
   

2. Is it possible to use RANS simualtion carried in NEK5000/ OpenFOAM as initial conditions in NEKO? If yes, can you direct me to the page/ tell me the procedure?

3. There is no definitive procedure of post processing possibilities in NEKO. Can you touch on how to utilise Python for reading the field files? I would like to know if KTH toolbox can be utilised for carrying out the postprocessing on the simualtion? Can you suggest resources or the path to see the use of this toolbox? This is needed as I need to extract Cp, Cl, Cd, other parameters from the simulation.

4. Do you have any more info on the use of statistics, simulation components, probes etc. It is not clear in the documentation on how to identify the GLL points using .csv file etc. Also when I tried to use it in the script as shown below the simualtion stalled after the gather scatter o/p.
   

5. Is the targetcfl update feature not implemented in the 0.7.2 version? I donot see the time step update accordingly in the output log.

6. Does the code doesnot consider the free stream velocity in mu calculation? i.e I have given the initial velocity as 2m/s with Re: 10^4, l=1, it shows mu as 10^-4. So irrespective of velocity input, it considers v=1m/s in calculating mu?

Thanks again for your patience in responding to all the queries :)

[timofeymukha]

This is indeed quite a bit so let me give some quick answers to some of the stuff.

1. There should be 1 MPI rank per GPU and you need to map the MPI ranks to the GPUs properly. For 14000 elements 1-2 GPU is probably enough.
2. No
3. You cannot use the Toolbox. But there is https://github.com/ExtremeFLOW/pyNekTools , under development.
4. You should ask more concretly what you are missing in the documentation and then we can try to answer and add it. "It is not clear in the documentation on how to identify the GLL points using .csv file etc." I don't quite understand. The probe points don't have to coincide with GLL nodes.
5. No.
6. Correct. If you use non-dimensional quantities to setup the material properties, the reference scales are considered to be unity.

[padakpra]

Thank you for the replies. I was not expecting an immediate answer :)

3. Regarding the postprocess - as mentioned I need to plot Cp, Cl, Cd plots and need some info about how to extract them from the field data. I tried using paraview but was unsuccessful. I you can provide some info on priciples of the PyNekTools, I can try to code and work with it. I have a deadline for my thesis in July and need these parameters.

4. Lets say for the probes, it is described as follows with points.csv. I am not sure how the coordinates of these points work.
   

[timfelle]

In regards to 4.

The file probes.csv should contain a list of locations in the domain for which you wish to sample the selected field. Each line in probes.csv should contain x, y and z coordinates in the global domain.

The fields probed in the example you linked are the z-component of the velocity (w) and the scalar field (s) representing temperature quite often.

The output probe file is structured as follows:

• First a list of all points probed (so you don't have to read the input file to get the output).
• Then all the probed fields for each point for each timestep.

p0_x, p0_y, p0_z
p1_x, p1_y, p1_z
p2_x, p2_y, p2_z
t0, w(p0), s(p0)
t0, w(p1), s(p1)
t0, w(p2), s(p2)
t1, w(p0), s(p0)
t1, w(p1), s(p1)
t1, w(p2), s(p2)

[padakpra]

Thank you for the info @timfelle, I will try to implement it during my run.

@timofeymukha, Can I know if @MartinKarp is available and would like to discuss the possibility of evaluating the Cp, Cl, Cd from the field files.

Also, I am trying to simulate the inlet velocity at an angle as the airfoil is an angle but I see my simualtion diverging. I assigned the top and bottom walls as periodic with inlet vector in x and y direction. Another case it was as shown in picture. Do let me know if there is any mistake or any way around it.

exaggerated

[padakpra]

Hello @timofeymukha and @MartinKarp

I am trying to execute the following code and it seems to work when I remove the statistics and simulation components block. I need to write the output from 10sec sampling at every 0.5sec and writing for every 1sec.

{
"version": 1.0,
"case":
{
"mesh_file": "75000_v3_z_5_sim.nmsh",
"output_at_end": true,
"output_boundary": true,
"end_time": 20.0,
"timestep": 1e-6,

"numerics": 
  {
    "time_order": 3,
    "polynomial_order": 5,
    "dealias": true
  },

"fluid": 
  {
  "scheme": "pnpn",
  "Re": 75000,

  "initial_condition": 
  {
    "type": "uniform",
    "value": [1.0, 0.0, 0.0]
  },
  
  "inflow_condition": 
  {
      "type": "uniform",
      "value": [1.0, 0.0, 0.0]
  },


  "velocity_solver": 
  {
    "type": "cg",
    "preconditioner": "jacobi",
    "projection_space_size": 0,
    "absolute_tolerance": 1e-9,
    "max_iterations": 800
  },

  "pressure_solver": 
  {
    "type": "gmres",
    "preconditioner": "hsmg",
    "projection_space_size": 0,
    "absolute_tolerance": 1e-4,
    "max_iterations": 800
  },
  
  "boundary_types": ["v", "o+dong", "sym", "sym", "w", "sym", "sym"],
  "output_control": "nsamples",
  "output_value": 20,
  
  "output_checkpoints": true,
  "checkpoint_control": "nsamples",
  "checkpoint_value": 4,
  
  "variable_timestep": true,
  "max_timestep": 1e-4,
  "target_cfl": 0.4
  
 },
 
 "statistics":
 {
      "enabled": true,
      "output_control": "simulationtime",
      "output_value": 1.0,
      "start_time": 10,
      "sampling_interval": 500000
 },  
 
 "simulation_components":
 [
    {
      "type": "vorticity",
      "compute_control": "tsteps",
      "compute_value": 1000000,
      "start_time": 10     
    },

    {
      "type": "lambda2",
      "compute_control": "tsteps",
      "compute_value": 1000000,
      "start_time": 10
    }
 ]

}
}

Am i declaring some of the parameters wrongly?

Thank you
Pranay

[timofeymukha]

  "output_checkpoints": true,
  "checkpoint_control": "nsamples",
  "checkpoint_value": 4,

  "variable_timestep": true,
  "max_timestep": 1e-4,
  "target_cfl": 0.4

This stuff should not be inside fluid, but rather directly under case.

Sampling statistics so seldom makes pretty much no sense. Sample every 10 timesteps or so. Moreover for dt=1e-4 500000 iteroation is 50 s, not 5. So you sample more seldom than you write.

I don't think simulation components have the concept of start_time, so that is probably just ignored. With your compute_value your interval is 100 s, I guess you are off by a 0 here as well. That may be causing the crash. But I am not sure.

[padakpra]

dt here should be 1e-6 not the max timestep 1e-4 @timofeymukha ?

[timofeymukha]

If you use the adaptivity, it will be <= 1e-4, right?

[padakpra]

True, so my calculations need to be around this value but I haven't seen this being executed till now. Maybe my declaration is wrong and I will push the block to the case level.

[padakpra]

Hello @timofeymukha , I found that the script works well if I remove the statistics block. Is the sampling frequency too low or the options being declared twice, could you clarify? The same block worked well when I implemented in 2d_cylinder example.

{
"version": 1.0,
"case":
{
"mesh_file": "75000_v3_z_5_sim.nmsh",
"output_at_end": true,
"output_boundary": true,
"end_time": 20.0,
"timestep": 1e-6,
"output_checkpoints": true,
"checkpoint_control": "nsamples",
"checkpoint_value": 4,
"variable_timestep": true,
"max_timestep": 1e-4,
"target_cfl": 0.4,

"numerics": 
{
    "time_order": 3,
    "polynomial_order": 5,
    "dealias": true
},

"fluid": 
{
  "scheme": "pnpn",
  "Re": 75000,

  "initial_condition": 
  {
    "type": "uniform",
    "value": [1.0, 0.0, 0.0]
  },
  
  "inflow_condition": 
  {
    "type": "uniform",
    "value": [1.0, 0.0, 0.0]
  },

  "velocity_solver": 
  {
    "type": "cg",
    "preconditioner": "jacobi",
    "projection_space_size": 0,
    "absolute_tolerance": 1e-9,
    "max_iterations": 800
  },

  "pressure_solver": 
  {
    "type": "gmres",
    "preconditioner": "hsmg",
    "projection_space_size": 0,
    "absolute_tolerance": 1e-4,
    "max_iterations": 800
  }, 

  "boundary_types": ["v", "o+dong", "sym", "sym", "w", "sym", "sym"],
  "output_control": "nsamples",
  "output_value": 20     
 
},
 
"simulation_components":
[
    {
      "type": "vorticity",
      "compute_control": "tsteps",
      "compute_value": 1  
    },

    {
      "type": "lambda2",
      "compute_control": "tsteps",
      "compute_value": 1
    }
],
 
"statistics":
{
    "enabled": true,
    "output_control": "simulationtime",
    "start_time": 10,
    "output_value": 1.0,
    "sampling_interval": 50
}

}
}

Thank you
Pranay

[MartinKarp]

and do you restart the simulation or run from 0. Statistics also adds quite some memory requirements, which might also be something to keep in mind.

[padakpra]

The simulation is from the start and my code is same as turbchannel just that the value of start time and output value is adjusted accordingly. I will increasing sampling interval.

[padakpra]

can we have a quick meeting tomorrow to debug what is the issue? I guess its the RAM but I have increased the sampling as well as the writing to 1sec but its still the same

[adperezm]

Perhaps not your problem, but changing the sampling time will not change the memory requirement of the statistics. A bunch of extra fields are allocated for the statistics (around 50!) but they are updated in place. So the extra requirement will remain 50ish fields whether you sample every timestep or every 100th.

Perhaps this is a good point on having some "light" statistics, where the overhead is smaller. This has been discussed but never executed.

[padakpra]

Yes, that's true I wanted to see the memory requirements and tested out with a smaller case and it ran as the memory requirements were met.

Thanks again

[padakpra]

@timofeymukha and @MartinKarp, can I know if there is a way to simulate turbulence intensity for the inflow wind in NEKO?

[timofeymukha]

You will need to prescribe some resolved turbulent structures at the inlet. Currently, we don't have any synthetic turbulence generators in the code, so the only option would be to implement it yourself, but is not so trivial, of course.

Perhaps an alternative would be to prescribe some volumetric forcing upstream of the foil to trigger some perturbations. But again, this is not super-trivial.

[padakpra]

@MartinKarp and @timofeymukha

I am trying to execute the same case on the CPU nodes but I have this error when i have the case running. Can I know what is the source of the problem?

Is it missing some requirements or my execution is wrong.

{
"version": 1.0,
"case":
{
"mesh_file": "75000_v3_z_5_sim.nmsh",
"output_at_end": true,
"output_boundary": true,
"end_time": 20.0,
"timestep": 1e-6,

"numerics":
  {
    "time_order": 3,
    "polynomial_order": 5,
    "dealias": true
  },

"fluid":
  {
  "scheme": "pnpn",
  "Re": 75000,

  "initial_condition":
  {
    "type": "uniform",
    "value": [1.0, 0.0, 0.0]
  },

  "inflow_condition":
  {
      "type": "uniform",
      "value": [1.0, 0.0, 0.0]
  },


  "velocity_solver":
  {
    "type": "cg",
    "preconditioner": "jacobi",
    "projection_space_size": 0,
    "absolute_tolerance": 1e-7,
    "max_iterations": 800
  },

  "pressure_solver":
  {
    "type": "gmres",
    "preconditioner": "hsmg",
    "projection_space_size": 0,
    "absolute_tolerance": 1e-4,
    "max_iterations": 800
  },


  "boundary_types": ["v", "o+dong", "sym", "sym", "w", "sym", "sym"],
  "output_control": "nsamples",
  "output_value": 400,

  "output_checkpoints": true,
  "checkpoint_control": "nsamples",
  "checkpoint_value": 4,

  "variable_timestep": true,
  "max_timestep": 1e-2,
  "target_cfl": 0.4

 },

 "simulation_components":
 [
    {
      "type": "vorticity",
      "compute_control": "tsteps",
      "compute_value": 1
    },

    {
      "type": "lambda2",
      "compute_control": "tsteps",
      "compute_value": 1
    }
 ],


"statistics": {
    "enabled": true,
    "output_control": "simulationtime",
    "output_value": 0.05,
    "start_time": 10.0,
    "sampling_interval": 1000

}

}
}

output:

[njansson]

@timfelle and @njansson , the code has been recompiled without CUDA and case initiated but stuck at gather scatter step. Is it because of too much of memory allocation? Allocated cores: 256, 2 cores

Thank you

Note, if the log has printed the gather-scatter statistics, it has passed the initialisations of gs. Thus it is not stuck at gs, rather it is busy going through the large amount of dofs elsewhere in the code

[padakpra]

How much memory do you have available per core? And can you try to build a toy problem with the same setup, if you can run an example with say 10 elements per core. Even if it diverges on the first timestep, we can narrow the problem down to be related to mesh size.

Btw, @timofeymukha I remember you mentioning that there ius a lower bound for how many elements must be available per core, but 10 should be fine right?

The ram per core is around 6 GB, so I guess I allocated few cores for the problem.

[padakpra]

@njansson , @timfelle and @timofeymukha - I ran the cyl_bl case on the cluster using different cores of 1,2,3 and 8. It started time stepping for cores 1 and 2 but halted at gather scatter from 3 cores. The case is of 7000 elements.

how long will I have to wait for iterations to start?

[padakpra]

Apologies for the previous comments but now I understood that the halt at gather scatter is not an error and it starts iterating after a while and I got my run executed on 6 nodes (6x128 cores). Generally, what is the min. RAM requirement per one element? I understand it depends on the BC and the statistics but just to get a rough idea of estimating no of cores needed for no of elements

[timofeymukha]

I mean, you have (p+1)^3 reals per element per field. In double precision, each real is 8 bytes, right? From this you can try to get an estimate, but of course there is more stuff like the mesh, temporary fields appearing here and there, etc.

[padakpra]

@njansson, are you free for a meet sometime this week? It is going to be a short discussion with my colleague, and we can probably look at the execution of NEKO using multi-GPU's and on CPU

Thank you
Pranay

[padakpra]

Hello @timofeymukha / @njansson / @MartinKarp / @timfelle ,

A general question on how much aspect ratio NEKO can handle, in my case as I am resolving the boundary layer for y+<1 , I am getting high aspect ratios near the wall of the airfoil.

What is the status of the o+dong bc, can it be implemented as I ran my case with o+dong and it diverged after 30% completion

Thank you
Pranay

[timofeymukha]

There is no exact answer regarding the ration. You can certainly handle say 50.

The issue is probably not with o+dong per se but with our boundary conditions now working with certain mesh topologies. This has not been fixed yet.

[padakpra]

I mean getting 50 is difficult especially near the walls :)

[timofeymukha]

Keep in mind that you are doing a DNS, so even x+ 50 is coarse. Even for LES that would be pushing it.

[padakpra]

Hello all.

I have few queries below:

1)Can I know if there is a way to capture the output exactly where the solver diverges as its difficult in predicting the frequency and re-running the simulation or starting from a check point

2. is there a way to run 2D simulations in NEKO?

Thanks

[timofeymukha]

Hi!

1. No automatic way at the moment.
2. Not true 2D, but you can have a single element in the 3D direction and set periodic boundaries on the corresponding sides. That mimics 2D pretty well.

[padakpra]

Hello all,

I have been trying to simulate the Re 75000 around an airfoil and following ppt describes the mesh configurations and Reynolds number achieved. I could manage to run Re 10000 and one key factor observed is whenever there are iterations in the z direction the pressure iterations ramp up and there is divergence.

Assessing the mesh quality parameters is a big gap in spectral element mesh and how do I ensure that the z direction has enough length and elements for characteristics to be captured? Also, is there a way for better debugging i.e. knowing where the error is?

Attached a presentation, excel and a simulation test file for verification. Please let me know if you need any more information. I know its too many questions, maybe a short meeting would be helpful.

Thank you
Pranay

NEKO Queries.pptx
Run_analysis.xlsx
New folder.zip

Results: https://fraunhofer-my.sharepoint.com/:f:/g/personal/sai_pranay_padakanti_iwes_fraunhofer_de/ElHWz3y_7VtKqgzgB8e31hEBqfFi-XZUtCKBT_unrcQnig?e=5egBJt

[timofeymukha]

Hi! Generally, the meshes for which one sees the details in your ppt suffer from very abrupt jump in element sizes. This is generally not suitable for high-fidelity CFD, and I would say bad practice for CFD in general. A jump i size larger than a factor of 1.5 is to be avoided, unless you are really somewhere in the free stream where the solution is constant. For SEM, one should be extra careful with this stuff.

The fact that checkMesh (right?) fails the well posedness check for some cells also looks very scary, did you check that?

Keep in mind that what you are trying to run is a DNS. For this type of simulation the mesh should generally be perfect, and well refined near the foil in terms of inner scales of the boundary layer.

Simulating in 2D could yield something at low Re number when the flow is still laminar, but for turbulent flow they are not valid. Turbulence is always 3D.

[padakpra]

I believe the abrupt changes in mesh is valid only for the Executed Run I mesh which was one of the first meshes I ran for testing though this ran for Re 75000 till 23%. Rest of the meshes specially configuration I is fine I guess in terms of jump in sizes

checkMesh is fine for configuration II except for determinant which is in the far field region.

Can you comment on the domain sizes and their influence especially in the z direction?

Thank you

[padakpra]

Just to add to that, I tried to reduce the z dimension to 0.1 and get the z scales to get kolmogorov lengths z+ ~20 but still the solution diverged

[padakpra]

Hello all,

Can I know if the variable timestep feature is being implemented as I have not seen implemented though I have called it during the run. Also, how to decide on the dt for the simulation, my simulation crashed around 0.13% when I went with order 6 and dt 1E-5 whereas it ran for around 10% when I changed the polynomial order to 3 and dt to 1E-4.

I know that it depends on the cfl and the velocity of the flow, but how to really determine the time step depending on the polynomial order?

Thanks

[timofeymukha]

Hi! Sorry, I simply have no time at the moment, and next week I will be travelling. We are all participating in the CTR Summer program, which is very intense.

[timofeymukha]

For the variable timestep, if it does start varying, you should be in business :-)

[pschlatt1]

We will schedule a meeting with Siavash the coming days, addressing some of the issues.

[Shiyu-Sandy-Du]

Can you analyze what is the error in this case which ran for 10%.

75000_per.zip

In this zip file, you don't activate the variable time stepping, do you? Could you share your log file in which you have variable time stepping?

[padakpra]

I haven't pulled in the develop version of NEKO which supports variable time stepping, currently using 0.7.2 version. I will try to execute the runs on the develop version now using variable time stepping.

Regarding the o/p log here you go:
New folder.zip

This was executed with p =3 and Δt = 1e-4 sec and with an increase in p the solution diverged quickly

Thank you

[padakpra]

Hello can I know when will @timofeymukha and @MartinKarp will return from vacation?
We have been trying to execute a run for Re 75000 case on guidance from Mr. Siavash Toosi but there has been no positives in getting the simulation running. We managed to run the 10000 case but not sure where the problem is - mesh (mesh requirements, gmsh, bc), case file, NEKO setup in the cluster etc.

Thanks
Pranay

[padakpra]

Hello @MartinKarp

Regarding the post processing, I wanted to know what the fields in stats meant?

⟨pp⟩ is the temporal average of pressure^2 or is it the fluctuations? same question for u, v and w indicated in the picture below

also is it the same here, the temporal average of pressure*du/dy and so on?

Just an update on the run, it is running around 60% and will keep you updated :)

Thanks
Pranay

[MartinKarp]

It should be at the top of the documentation
https://neko.cfd/docs/release/df/d8f/statistics-guide.html

u, the instantaneous field.
⟨u⟩t, the temporal average of u.
u=⟨u⟩+u′, the Reynolds decomposition of u, where u′ is the fluctuation of u around the mean field.

So not fluctuations, but the actual field.

[pschlatt1]

This follows standard notation within turbulence, i.e. the brackets correspond to a time (and potentially spatial) average, and primes indicate fluctuating quantities. You can then construct the fluctuations (variance) afterwards as <u'u'> = - .

There is a recent journal paper on the KTH framework for Nek5000, which describes this quite well: https://www.sciencedirect.com/science/article/pii/S0010465524001723

[padakpra]

I understand the notation, but I was wondering about the application of fields as we already have

. I will look into the paper for my understanding.

Thank you

[pschlatt1]

prms^2 = <p'^2> = <(p-

)^2> = <p^2> +

^2 -<2p

> = <p^2>-

^2. So pp is <p^2>.

[padakpra]

Hello @MartinKarp

Can I know what is meant by zone number in the script of contrib/calc_lift_from_field/calc_lift_from_field.f90 that we need to pass as arguments?
is it the boundary ID?

Thank you

Thank you

[padakpra]

Will it work if I do it with 0.7.2 or 0.8.0? I did not run with the latest version as the code had commit that it was updated last year and thought it would be the same

[padakpra]

@MartinKarp is the problem in the rea2nbin version or calc_lift_from_field version?

[MartinKarp]

rea2nbin

[padakpra]

can I work with the develop version as I see changes being made in the gmsh2nek or should I proceed with the 0.8 version?

[padakpra]

@MartinKarp , I got the ID's in the new version, dont we get the periodic BC as an output

[padakpra]

Hello all,

In my cases, I tried implementing the statistics from t=10 sec till t = 20sec and I get mean fields every 0.05 sec. I implemented the same in case file below

I used the average_fields_in_time routine to average out the mean_fields by implementing the command

~/Desktop/NEKO/source/neko_0.7.2/contrib/average_fields_in_time/average_fields_in_time mean_field0.fld 10 mean_field0_avg.fld

When I observe the mean_field_avg the domain velocity is around 0.28 but the inflow velocity is around 1.0 m/s. When a temporal average is done, it should be around 1.0 m/s right?
Instantaneous field at t=10sec:

mean_field_avg

I am not sure if my understanding is incorrect/ mean fields average consider initial fields as well i.e. from t=0/ T0 is incorrect. Please comment.

Thank you
Pranay

[MartinKarp]

Check for example the turb_channel example for some idea on parameters. Then we would also often run for quite some time (I guess around 100 time units based on chord length and free stream velocity might be reasonable, but I haven't really run many wings).

[padakpra]

Thanks will look into it..

[timofeymukha]

No-no, @MartinKarp , you don't need that much averaging time for an airfoil. Even 5 c / U_inf will give you a very good mean if the transient has been removed. Channels are kind of horrible from this perspective because of periodicity.

[timofeymukha]

(but as usual, it is better to look at local delta/utau, in attached flow regions)

[timofeymukha]

So, if @padakpra if you have c=1 and Uinf=1, then starting at 10 and ending at 20 is quite reasonable.