Axisymmetric, hydromagnetic, modulated Couette flow between infinite or finite
cylinders in the small Prandtl number limit. Time-stepping is via 2nd order
accurate implicit Crank-Nicolson for the linear terms and 2nd order accurate
explicit Adams-Bashforth for the non-linear terms. For the diffusive equations
the code uses operator factorisation to allow a tridiagonal system to be
solved; the linear algebra package LAPACK is used to solve the Poisson
equations associated with the stream function, current and magnetic field. The
spatial discretisation is via 2nd order accurate centred finite differences.
It includes a 'homotopy' parameter,
tau, to continuously deform the
boundaries from the infinite cylinder case (
tau = 0) to the finite cylinder
tau = 1). A basic particle path subroutine to track the trajectory of
a particle is also included. Spatial modulation of the inner and/or outer
Reynolds numbers in the axial direction is possible to mimic a wavy cylinder
- LAPACK - the linear algebra package.
- BLAS - the basic linear algebra subprograms.
- A Fortran 90/95 compiler.
- couette_mod.f90 - Main program file.
- current.f90 - Routines to solve the Poisson equation for the azimuthal current.
- derivs.f90 - Routines to calculate derivatives.
- ic_bc.f90 - Routines to set up initial and boundary conditions.
- io.f90 - Routines to do with input/output.
- linear.f90 - Routines to set up the linear parts of the RHS of the azimuthal velocity and vorticity fields.
- magnetic.f90 - Routines to solve the Poisson equation for the azimuthal magnetic field.
- Makefile - Makefile to compile the code (see below).
- matrices.f90 - Routines to set up matrices involved in time-stepping and solving Poisson equations.
- nonlinear.f90 - As linear.f90 but for the nonlinear terms.
- parameters.f90 - Parameters to set (see below).
- README.md - This file.
- setup - Setup script to compile and set up ready to run in separate directory.
- solve.f90 - Routines to solve for the azimuthal velocity and vorticity fields, as well as the Thomas algorithm.
- stream.f90 - Routines to solve Poisson equation for stream function.
- variables.f90 - Derived types and routines to do with variables.
alpha - Wavenumber (
2pi/wavelength) in the infinite case. Set equal zero if using finite cylinders.
gamma - Aspect ratio (
height/gap) for the finite case. Set equal to
2pi/alphaif using infinite cylinders.
eta - Radius ratio (
- Q - Chandrasekhar number. Measure of the imposed magnetic field.
Re1,2 - Reynolds number in
- Re1,2_mod - Modulation amplitude as above.
- om1,2 - Frequency of modulation as above.
Re_incr - How much
Re1should be incremented or decremented when searching for critical values.
growth_tol - How close two successive growth rates should be before
Re1is altered in searching for critical values.
dt - Time step. In general
10^-4is a good choice. Once the Reynolds numbers are large (>700, say) and/or the spatial resolution is increased significantly (>80 points in z) then
dt = 10^-5is a better choice for stability.
seed - Initial seed for initial conditions. Set small
O(10^-10)for calculating linear growth rates. Set
O(10^-3)for non-linear saturation. In practice this can be set to zero, since the boundary conditions of the azimuthal velocity at the cylinder wall(s) can start the flow.
- end_time - Final viscous diffusion time.
tau_init - Initial value of homotopy parameter, tau.
tau = 0=> infinite cylinders,
tau = 1=> finite cylinders.
- tau_step - For steady case, how much tau should be increased after saturation at each tau.
- tau_end - For steady case, the final value of tau desired.
- nx_init - Number of radial grid points.
- nt_init - Number of azimuthal grid points for use when a 3D OpenDX isosurface is desired. This is not actually used in any computation in the code other than for the isosurface plots.
nz_init - Number of axial grid points. For infinite
2*nxis sufficient. For finite
gamma*nxis usually required for full resolution.
- save_rate - After how many time steps should velocities be saved?
- save_rate_2 - After how many time steps should cross-sections be saved?
- xsect_save - Should cross-sections (surfaces) of fields be saved?
- save3d - Should a 3D isosurface be saved (OpenDX)?
iso_hel - If
save3d = .true.should the isosurface be of the helicity? If
iso_hel = .false.then the stream function is saved.
restart - Should we restart from a previous run? If
end_state.datshould be copied to restart directory.
auto_tau - For steady case, should tau be automatically stepped after saturation at each tau. Set in conjunction with
auto_Re - For steady case, should
Re1be stepped automatically?
dec_Re - Can the critical value for
Re1only be found by a quasi-static decrease of
hyst_Re - Is the critical value of
Re1in a hysteresis region? (Specifically for 1- and 2-cell flows at very small aspect ratio).
Rarely used parameters
eps1,2 - Amplitude of oscillation of
Re_1,2(t,z)in axial direction.
freq1,2 - Frequency of oscillation of
Re_1,2(t,z)in axial direction.
- x_par_pos - Initial radial position of a particle in the fluid as a fraction of gap width.
z_par_pos - Initial axial position of a particle in the fluid as a fraction of
- save_part - Should a particle path be saved?
Use the makefile to compile the code and handle module dependencies.
- OBJECT - The compiled executable name.
- OBJS - The object files that should be linked.
- FC - The Fortran compiler.
- FFLAGS - Compiler flags.
- LDFLAGS - Any extra flags required by the linker.
parameters.f90 then run
source setup <directory>
which compiles code and copies
parameters.f90 and moves
./couette_mod runs code.
If restarting from a previous run then be sure to set
restart = .true. in
parameters.f90, and have the file
end_state.dat from the preceding run in
the run/submit directory.
Errors are output if either:
restart = .false.or
restart = .true.but
end_state.datdoes not exist.
If, at any time during a run, the file
SAVE is present in the run directory
(e.g. by giving the command
touch SAVE), the program will note this and
output data files for contours and surfaces at that time. Once the data files
SAVE is removed. This can be done at any time and as many times
as desired during the run.
Any run should be cleanly stopped by removing the empty file
RUNNING which is
present in the run directory once the job has started. This avoids having to
Ctrl-C to interrupt the program. Using this procedure will ensure that
all run-time files are cleaned up and unneeded directories removed.
- end_state.dat - Saves time index, p, time-step, dt, and fields u, Z, psi, current and magnetic field for restart.
- energy.dat - Saves the kinetic energy due to CCF and the total kinetic energy including CCF at each time.
- max_psi.dat - Saves maximum value of psi (and vr, vz) over whole field. Mainly for use with IDL plots.
- particle.dat - Saves the radial and axial position of a particle at each time.
time_tau.dat - If
auto_tau = .true.saves time at which each step in tau occured.
- torque.dat - Saves the torques on the inner (G1) and outer cylinders (G2), the torque due to CCF (Gc), and the ratio G1/Gc. For steady states G1=G2.
u_growth.dat - Saves time, radial velocity (outflow, inflow), growth rate, axial velocity, stream function, azimuthal velocity, vorticity, azimuthal current and magnetic field, and Reynolds number. Velocities are saved at the points defined in
p1234567.dat - Stream function field saved at the time defined by
(time-step)*(number following 'p')in filename. (Only if
xsect_save == .true.).
- u1234567.dat - As above but for azimuthal velocity field.
- z1234567.dat - As above but for vorticity field.
- vr1234567.dat - As above but for radial velocity field.
- vz1234567.dat - As above but for axial velocity field.
- j1234567.dat - As above but for azimuthal current field.
- b1234567.dat - As above but for azimuthal magnetic field.
xsect1234567.dat - Cross-sections of all fields except vorticity saved at the time defined as above. (Only if
xsect_save == .true.).