====== Channelflow Tutorial ====== So you've installed channelflow. Now what? Well, computational fluid dynamics is a pretty big field, and there's no telling what kinds of ideas you might want to explore. For this reason Channelflow was developed primarily as a *programming language*. If you're using Channelflow for research towards a Ph.D. thesis, you will probably eventually end up writing your own channelflow programs, for example, modifying existing programs to integrate flows and analyze data. Or you might need to modify the time-integration routines to incorporate additional physics (polymer additives, bubbles, etc.) However, channelflow also includes a number of predefined utility programs that a basic set of important computations, such as time-integration of plane Couette and channel flows and measuring properties of velocity fields. These utilities suffice for the bulk of my own research. Probably the best way to get started with channelflow is to step through a few examples of run-of-the-mill calculations using these utilities. Please refer to [[:docs#utilities|Utilities]] and [[docs:utils:options|Utility Options]] for a detailed guide of individual utilities and their options. You can also run any utility with a ''-h'' or ''--help'' option to get a brief description of the utility's purpose and options, e.g couette --help ====== Example Calculations ====== ===== Making a movie ===== ===1. Generate an initial condition and examine its properties === gibson@akbar$ randomfield -Nx 48 -Ny 33 -Nz 48 -lx 2 -lz 1 -m 0.20 u0.ff This command generates a no-slip, divergence-free velocity field with random spectral coefficients on a 48 x 33 x 48 grid, on [0, 2pi] x [-1, 1] x [0, pi], with magnitude 1/V \integral_V |u|^2 dx = 0.2. The field is a perturbation from the laminar flow --by default, velocity fields in channelflow are differences from laminar. The spectral coefficients are random within an exponentially decaying envelope, roughly similar to turbulent fields. The velocity field is saved to disk in binary format in file u0.ff (the .ff stands for "FlowField", the name of the C++ class for velocity, pressure and tensor fields in channelflow). The channelflow binary format stores the spectral coefficients, the geometry, and all discretization information so the field can be reconstructed entirely from the data in the file. You can list this information and some dynamical properties of the field by running gibson@akbar$ fieldprops u0.ff === 2. Integrate a flow in time, saving the results to disk === gibson@akbar$ couette -T0 0 -T1 200 -l2 -o data u0.ff This command load the velocity field u0.ff from disk and integrates it in time (using the default integration algorithm and parameters) from t=T0=0 to t=T1=200, and saves the time varying velocity field to disk at the interval dT=1.0 (the default save interval) in a directory named data/. The -l2 option prints out the L2 norm of u as well as the Chebyshev-weighted norm. After this command finishes, look in the data/ directory, and you will see u0.ff u1.ff u2.ff etc. The integer label is the time (remember the save interval is dT=1.0). If you choose a noninteger save interval, the filenames will be something like u0.000.ff u0.975.ff etc. === 3. Extract data from the sequence of stored velocity fields for plotting === gibson@akbar$ movieframes -T0 0 -T1 200 -d data -o frames The movieframes program reads in the series of files data/u0.ff, data/u1.ff, etc. and extracts a number of 2D slices of the 3D fields that are good for visualizing the flow. These 2D slices are stored in the frames/ directory with filenames like u0_yz_slice.asc. === 4. Make a movie from the extracted data === To make a movie using channelflow's existing visualization tools, you need Matlab. (If you would like to write similar tools for another visualization package, please do so and send them to me!). Start up Matlab. Get all the scripts in channelflow/matlab into your Matlab path. Do this either by copying the scripts into the current directory, by copying them to wherever you store your Matlab scripts, or by putting channelflow/matlab in your Matlab path. You can do that by typing 'addpath /home/larry/channelflow-1.3.2/matlab' at the Matlab prompt (changing the directory as appropriate). Then, within Matlab change to the directory that where you ran the couette programs, the one with data/ and frames/ subdirectories. Within Matlab run makemovie(0,1,200,0,1,10, 'couette.avi'); This will construct a movie of the 3D velocity field as it evolves in time and store the result in file couette.avi, in AVI format. Running 'help makemovie' will give you a help string about the makemovie script and its arguments; briefly, here the arguments are 0 starting frame number int 1 frame interval int 200 ending frame number int 0 starting time float 1 time interval float 10 frames per second int 'couette.avi' output filename string further optional arguments are title printed this in the upper-left corner string credit printed this in the lower-right corner string xstride plot every xstride-th gridpoint int ystride plot every ystride-th gridpoint int zstride plot every zstride-th gridpoint int perspect do a perspective plot 0 or 1 (false or true) framedir directory containing frame data string (default='frames') Note: The Matlab scripts provided with channelflow are kludgy. I cobbled them together in order to get the plots I want. Some things, like the position of the title and credit strings, must be positioned manually by editing values in the script files. Improvements in the scripts and alternatives for systems other than matlab are welcome. === 5. Convert the AVI file === Matlab produces only uncompressed AVI files on Linux. You will probably want to compress the AVI file and convert it to another format. On Linux you can do this with ''mencoder'', which is part of the MPlayer package. For example, this command will convert ''couette.avi'' file to a flash video file ''couette.flv''. gibson@akbar$ mencoder couette.avi -nosound -of lavf -lavopts format=flv -ovc lavc -lavcopts vcodec=flv:vmax_b_frames=0:vbitrate=1600 -o couette.flv Adjust the bitrate to balance filesize and video quality. ===== Computing a 1d unstable manifold ===== The Nagata (1990) "lower-branch" equilibrium has a one-dimensional unstable manifold. Here we compute the unstable manifold by integrating two 1d trajectories u_{\pm}(x,t) = f^t(u_{LB} \pm v_{LB}), t \in [0, \infty] using several channelflow utilities: * ''fieldprops'' * ''arnoldi'' * ''addfields'' * ''couette'' * ''seriesprops'' * ''makebasis'' * ''projectseries'' === 1. Download the Nagata lower-branch solution === ...from the [[http://www.channelflow.org/database|channelflow database]]. ''LB'' stands for 'lower-branch'. wget http://channelflow.org/database/a1.14_g2.5_Re400/LB.ff === 2. Examine the solution's properties === The ''fieldprops'' utility will print out basic information about the field. For example, fieldprops -g LB prints out the field's geometrical properties: cell size, grid size, etc. Try ''--help'' to get a list of other options. Channelflow adds a ''.ff'' file extension to ''LB'' if you leave it off. === 3. Plot the solution === Visualization of fluid velocity fields is an art in itself. Channelflow provides a few scripts for plotting the velocity field on certain slices of the rectangular domain. I've found these plots useful, but if you have better ideas please adapt the scripts accordingly. Plotting take two steps. First you extract some 2D slices from the 3D field with a channelflow utility, like this fieldplots -o plot LB That saves the 2D slices as ASCII data files in a plot/ directory. Then within Matlab, go to the plot/ driectory and run plotbox('LB') The matlab ''plotbox'' script has a number of default parameters that you can change. Try ''help plotbox'' within Matlab for more information. === 4. Compute the eigenfunctions === The Nagata lower-branch solution is an equilibrium of plane Couette dynamics. You can compute the eigenvalues and eigenfunctions of the linearized dynamics about the equilbrium with the ''arnoldi'' utility. (Will write documentation on Arnoldi iteration later). arnoldi --flow LB.ff This produces a set of (approximate) eigenfunctions ''ef1.ff, ef2.ff, ...'' and a file of eigenvalues ''lambda.asc''. === 5. Perturb along the unstable manifold === The Nagata lower branch has a single unstable eigenvalue, so its unstable manifold is 1d and can be computed as a trajectory initiated with small perturbations in the +/- directions of the unstable eigenvector/eigenfunction. The following calculates LB +/- 0.01 ef1 and saves the results into files LBp01ef1 and LBm01ef1 addfields 1 LB 0.01 ef1 LBp01ef1 addfields 1 LB -0.01 ef1 LBm01ef1 === 6. Integrate the perturbations === couette -T0 0 -T1 400 -o data-LBp01 LBp01ef1 couette -T0 0 -T1 400 -o data-LBm01 LBm01ef1 === 7. Produce input vs dissipation curves === The ''seriesprops'' utility computes a few quantities like energy dissipation D and wall shear I for a time series of stored velocity fields seriesprops -T0 0 -T1 400 -d data-LBp01ef1 -o props-LBp01ef1 seriesprops -T0 0 -T1 400 -d data-LBm01ef1 -o props-LBm01ef1 The results will be stored in files in props-LBp01ef1/ and props-LBm01ef1/ directories === 8. Make movies === movieframes -T0 0 -T1 100 -d data-LBp01ef1 -o frames-LBp01ef1 movieframes -T0 0 -T1 100 -d data-LBm01ef1 -o frames-LBm01ef1 From here you can adapt the [[#make_a_movie_from_extracted_data|movie-making instructions]] from above. ===== Project movie data onto state-space coordinates ===== It can be useful to look at the temporal evolution of a fluid as a trajectory in state space. The number of degrees of freedom in a fluid simulation is very high (e.g. 10^5), so it is necessary to project the fields into a low-dimensional basis in order to plot the trajectory and look at it. We have found that good projection bases can be constructed from the "group orbits" of equilibria under the symmetries of plane Couette flow. In simple language, we take linear combinations of equilibria and their translations in x,z to form orthonormal basis sets. For a more detailed description of the logic and mathematics of this approach, see [[:references|Gibson et al (2007) JFM 611]]. Here we will just outline how the computation is done using channelflow. === 1. Make a low-d basis === Make a subdirectory and descend into it, so that the following steps don't pollute the current directory with a bunch of extraneous files mkdir basis-UBtrans cd basis-UBtrans Download an equilibrium solution of plane Couette flow from the channelflow website, one that is compatible in geometry and discretization. For example, you can get the Nagata upper-branch equilibrium (UB) with the Unix "wget" utility. wget http://www.channelflow.org/database/a1.14_g2.5_Re400/UB.ff Compute the half-cell translations of UB in x, in z, and in x,z with the channelflow [[:docs:utils:symmetryop]] utility: symmetryop -ax 0.5 UB UBx symmetryop -az 0.5 UB UBz symmetryop -ax 0.5 -az 0.5 UB UBxz Briefly, symmetryop constructs a symmetry σ parameterized by the options, applies it to the first FlowField argument, and saves the result to the second FlowField argument, according to the symmetry parameterization described in [[:docs:math:symmetry]]. Let if τx be translation by Lx/2, etc. Then the above lines compute τx UB, τz, and τxz respectively. Now construct the following orthogonal linear combinations of the above fields $ \begin{align*} UB_{pppp} = UB + \tau_x UB + \tau_z UB + \tau_{xz} UB \\ UB_{ppmm} = UB + \tau_x UB - \tau_z UB - \tau_{xz} UB \\ UB_{pmpm} = UB - \tau_x UB + \tau_z UB - \tau_{xz} UB \\ UB_{pmmp} = UB - \tau_x UB - \tau_z UB + \tau_{xz} UB \end{align*} $ with the channelflow [[:docs:utils:addfields]] utility: addfields 1 UB 1 UBx 1 UBz 1 UBxz UBpppp addfields 1 UB 1 UBx -1 UBz -1 UBxz UBppmm addfields 1 UB -1 UBx 1 UBz -1 UBxz UBpmpm addfields 1 UB -1 UBx -1 UBz 1 UBxz UBpmmp Finally, use the channelflow [[:docs:utils:makebasis]] utility to apply Gram-Schmidt orthogonalization on those fields and form an orthonormal basis set: makebasis UBpppp UBppmm UBpmpm UBpmmp The output of "makebasis" will be four orthonormal basis elements e0.ff, e1.ff, e2.ff, and e3.ff saved to disk. In this case the input fields are already orthogonal and all "makebasis" does is normalize. Now pop out of the basis-UBtrans subdirectory cd .. === 2. Project a series of fields onto the basis === Ok. Suppose you have a series of velocity fields u0.ff, u1.ff, etc for t=0,1,2,...1000 in a data/ directory and a set of basis elements e0.ff, e1.ff, e2.ff, e3.ff in a basis-UBtrans/ directory. To project the fields onto the basis, run projectseries -T0 0 -T1 1000 -d data -b basis-UBtrans -Nb 4 -o a.asc That will produce an ASCII file a.asc with 4 columns and 1001 rows. The t-th row and jth column is the value of (u(t), ej), where ( , ) signifies the L2 inner product (f,g) = 1/V \int_V f \cdot g dx dy dz