User Tools

Site Tools


docs:classes:flowfield

This is an old revision of the document!


A PCRE internal error occured. This might be caused by a faulty plugin

====== FlowField ====== The FlowField class represents vector-valued fields of the form <latex> {\bf u}({\bf x}) &= \sum_{k_x,k_y,k_z,j} \hat{u}_{k_x k_y k_z j}\bar{T}_{k_y}(y) \; e^{2 \pi i (k_x x/L_x + k_z z/L_z)} {\bf e}_j </latex> 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> \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. are stored as variables of type FlowField. The main functionality of the FlowField class is * algebraic, differential, 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>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) * 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, or assigned from computations. Examples: <code c++> FlowField f; // null value, 0-d field on 0x0x0 grid 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 h("h"); // read from file "h.ff" FlowField omega = curl(u); </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.1234888606.txt.gz · Last modified: 2009/02/17 08:36 by gibson