CRUNCH3D

A 3D Fourier Collocation Code
for Dissipative, Compressible Magnetohydrodynamics

GENERAL DESCRIPTION

• Numerical 3D MHD
• Finite Viscosity and Resistivity
• Periodic Cartesian Geometry
• Parallelized FORTRAN Using MPI

NUMERICAL TECHNIQUES

• Fourier Collocation Discretization
• Second Order Runge-Kutta Time Integration

TYPICAL APPLICATIONS

• MHD Turbulence

• Magnetic Fluxtube Reconnection

• Where to Find Source Code

• How to Compile & Run CRUNCH3D

• CRUNCH3D Users Guide

• Subroutines in the CRUNCH Code

• Platforms Supported & Code Performance

• References

• Credits

• Point of Contact

WHERE TO FIND SOURCE CODE

The code can be found at the anonymous ftp site at lcp.nrl.navy.mil (enter anonymous as your username). Once you are in, cd to /pub/hpcc-ess. Everything is there in the compressed tar file crunch3d.t3d.tar.Z.

HOW TO COMPILE & RUN CRUNCH3D

Here we present a simplified description of how to run the CRUNCH code on the Cray T3D. The user first makes any desired chances to crunch.f, e.g., changing the initial conditions, changing the target directory for the output files (write large files to /scr/username ). The desired mesh size and number of processors must then be specified in param.inc. The code then is compiled and loaded using the Makefile, i.e., type "make". An executable module called crunch will be created. The user specifies desired characteristics for the run, e.g., number of processors and run time in the files zsub. To submit a job to the T3D type "qsub zsub" .

CRUNCH3D USERS GUIDE

FILES

crunch.f

This is the basic code for time advancing the compressible MHD equations. A full description is found in the SUBROUTINES IN THE CRUNCH CODE subsection.

param.inc

Use this file to set information about size characteristics of the simulation. The size of the grid can be changed by specifying new values for nx, ny, and nz, the sizes in the x, y, and z directions, respectively. The number of processors is set by nprocs: note that the same value for the number of processors must also be specified in the file zsub.

RUNDAT

Use this file to set information about characteristics of the simulation. The parameters that are set in RUNDAT are:

first row:

istart = 1: initiate a new simulation

istart = 2: restart an old simulation

nd = number of time steps for the simulation

kon = frequency of calls to subroutine output

second row:

rmu = value of viscosity

eta = value of resistivity

pr = magnetic prandtl number

third row:

rho0 = base value of mass density

e0 = base value of total energy density

fourth row:

dtmax = desired timestep for the numerical simulation

Makefile

Use this file to compile and load the code and create an executable module.

zsub

Use this file to submit a run to the T3D. The user should set the number of processors requested to be the same as the number specified in param.inc (at present this number is set to 2). The user should also change the cd command.

mpif.h

This is an MPI Fortran interface header file.

SUBROUTINES IN THE CRUNCH3D CODE

C3D

This subroutine contains both the initialization and the master time loop.

OUTPUTKCOL

This subroutine prints out a subset of the data at a frequency prescribed by the "kon" variable in the file INPUT. A sample of the fields being computed at a particular value of x and y is printed out as a function of z.

EIKPI

This subroutine initializes the fast Fourier transforms by executing calls to the Cray Mathlib subroutine scfft. Scale factors for the fft's are also precomputed here.

PDX

This subroutine computes the x-derivative of the function f(x,y,z), using a Fourier collocation derivative (for a discussion of the Fourier collocation derivative, see Canuto et al.(1984)). In sum, first the data is Fourier transformed on the first index (the x-direction), the result is multiplied by the appropriate wavenumbers, and the result is then inverse transformed. For forward and inverse fast Fourier transforms the subroutine uses the Cray Mathlib subroutines scfft and csfft. These subroutines are described on Cray man pages.

PDY

This subroutine computes the y-derivative of the function f(x,y,z), using a Fourier collocation derivative and a transposition in processor of the data. In sum, first the data is transposed to exchange the second array index with the first, then the result is Fourier on the first index of the array f, the result is multiplied by the appropriate wavenumbers, the result is then inverse transformed, and a final transposition is performed to return the first and second array indices to their original locations. For forward and inverse fast Fourier transforms the subroutine uses the Cray Mathlib subroutines scfft and csfft. These subroutines are described in the Cray man pages.

BFLD

This subroutine computes the magnetic field by computing the curl of an input magnetic vector potential.

STRESS

This subroutine computes the viscous stress tensor for use in the momenta equations and the equation for the total energy density. When called the updated value for the velocities should be in the u, v, and w arrays.

GRAD

This subroutine computes the gradient of an input scalar field. This routine is used only for computing temperature gradients needed for time-advancing the total energy density.

CD

This subroutine computes the negative of the electric current density from an input magnetic field. When it is called the updated value for the magnetic field should be in the bx, by and bz arrays. The negative of the curl of the magnetic field is returned.

BIGC

This subroutine forms the right-hand side of the continuity equations for the mass density, momenta, and total energy.

PDZ

This subroutine computes the z-derivative of the function f(x,y,z), using a Fourier collocation derivative and a transposition across processors of the data. In sum, first the data is transposed to exchange the third array index with the first, then the result is Fourier on the first index of the array f, the result is multiplied by the appropriate wavenumbers, the result is then inverse transformed, and a final transposition is performed to return the first and third array indices to their original locations. As mentioned previously, the transposes are performed across processors using MPI. For forward and inverse fast Fourier transforms the subroutine uses the Cray Mathlib subroutines scfft and csfft. These subroutines are described in the Cray man pages.

INPUT5

This subroutine reads in the contents of the file RUNDAT. The file RUNDAT contains information which determines computational and physical characteristics of the simulation.

WRITEDIOBINARY

This subroutine writes out files at the termination of the simulation. These files can be used to restart the run or for plotting or analysis. At present the code is set up to write to a file named AF. The user should create a directory called /scr/username and then change the call to WRITEDIOBINARY to reflect this change. This will be necessary for large files.

READDIOBINARY

This file reads in files to restart a simulation. At present the code is set up to read from a file named AN. The user should create a directory called /scr/username and then change the call to READDIOBINARY to reflect this change. This will be necessary for large files.

PLATFORMS SUPPORTED & CODE PERFORMANCE

This distribution of the code runs on the Cray T3D. Floprates can be found in the CRUNCH3D code performance report.

REFERENCES

REFERENCES - CRUNCH ALGORITHM

"Evolution of the Orszag-Tang Vortex System in a Compressible Medium. I. Initial Average Subsonic Flow" R. B. Dahlburg and J. M. Picone, Phys. Fluids B 1, 2153 (1989).

"Evolution of the Orszag-Tang Vortex System in a Compressible Medium. II. Supersonic Flow" J. M. Picone and R. B. Dahlburg, Phys. Fluids B 3, 29 (1981).

"The Kelvin-Helmholtz Instability in Photospheric Flows: Effects on Coronal Heating and Structure", J. T. Karpen, S. K. Antiochos, R. B. Dahlburg and D. S. Spicer, Astrophys. J. 403, 769 (1993).

"The Effects of Kelvin-Helmholtz Instability on Resonance Absorption Layers in Coronal Loops", J. T. Karpen, R. B. Dahlburg and J. M. Davila, Astrophys. J. 421, 372 (1994).

"Reconnection of Antiparallel Magnetic Fluxtubes," R. B. Dahlburg and S. K. Antiochos, J. Geophys. Res. 100, 23489 (1995).

"Parallel Computation of Magnetic Fluxtube Reconnection," R. B. Dahlburg and D. Norton, in "Small-Scale Structures in Three-Dimensional Hydro and Magnetohydrodynamic Turbulence," eds. M. Meneguzzi, A. Pouquet and P. -L. Sulem, Springer-Verlag, Heidelberg (1995), pp 317-325.

REFERENCES - SPECTRAL METHODS

"Numerical Analysis of Spectral Methods: Theory and Applications," D. Gottlieb and S. A. Orszag, SIAM, Philadelphia (1977).

"Spectral Methods for Partial Differential Equations," eds. R. Voight, D. Gottlieb and M. Y. Hussaini, SIAM, Philadelphia (1984).

"Spectral Methods in Fluid Mechanics,"C. Canuto, M. Y. Hussaini, A. Quarteroni and T. A. Zang, Springer-Verlag, New York (1988).

"Chebyshev and Fourier Spectral Methods," J. P. Boyd, Springer-Verlag, New York (1989).

CREDITS

The basic algorithm was put together by Russ Dahlburg. Contributors to the development of the algorithm over the years include Mike Picone, Jill Dahlburg, Judy Karpen and Dave Norton. This work is documented in the references listed at the end of this write-up. The lion's share of the MPI interface in the present code is the work of Edith Huang at JPL. Tom Clune of Cray Research provided valuable suggestion for increasing the floprate of the code.

This work was supported by the NASA HPCC/ESS program under Cooperative Agreement Notice 21425/041.

POINT OF CONTACT