docs:tutorial

# Differences

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

 docs:tutorial [2008/12/22 08:35]gibson created docs:tutorial [2009/02/25 15:23]dspieker Fixed small typo that was bothering me. 2009/02/25 15:23 dspieker Fixed small typo that was bothering me.2009/02/13 07:59 gibson 2009/02/13 07:58 gibson 2009/02/13 07:58 gibson 2009/02/13 07:28 gibson 2009/02/11 09:33 gibson 2009/02/10 19:08 wikiadmin 2009/02/10 16:57 wikiadmin 2009/01/31 17:37 john_h Changed "Arnolid" to "Arnoldi2009/01/30 07:51 gibson 2009/01/29 12:39 gibson 2009/01/28 10:06 wikiadmin 2008/12/22 08:35 gibson created Next revision Previous revision 2009/02/25 15:23 dspieker Fixed small typo that was bothering me.2009/02/13 07:59 gibson 2009/02/13 07:58 gibson 2009/02/13 07:58 gibson 2009/02/13 07:28 gibson 2009/02/11 09:33 gibson 2009/02/10 19:08 wikiadmin 2009/02/10 16:57 wikiadmin 2009/01/31 17:37 john_h Changed "Arnolid" to "Arnoldi2009/01/30 07:51 gibson 2009/01/29 12:39 gibson 2009/01/28 10:06 wikiadmin 2008/12/22 08:35 gibson created Line 1: Line 1: ====== Channelflow Tutorial ====== ====== Channelflow Tutorial ====== - - ===== Intro ===== So you've installed channelflow. Now what? Well, computational fluid So you've installed channelflow. Now what? Well, computational fluid Line 19: Line 17: of my own research. Probably the best way to get started with 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 channelflow is to step through a few examples of run-of-the-mill - calculations using these utilities. If you want to get right to these + calculations using these utilities. ​ - examples, skip to Section 3. + + 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 - ===== Overview of channelflow utility programs ===== + <​code>​ + couette --help + ​ - Here's a list of current channelflow utilities. The first three - are taken out of alphabetical order because they'​re featured in - Section 3, Example Calculations. - ^ program name ^ purpose ^ + ====== Example Calculations ====== - | randomfield ​  | build a random initial velocity field, save to disk | + ===== Making ​a movie ===== - | couette ​      | integrate an initial condition, save results to disk | + - | fieldprops ​   | print out norms, symmetries, geometrical data of a stored field | + - | makemovie ​    | extract slices of fields in order to make a movie | + - | addfields ​    | compute sum a_n u_n and store result to disk | + - | arnoldi ​      | compute the eigenvalues and eigenfunctions of eqbs and orbits | + - | ascii2field ​  | convert a file of ASCII data to a channelflow FlowField | + - | changegrid ​   | change the discretization or box size of a field | + - | field2ascii ​  | convert a channelflow FlowField to a file of ASCII data | + - | fieldplots ​   | extract a number of 2D slices of the 3D field, good for plots | + - | findorbit ​    | compute an equilibrium or periodic orbit of plane Couette | + - | L2Dist ​       | compute the L2 distance between two fields | + - | L2IP          | compute the L2 inner product | + - | makebasis ​    | construct an orthonormal basis from a set of fields | + - | makeheatmode ​ | construct a field that decays in time according to Laplace eqn | + - | makestokesmode | construct a stokes eigenfunction of laminar equilibrium | + - | perturbfield ​ | add random perturbations to a given field | + - | projectfields | project a set of fields onto a given basis | + - | projectseries | project a sequence of fields onto a given basis | + - | seriesprops ​  | compute statistics on a sequence of data | + - | symmetrize ​   | find the phase shift of a field that optimizes symmetries | + - The utilities are stand-alone command-line programs that are run from + ===1. Generate an initial condition and examine its properties === - the Unix shell. You can get brief built-in help information on each + - utility by running it with a -h or --help option. For example, running + - "​couette --help"​ produces + - + - gibson@akbar$couette --help + - couette : + - integrate an initial condition and save velocity fields to disk. + - + - options : + - -T0 ​--T0 ​ <​real> ​ ​default ​== 0 start time + - -T1 ​--T1 ​ <​real> ​ ​default == 100 end time + - -vdt --variabledt ​ ​adjust dt for CFL + - -dt ​--dt ​ <​real> ​ ​default == 0.03125 timestep + - -dtmin ​ --dtmin ​ <​real> ​ ​default == 0.001 ​minimum time step + - -dtmax ​ --dtmax ​ <​real> ​ ​default == 0.05 maximum time step + - -dT ​--dT ​ <​real> ​ ​default ​== 1 save interval + - -CFLmin ​ ​--CFLmin ​ <​real> ​ ​default == 0.4 ​minimum CFL number + - -CFLmax ​ ​--CFLmax ​ <​real> ​ ​default == 0.6 ​maximum CFL number + - -ts ​--timestepping <​string>​ default == sbdf3 ​timestepping algorithm + - ... + - ​-p ​ --pressure ​ print pressure grad + - <​flowfield> ​ (trailing arg 1) initial condition + - + - The built-in help gives a brief description of each utility'​s purpose + - + - and a list of its command-line options and arguments. Channelflow + - utilities are invoked at the command line with syntax like + - + - ​utility -opt1 value -opt2 value -flag1 arg3 arg2 arg1 + - + - or concretely + - + - ​couette -T0 0 -T1 -vdt -dt 0.02 -ts sbdf4 u0.ff + - + - + - "​Options"​ (e.g. -opt1 value) are used to reset default values + - of parameters. For options, the first two columns in the built-in + - help give the short and long form of the option (e.g. -ts and + - --timestepping),​ the third column indicates the type of parameter + - expected (e.g. real, int, bool, string), and the fourth gives the + - the default value. For example, "​couette -dt 0.02 -ts cnab2" sets + - the time stepping method to 2nd order Crank-Nicolson Adams-Bashforth + - with dt=0.02. + - + - "​Flags"​ simply turn on boolean options that would otherwise be set + - to false. For example, calling "​couette -vdt" turns on variable-dt + - timestepping,​ which adjusts dt at fixed intervals to keep the CFL + - number within bounds. For flags the third and fourth columns of + - built-in help are left blank. + - + - "​Arguments"​ always come after all options and flags. Arguments usually + - specify the filenames of binary velocity fields that the utility will + - load and operate on. Most channelflow programs have one required + - argument (e.g. "​couette u0.ff"​) some two (e.g. "​L2Dist u0.ff u2.ff"​). + - Others take a variable number of arguments (e.g. makebasis u0 u1 u2"​). + - Unfortunately it's difficult to document variable-number arguments + - properly in the four-column option system, so variable-number arguments + - are usually documented with a "​usage:​ line right after the description + - of the utility'​s purpose. + - + - So, as you read work through the Example Calculations,​ you can run the + - suggested command with a --help option to clarify what the options are + - doing and what other options are possible. + - + - + - ===== Example Calculations ===== + - + - + - ==== How to make a movie of a flow ==== + - + - + - === Generate an initial condition and examine its properties === + gibson@akbar$ randomfield -Nx 48 -Ny 35 -Nz 48 -lx 0.875 -lz 0.6 -m 0.20 u0.ff gibson@akbar$randomfield -Nx 48 -Ny 35 -Nz 48 -lx 0.875 -lz 0.6 -m 0.20 u0.ff Line 146: Line 54: - === Integrate a flow in time, saving the results to disk === + === 2. Integrate a flow in time, saving the results to disk === gibson@akbar$ couette -T0 0 -T1 200 -l2 -o data u0.ff gibson@akbar$couette -T0 0 -T1 200 -l2 -o data u0.ff Line 163: Line 71: - === Extract data from the sequence of stored velocity fields for plotting === + === 3. Extract data from the sequence of stored velocity fields for plotting === gibson@akbar$ movieframes -T0 0 -T1 200 -d data -o frames gibson@akbar$movieframes -T0 0 -T1 200 -d data -o frames Line 169: Line 77: The movieframes program reads in the series of files data/u0.ff, 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 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 that are good for visualizing the flow. These 2D slices are stored in the frames/ directory with filenames like u0_yz_slice.asc. the frames/ directory with filenames like u0_yz_slice.asc. - + === 4. Make a movie from the extracted data === - === Make a movie from the extracted data === + To make a movie using channelflow'​s existing visualization tools, you To make a movie using channelflow'​s existing visualization tools, you Line 197: Line 103: script and its arguments; briefly, here the arguments are script and its arguments; briefly, here the arguments are - 0 starting frame number + 0 starting frame number ​ int - 1 frame interval + 1 frame interval ​ int - 200 ending frame number + 200 ending frame number ​ int - 0 starting time (t=0.0) + 0 starting time float - 1 time interval ​(dT=1.0) + 1 time interval ​ float - 10 ​frames per second + 10 ​frames per second ​ ​int ​ - '​couette.avi'​ output filename + '​couette.avi'​ output filename ​ string - I have a number of scripts that will convert these AVI animations to + further optional arguments ​are - other video formats and compress them along the way. These conversion + - scripts ​are not yet part of the channelflow distribution package, but + - I will try to work them in in the near future. Sorry for the delay, + - but these things take time. + - ==== Project movie data onto state-space coordinates ==== + 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 + + <​latex>​ + 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 It can be useful to look at the temporal evolution of a fluid as Line 221: Line 250: equilibria under the symmetries of plane Couette flow. In simple equilibria under the symmetries of plane Couette flow. In simple language, we take linear combinations of equilibria and their language, we take linear combinations of equilibria and their - translations in x,z to form orthonormal basis sets. + 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. - === Make a low-d basis === + === 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 Download an equilibrium solution of plane Couette flow from the channelflow website, one that is compatible in geometry and channelflow website, one that is compatible in geometry and - discretization. + 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 τ<​sub>​x​ be translation by Lx/2, etc. + Then the above lines compute τ<​sub>​x​ UB, τ<​sub>​z,​ and + τ<​sub>​xz​ respectively. + + Now construct the following orthogonal linear combinations of the above fields + + <​latex>​ \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 + + <​latex>​ + (f,g) = 1/V \int_V f \cdot g dx dy dz + ​ - (to be continued...) + ​ - + ​