CFD: First Steps for Newbies - Pooven

[Pooven]

Hi there, I'm fairly new to Computation Fluid Dynamics and the engineering concepts as a whole. I believe CFD is generally a Chemical and Mechanical Engineering interest but I thought I'd post in this sub forum. I'd like to share what little I know and learn from others.

I understand that the basic idea is to model the continuous domain into a discrete domain usually by using a grid or mesh. I'm particularly interesting in the 3D discretization where the mesh is made up of http://en.wikipedia.org/wiki/Hexahedral" (and perhaps other) geometries. So the first step is to define the geometries to be used.

I think, that the next step is to simplify the https://www.physicsforums.com/showthread.php?t=6148" on Wikipedia for brief overview. Wikipedia refers to these as discretization methods.
I also read this tutorial which suggested the conversion from a differential equation to an algebraic equation (but it was within a 1D domain): http://courses.cit.cornell.edu/fluent/cfd/intro.pdf"

There is also a http://www.cfd-online.com/Wiki/Turbulence_modeling" that helps to describe the seemingly 'random' behaviour of particles.

What I'm a bit confused about is how the equations are solved over the geometries of the mesh. Suppose I have a hexahedral celled mesh... would the equations be different from a mesh with tetrahedral cells? Is it the discretization method that adapts the equations to the chosen geometries? The equations I mention here refer to the governing mass and momentum equations for the flow.

Any suggestions or insight would be deeply appreciated.
Kind regards,
Pooven

 

[minger]

Firstly, my adviser would frown for you calling turbulence random. Many of us prefer chaotic, but that's another story. You have the steps fairly understood. My area of background is in fully structured, finite difference codes. Structured means that essentially any grid point in the domain can be referenced by a single index for each dimension. Think of an evenly spaced cartesian grid. I can reference any point by a simple (i,j,k) pointer.

With grids like this, the discretization can be fairly easy. Whereas finite volume is an integral approach (summing fluxes in and out of an element), finite difference takes derivatives at grid points to determine the time derivatives.

In finite differencing there is an additional step between meshing and actually solving. Analysis of the Taylor series derivative approximations show that they are very unstable in two or more dimensions. To avoid the unnessesary terms, we "convert" the grid to a "perfect" evenly spaced cartesian grid. This is done by taking derivatives of the grid itself to form the Jacobian grid metrics. Once the grid is converted into what we refer to as computational space, we can start taking derivatives.

You did not mention one, and possible the most important part of the process: boundary conditions. Boundary conditions specify what happens at walls and the ends of the domain. For years, these have been a pain the ***, and just as much work goes into avoiding the use of them as does the refinment of them.

As far as solving between different mesh shapes, the answer is yes and no. In the most general case, we are solving the same equations
$$\frac{\partial \vec{Q}}{\partial t} + \frac{\partial \vec{E}}{\partial x}+\frac{\partial \vec{F}}{\partial y}+\frac{\partial \vec{G}}{\partial z} = S$$
However, what changes is, as you mentioned, the discretization methods. The methods will be difference for the different shapes. Further more, there are different stencils for each shape. In our world of finite differences, we have one-sided schemes such as:
First-Order Backwards
$$\frac{du}{dx} = \frac{u_{i+1}-u_{i}}{\Delta x}$$
Two-sided schemes such as the second-order central
$$\frac{du}{dx} = \frac{u_{i+1}-u_{i-1}}{2\Delta x}$$
Additionally there are high order DRP schemes which play with the coefficients to get high accuracy at lower grid spacing, there are compact schemes which use derivatives themselves to obtain the derivatives. It can get confusing.

 

[Pooven]

Hi minger , thank you for reply! I actually agree with your advisor about turbulence not being random. I think there are cumulative forces in nature that act together and can perfectly explain phenomena such as turbulence - I was making reference to the explanation given on the http://www.cfd-online.com/Wiki/Intr...e/Nature_of_turbulence#What_is_turbulence.3F" website.

I do remember reading about structured grids and how they are stored in memory (using a 3D array if I recall properly - but it corresponds to what you just said anyway). On that note, under the FAQ section of http://www.cfd-online.com/Wiki/General_CFD_FAQ" that unstructured meshes over a terrain (for use with air movement) has proven successful - though I'm not sure how much more successful over a structured mesh.

In the tutorial I mentioned, they did explain the use of the Tylor series in the infinite difference method - and about the truncation error. Though the extra step seems to need a good understanding of the physics. And then I'm confronted again with the question of the more appropriate method to use. Albeit I haven't read up on this enough, all I know is that Wikipedia states that the finite volume method is the standard or classical approach.

I did read about the boundary conditions - and failed to mention them because I was still unable to understand how the equations governing motion was applied to the primitive elements of the different meshes. I'm actually hoping that I can ignore the boundary conditions but it's quite possible that I don't fully understand them; in my case, if I had a 90mx90m area of land, would the boundary conditions specify what happens at the cells adjacent to the perimeter of my land? If so, can I not solve for the other cells ignoring those at the boundary (since I'm not that interested about the air movement in those cells)? In the literature I read, a boundary condition similar to the work by Richards and Hoxey (2003) was used - though I'm not sure what that means! You did mention that effort was placed in trying to avoid using boundary conditions... so can it be avoided?

So based on the geometry used, I will have fewer discretization methods that are applicable? I am at a loss with the mention of different stencils for each shapes... I've started reading ' by John David Anderson' (recommended via the book guide on cfd-online.com), so I'll hopefully understand what you spoke about soon Is there a cheat-sheet of sorts that provides a guideline that shows which methods are applicable to which mesh geometries?

 

[minger]

OK, let's try to answer some of these for you:

'What are the advantages/disadvantages of structural and unstructured Mesh?

The first main advantage of unstructured meshes are the meshing process themself. Unstructured meshes are extremely easy to make, and computers are very efficient at creating them. Furthermore, they are easy to refine as well. If I have a group of triangular elements, I can refine that area by simply splitting a group of triangles into more, and adjusted the connectivity.

With structured meshes, because the gridlines extend all the way to the boundary, refining meshes can not only be a tedious process, but it also typically forces the grid refinement all the way to the boundary, which is computationally expensive. Not good.

The advantage of structured meshes, aside from easy storage is the derivative stencils themselves, with particular importance to high-order schemes. Because structured grids have grid lines which extend all the way to the boundaries, we can develop derivative schemes that employ lots of data (6th-order DRP scheme is an 11-point stencil!). This gives the derivative schemes high order of accuracy. DO NOT confuse order of accuracy with accuracy. The order of a scheme dictates how the numerical solution approaches the actual solution as the grid becomes more fine.

In my area of aero-acoustics, we are trying to resolve waves that are orders of magnitude smaller than the mean flow. That means we need great accuracy. Let's examine two schemes to look at number of grid points vs accuracy.

A first-order backwards scheme obviously has an order of 1. This means that if I double the grid, the error decreases by a factor of:
$$\Delta\epsilon = 2^n = 2^1 = 2$$
However, a sixth order scheme is:
$$\Delta\epsilon = 2^n = 2^6 = 32$$
We see that the the high-order scheme approaches the actual solution 16 times faster than the first-order scheme. At this point it's important to say that the low order scheme might perform better than the high order scheme with a coarse grid. However, if one wants a certain small accuracy, then using a high-order scheme will produce better results with lower number of grid points. Using low-order schemes is possible to get CAA results, but typically the amount of grid needed to resolve high gradient areas is staggering.

boundary conditions

You are correct about the boundary conditions. However, no matter how far away you place them, any numerical run is finite and needs them. The boundary conditions essentially try to replicate the real world. You want the code to think that they're not there. However, they can cause all sorts of problems.

Some of the first types of boundary conditions are the Thompson-style BCs. They are what's called characteristic style. They specifiy the characteristic waves as they leave the domain. A common implementation is to specifiy a constant stagnation pressure. Think of a wave leaving the domain. In order to maintain a constant stagnation pressure at the outlet, the BC will send an incoming wave back into the domain in order to maintain this pressure.

The problem is that you now have reflecting waves bouncing all over the place and its hard to determine what is real, and what's not. Some of the first attempts to get around this was grid stretching. Essentially if you dig a little deeper into derivative stencils, you'll see that the dissipation of schemes is related to the grid spacing. If I stretch the grid, the grid is no longer able to resolve the propagating waves, and the wave will "die". If these waves are "dead" before they get to the boundary, then the BC has no work to do; the BC can't screw anything up.

You have a good start with Anderson; anything that he writes is worth reading. There is no "cheat-sheet" of sorts, and many people in the craft will say that CFD even now is quite the "gray" art. Two different CFD'ers can get quite different results. If you're using commercial software, then my strategy is just to try and get brick elements, otherwise just let it do what it wants, you typically can't tell it what to do anyways, nor do you know how its solving.

So hopefully some of this helps.

 

[Pooven]

Hi again minger. Thank you for explaining some of the advantages and disadvantages of the two methods... I'm not sure how to phrase this, but I'd assume that the mesh representing the object should as close as possible of a match to the geometry of the object; that is, for my situation (where a digital elevation model is used), the larger the gap between the generated mesh and the actual topography of the land, the less accurate the results might be? Though I suppose this too could be negligible considering the many other factors influence the flow of air over land (such as terrain cover - grass, forest, etc) and such properties will not be accurately represented.
I suppose what I really asking is, how important is grid refinement? I hope I'm not misunderstand grid refinement > I think it's the process of making the grid match the geometry with greater precision? In industrial processes, I suppose grid refinement is very important especially since the flow, say through a tube, is controlled (since it's free from the effects of say ambient wind).
I suppose we should always aim for a complete as possible system and so I think in my situation where there are complex geometries, an unstructured mesh would be better suited not only because it's easier to implement but because it's easier to refine?

I was reading up on stencils - http://en.wikipedia.org/wiki/Five-point_stencil" )? And as you said, different shapes have different stencils, which now makes sense A quick question, do all discretization methods use stencils?

Was was going to ask if higher order is in reference to the truncation of the Taylor series? I suppose in a way; I read, 'High-Order Compact (HOC) Schemes, in the context of my research, are numerical approximations which achieve high-order accuracy on a http://www.cisl.ucar.edu/css/staff/spotz/research/hoc.html#compact"
This concept is still a bit abstract to me since I haven't done any applications... is the order of accuracy directly related to which term the http://en.wikipedia.org/wiki/Taylor_series" truncation occurred? And this is different from a higher-order scheme? I think I should read more before I ask these question so let's put that on hold.

As for the boundary conditions, in the Dr. Anderson's book it reads, 'Finally, the physical boundary conditions and their appropriate mathematical statements will be developed. The governing equations must be solved subject to these boundary conditions.' - Ch2. Not that I was doubting you, just reaffirming what you said. Thank you for the example, that added much clarity... the events that occur at the boundary could impact the overall fluid flow and this can't be ignored anyway.

I find that this discussion has been quite useful; especially with my overall understand of the subject matter - as I am just a beginner. So I'm very grateful to you minger! Thank you ever so much :D

 

[Pooven]

Salome

I'm trying to generate a mesh over my terrain but I'm not sure if my approach is correct. I'm using Salome and I've failed to find any tutorials online... If anyone knows of such tutorials, please let me know. I did find this one: http://www.caelinux.org/wiki/index.php/Doc:CAETutorials
Which was really cool since I want to use the mesh in OpenFOAM. My problem arises because of the specific nature of the tutorial.

I basically want to generate a mesh (of some fixed height) over a digital elevation model (DEM). The DEM is made up of a matrix of points (in essence) and so to test my ability to simulate fluid (in my case air), I fabricated points for the DEM. Here's what I've done in Salome:

Once I had a grid of points (on the floor of the DEM and parallel above the floor to form a 'roof') I used the New Entity > Blocks > Quadrangle Face option to join four adjacent points together to form tiles until eventually I had a closed object.

I used the Boolean > Fuse option to combine the quadrangle shapes together (which took allot of effort since you can only fuse two objects together at a time?) and finally the New Entity > Group > Create option to define walls.​

Am I going about this the wrong way? When I try to compute the mesh I end up with an error
I'm hoping there's some with experience in Salome (I'm using 4.1.4) that has some suggestions/hints for me? Hmmm, I wonder if I should start this as a new topic

 

[minger]

I'd assume that the mesh representing the object should as close as possible of a match to the geometry of the object;

Yes

that is, for my situation (where a digital elevation model is used), the larger the gap between the generated mesh and the actual topography of the land, the less accurate the results might be? Though I suppose this too could be negligible considering the many other factors influence the flow of air over land (such as terrain cover - grass, forest, etc) and such properties will not be accurately represented.

There should be no gap 'ideally' between your mesh and physical model. When considering terrain, I would assume that the cover is part of the boundary conditions. When we model something such as pipe flow, there is inherent roughness in the pipes. Of course we can't actually model that roughness, our grid spacing would be thousandths of an inch. So, we model some nominal diameter of the pipe and then the pipe roughness gets factored into the boundary conditions.

I suppose what I really asking is, how important is grid refinement? I hope I'm not misunderstand grid refinement > I think it's the process of making the grid match the geometry with greater precision?

Grid refinment is very important. Remember that you will end with a discrete set of data. That means that you don't really know what's happening between grid points. So, if something is change a lot in an area with little gridpoints, then you aren't "resolving" what's there. Imagine trying to draw a sine wave with only three points. Imagine how it would look.

A quick question, do all discretization methods use stencils?

Discretization is simply the process by which a continuous distribution of data is converted into a finite set of data. Because there won't be a closed form of the Navier-Stokes equations, we need to discretize in order to "solve" the equations.

is the order of accuracy directly related to which term the http://en.wikipedia.org/wiki/Taylor_series" truncation occurred? And this is different from a higher-order scheme? I think I should read more before I ask these question so let's put that on hold.

Your are correct. Using the Taylor series derived stencils, if we use more data (i.e. more grid points) we can knock more and more terms out. When we write the derivatives, there are actually an infinite number of terms in the equation. If I take two grid points (at i-1, and i+1), and subtract them from each other, and divide by the spacing, I end up with
$$\frac{u_{i+1}-u_{i-1}}{2\Delta x} = \left.\frac{du}{dx}\right\vert_i + ...$$
Now, it happens that this stencil is second order accurate. This means that the first term that shows up besides the derivative that we're looking for the second order term. The general form of the "next" term looks like:
$$\frac{\Delta x^n}{n!}\frac{d^n u}{dx^n}$$
Where n is the order. We then try and say that "hopefully" this term is small enough compared to the derivative that we can ignore it. By adding more data, we can algebraically eliminate these terms.

... the events that occur at the boundary could impact the overall fluid flow and this can't be ignored anyway.

Could and almost always certainly will.

edit: Also, I'm sorry but I'm really not familiar with the meshing package that you're using. I can't really help you with that.

 

[Pooven]

Thank you so very much for answering my questions minger! I truly understand the content so much more because of you. I especially like the way you explain things :)

Do you (or anyone else) have any Open Source recommendation to meshing solutions? I'm going to use the mesh tool provided with OpenFOAM but I'd really like to compare the results with another meshing tool.

 

[minger]

I'm not familiar with any other opensource meshing packages. When all else fails, try to do it yourself. You said that you have a "elevation" map. As a first cut, if these points were taken in a grid, try just connecting the grid into a structured grid. So, for example, let's say you have an array that is (n²,3) long: points(n²,3), so basically the person took n points, took a step over, then repeated n times, with the second indice representing x,y,z. You can do something like this:

DO i=1,n
 DO j=1,n
   ii = ii + 1
   x(i,j,1) = points(ii,1)
   y(i,j,1) = points(ii,2)
   z(i,j,1) = points(ii,3)
  END DO
END DO

!--"Fill" in the rest of the x,y locations
!--You now have your surface grid points 
vPoints = 10      !--vertical points
vMax    = 100.0  !--ceiling

DO k=1,vPoints
 x(:,:,k) = x(1,1,k)
 y(:,:,k) = y(1,1,k)
 z(:,:,k) = z(1,1,k)
END DO

!--evenly space grid points above surface points
counter = 0.0
DO i=1,n
 DO j=1,n
   int = (vMax - z(i,j,1)) / vPoints
  DO k=1,vPoints
   counter = counter + int
   z(i,j,k) = counter
  END DO
 counter = 0
 END DO
END DO

!--write out as plot3d file (hopefully your open source can read)
OPEN(11,file='grid.x',form='formatted')
WRITE(11,*) '1' !--number of blocks
WRITE(11,*) n,n,vPoints !--bounds of each block
WRITE(11,*) x,y,z
CLOSE(11)

Something as simple like this might get you going. After you then see where your problem areas are, manually insert points between points. Use a second order polynomial to insert the points, which should keep the first derivative continuous at those points.