User Tools

Site Tools


docs:classes:flowfield

Differences

This shows you the differences between two versions of the page.

Link to this comparison view

Both sides previous revision Previous revision
Next revision
Previous revision
docs:classes:flowfield [2009/02/16 10:52]
gibson
docs:classes:flowfield [2010/02/02 07:55] (current)
Line 9: Line 9:
 and also scalar and tensor fields with appropriate changes in the dimensionality of the coefficients. ​ and also scalar and tensor fields with appropriate changes in the dimensionality of the coefficients. ​
 The barred %%T%% function is a Chebyshev polynomial scaled to fit the domain y ∈ [a,b]. ((<​latex>​ The barred %%T%% function is a Chebyshev polynomial scaled to fit the domain y ∈ [a,b]. ((<​latex>​
-\bar{T}_{k_y}(y) = T_{k_y}\left(\frac{2}{b-a}(y - \frac{b+a}{2})\right)</​latex>​))+\bar{T}_{k_y}(y) = T_{k_y}\left(\frac{2}{b-a}(y - \frac{b+a}{2})\right)</​latex>​)) ​The spatial domain  
 +of a FlowField is  Ω = [0,Lx] x [a,b] x [0,Lz], with periodicity in x and z.
  
 In channelflow programming,​ fields such as velocity, pressure, stress tensors, vorticity, etc.  In channelflow programming,​ fields such as velocity, pressure, stress tensors, vorticity, etc. 
 are stored as variables of type FlowField. The main functionality of the FlowField class is are stored as variables of type FlowField. The main functionality of the FlowField class is
  
-   * algebraic ​and differential operations, +/-, +=, ∇, ∇<​sup>​2</​sup>,​ etc. +   * algebraicdifferential, and symmetry ​operations, +/-, +=, ∇, ∇<​sup>​2</​sup>​, norms, inner products, etc. 
-   * transforming back and forth between spectral coefficients <​latex>​ \hat{u}_{k_x k_y k_z j}</​latex>​ and gridpoint values <​latex>​{\bf u}(x_jy_k, z_{\ell})</​latex>​+   * transforming back and forth between spectral coefficients <​latex>​ \hat{u}_{k_x k_y k_z j}</​latex>​ and gridpoint values <​latex>​u_j (x_{n_x}, y_{n_y}, z_{n_z})</​latex>​
    * serving as input to DNS algorithms, which map velocity fields forward in time: u(x,t) → u(x, t+Δt)    * serving as input to DNS algorithms, which map velocity fields forward in time: u(x,t) → u(x, t+Δt)
 +   * setting and accessing scpetral coefficients and gridpoint values
 +   * reading and writing to disk
 +
 +For a complete description of FlowField functionality,​ see the header file {{:​librarycode:​flowfield.h}}.
 +
 +===== Constructors / Initialization =====
  
 FlowFields are initialized with gridsize and cellsize parameters, read from disk,  FlowFields are initialized with gridsize and cellsize parameters, read from disk, 
 or assigned from computations. Examples: or assigned from computations. Examples:
  
-<code c++>+<code c++> ​
    ​FlowField f;                                   // null value, 0-d field on 0x0x0 grid     ​FlowField f;                                   // null value, 0-d field on 0x0x0 grid 
-   ​FlowField u(Nx, Ny, Nz, Nd, Lx, Lz, a, b);     // Nd-dim field on Nx x Ny x Nz grid, [0,Lz]x[a,​b]x[0,​Lz]+   FlowFIeld f(g);                                // make a copy of g 
 +   FlowField u(Nx, Ny, Nz, Nd, Lx, Lz, a, b);     // Nd-dim field on Nx x Ny x Nz grid, [0,Lx]x[a,​b]x[0,​Lz]
    ​FlowField g(Nx, Ny, Nz, Nd, 2, Lx, Lz, a, b);  // Nd-dim 2-tensor  ​    ​FlowField g(Nx, Ny, Nz, Nd, 2, Lx, Lz, a, b);  // Nd-dim 2-tensor  ​
    ​FlowField h("​h"​); ​                             // read from file "​h.ff"​    ​FlowField h("​h"​); ​                             // read from file "​h.ff"​
-   FlogField ​omega = curl(u);+   FlowField ​omega = curl(u);
 </​code>​ </​code>​
 +===== Algebraic and differential operators =====
  
 +Assume f,g,h etc. are FlowField variables with compatible cell and grid sizes. Examples of 
 +possible operations
 +
 +<​code>​
 +   f *= 2.7;                // f = 2.7*f
 +   f += g;                  // f = f + g
 +   f -= g;                  // f = f - g
 +
 +   f = xdiff(g); ​           // f_i  = d g_i/dx
 +   f = ydiff(g); ​           // f_i  = d g_i/dy
 +   f = zdiff(g); ​           // f_i  = d g_i/dz
 +   f = diff(g, j, n);       // f_i  = d^n g_i/dx_j
 +   f = diff(g, j, n);       // f_i  = d^n g_i/dx_j
 +   f = grad(g); ​            // f_ij = dg_i/dx_j
 +   f = curl(g); ​
 +   f = lapl(g);
 +   f = div(g);
 +   f = cross(g,​h); ​     ​
 +  ​
 +   ​xdiff(g,​ dgdx); ​         // same as dgdx = xdiff(g), but often more efficient
 +   ​curl(g,​ curl_g); ​        // ditto
 +   ​lapl(g,​ lapl_g); ​        // ditto
 +   ...
 +
 +   Real c = L2IP(f,​g); ​     // L2 inner product of f,g
 +   Real n = L2Norm(u);
 +   Real D = dissipation(u);​
 +   Real E = energy(u); ​
 +   Real I = wallshear(u);​
 +</​code>​
 +
 +The latter functions are defined as 
 +
 +<​latex>​ $ \begin{align*}
 +L2IP(f,g) &= \frac{1}{L_x L_y L_z} \int_{\Omega} {\bf f} \cdot {\bf g} \,\, d{\bf x} \\
 +L2Norm(u) &= \left(\frac{1}{L_x L_y L_z} \int_{\Omega} \|{\bf u}\|^2\,\, d{\bf x} \right)^{1/​2}\\
 +     E(u) &= \frac{1}{2 L_x L_y L_z} \int_{\Omega} \|{\bf u}\|^2\,\, d{\bf x}  \\
 +     D(u) &= \frac{1}{L_x L_y L_z} \int_{\Omega} \|\nabla \times {\bf u}\|^2 \,\, d{\bf x} \\
 +     I(u) &= \frac{1}{L_x L_z} \int_{y=a,​b} \frac{\partial u}{\partial y} \, dx dz
 +\end{align*} $ </​latex>​
 +
 +Note that expressions such as %%f = g + h%% or %%f = 0.5*(g + h)%% are **not allowed**
 +on FlowFields, since these would generate temporary FlowField variables %%g + h%% and
 +%%0.5*(g+h)%% during expression evaluation. Instead, use sequences such as 
 +
 +<​code>​
 +  // A sequence that results in f = 0.5*(g + h);
 +  f = g;
 +  f += h;
 +  f *= 0.5;
 +</​code>​
 +  ​
 +As C++ objects, FlowFields are huge monsters. It is best to minimize the amount
 +of construction,​ copying, assignment of FlowFields by reusing temporaries and
 +figuring out the minimal sequence of operations to get the desired result.
 +
 +===== Symmetry operations =====
 +
 +The symmetry group of FlowFields is represented by the [[fieldsymmetry|FieldSymmetry]]
 +class. Briefly, the symmetries of 3D FlowFields can be parameterized as 
 +
 +<​latex>​ $ \begin{align*}
 +  \sigma &= (s_x, s_y, s_x, a_x, a_z, s)\\
 +  s_x, s_y, s_z, s &= \pm 1\\
 +  a_x, a_z &\in [-0.5, 0.5)
 +\end{align*} $ </​latex>​
 +
 +with the action of σ on a velocity field u as
 +
 +<​latex>​
 +  \sigma [u, v, w](x,y,z) = s (s_x u, s_y v, s_z w)(s_x x + a_x L_x, s_y y, s_z z + a_z L_z)
 +</​latex>​
 +
 +A FieldSymmetry can be constructed and applied to a FlowField as follows
 +
 +<code c++>
 + ​FieldSymmetry sigma(sx, sy, sz, ax, az, s);  // construct sigma = (sx,​sy,​sz,​ax,​az,​s)
 + ​FlowField sigma_u = sigma(u); ​               // apply symmetry sigma to u
 +</​code>​
 +
 +Or, the symmetric component of a field can be obtained by 
 +
 +<code c++>
 + ​FlowField Pu = u;  ​
 + Pu += sigma(u); ​   // Pu now equals u + sigma u
 + Pu *= 0.5;         // Pu now equals (u + sigma u)/2
 +</​code>​
 +
 +For more examples of FlowField and FieldSymmetry usages, see 
 +[[:​docs:​classes:​fieldsymmetry|the FieldSymmetry documentation]].
 +
 +===== Transforms and data access =====
 +
 +FlowField transforms are a complicated subject --there are transforms in x,y, and z and 
 +implicit symmetries in complex spectral coefficients due to the real-valuedness of the 
 +field, for instance. This section outlines the bare essentials of transforms and data
 +access methods. For further details see the {{docs:​chflowguide.pdf|Channelflow User Guide}}.
 +
 +The FlowField class has a large data array that contains the spectral coefficients of
 +the expansion listed at the top of this page. Most operations on FlowFields are 
 +calculated in terms of these spectral coefficients. But sometimes we need to know 
 +the value of the field at gridpoints. Rather than directly evaluating the expansion sum
 +for given values of (x,y,z) (which would be very slow), we use fast Fourier transforms
 +to transform the array of spectral coefficients into an array of gridpoint values. ​
 +The main FlowField class transform methods are
 +
 +<​code>​
 +   ​u.makePhysical();​ // transform spectral coeffs to gridpt values
 +   ​u.makeSpectral();​ // transform gridpt values to spectral coeffs
 +</​code>​
 +
 +Because the transforms change the meaning of the FlowField'​s internal data 
 +array, ***you need to make sure the FlowField is in the proper state before ​
 +trying to access either its spectral coefficients or its gridpoint values.**
 +
 +For example, to print out the entire set of gridpoint values of a FlowField,
 +you would want to make the field Physical first.
 +<​code>​
 +   char s = ' ';
 +   ​u.makePhysical();​
 +   for (int i=0; i<​u.Nd();​ ++i)
 +     for (int nx=0; nx<​u.Nx();​ ++nx)
 +        for (int ny=0; ny<​u.Ny();​ ++ny)
 +           for (int nz=0; nz<​u.Nz();​ ++nz)
 +              cout << nx <<​s<<​ ny <<​s<<​ nz <<​s<<​ u(nx,​ny,​nz,​i) << endl;
 +</​code>​
 +
 +To print out its spectral coefficients,​ you need to make it
 +Spectral first 
 +<​code>​
 +   ​u.makeSpectral();​
 +   for (int i=0; i<​u.Nd();​ ++i)
 +     for (int mx=0; mx<​u.Mx();​ ++mx) {
 +        int kx = u.kx(mx);
 +        for (int my=0; my<​u.My();​ ++my)
 +           for (int mz=0; mz<​u.Mz();​ ++mz) {
 +              int kz = u.kz(mz);
 +              cout << kx <<​s<<​ my <<​s<<​ kz <<​s<<​ u.cmplx(mx,​my,​mz,​i) << endl;
 +           }
 +     }
 +</​code>​
 +Note that the loop variables for mx,mz are //not// the same as the wavenumbers kx,​kz. ​
 +That's because the Fourier transforms leave the data in a peculiar order. Channelflow
 +tries to ease the pain of this difference by providing functions %%int kx = u.kx(mx)%%
 +and %%int mx = u.mx(kx)%% that translate between data ordering %%mx%% and Fourier ​
 +wavenumbers %%kx%%, and similarly for %%mz,kz%%. ((We could eliminate the issue entirely, ​
 +but at the cost of run-time efficiency)).
 +Note also that the data access method for spectral coefficients is the complex-valued ​
 +%%u.cmplx(mx,​my,​mz,​i)%%,​ compared to the real-valued gridpoint access method %%u(nx,​ny,​nz,​i)%%
 +and that the bounds of the indexing variables are different. ​
 +
 +If you really want to loop in %%kx,kz%% order (at the slight cost in efficiency),​ do this
 +<​code>​
 +   ​u.makeSpectral();​
 +   for (int i=0; i<​u.Nd();​ ++i)
 +     for (int kx=u.kxmin();​ kx<​u.kxmax();​ ++kx) {
 +        int mx = u.mx(kx);
 +        for (int my=0; my<​u.My();​ ++my)
 +           for (int kz=u.kzmin();​ kz<​u.kzmax();​ ++kz) {
 +              int mz = u.mz(kz);
 +              cout << kx <<​s<<​ my <<​s<<​ kz <<​s<<​ u.cmplx(mx,​my,​mz,​i) << endl;
 +           }
 +     }
 +</​code>​
 +
 +But in general it is better to use built-in FlowField operations such as %%curl%% and %%diff%%
 +than to loop over the data arrays, if you can.
 +
 +You can also perform the $x,z$ and the $y$ transforms independently. For example, if 
 +%%u%% is representing pure gridpoint values you could do this
 +
 +<code c++>
 +  // get a gridpoint value
 +  Real u_nxnynzi = u(nx,​ny,​nz,​i); ​
 +
 +  u.makeSpectral_xz();​
 +
 +  // get kx,kz Fourier coefficient at ny-th gridpoint in y
 +  Complex ukxnykzi = u.cmplx(u.mx(kx),​ ny, u.mz(kx), i) 
 +</​code>​
 +
 +The complete set of such transform functions is 
 +
 +<code c++>
 +  u.makeSpectral(); ​   // to pure spectral coeffs
 +  u.makePhysical(); ​   // to pure gridpoint values
 +  u.makeSpectral_xz();​ // to spectral coeffs in x,z
 +  u.makeSpectral_y(); ​ // to spectral coeffs in y
 +  u.makePhysical_xz();​ // to gridpoint values in x,z
 +  u.makePhysical_y(); ​ // to gridpoint values in y
 +  u.makeState(Physical,​ Spectral); // to x,z Physical and y Spectral
 +  u.makeState(...,​ ...); // and the other three combinations of (Physical,​Spectral);​
 +</​code>​
 +The FlowField keeps track of its spectral/​physical states in
 +x,z and y performs the desired transform only if it's in the
 +opposite state. You can query the state of a FlowField like this
 +
 +<code c++>
 +  fieldstate xzstate = u.xzstate();​
 +  if (xzstate == Physical)
 +    ....
 +
 +  fieldstate ystate = u.ystate();
 +  if (ystate == Spectral)
 +    ....
 +</​code>​
  
 +The FlowField class has quite a few other member functions and operators.
 +For a complete description,​ see the header file {{:​librarycode:​flowfield.h}}.
docs/classes/flowfield.1234810322.txt.gz · Last modified: 2009/02/16 10:52 by gibson