Clone wiki

gyre / Running GYRE


Introduction

These instructions assume you've already worked through the downloading, compiling and testing steps documented on the Getting Started page. If that's the case, then read on to learn how to run GYRE!

Making a Place to Work

When starting a new project, it's a good idea to create a dedicated directory to contain the various input and output files GYRE operates on. These commands will make a new subdirectory with the name $GYRE_DIR/work:

cd $GYRE_DIR
mkdir work

Grabbing a Stellar Model

The next step is to grab a stellar model which will form the basis for GYRE's calculation. There are a number of models provided beneath the $GYRE_DIR/models subdirectory; the following commands will copy a model for a 5M slowly pulsating B (SPB) star into the newly-created work area:

cp $GYRE_DIR/models/mesa/spb/spb.mesa $GYRE_DIR/work/

Creating a Namelist File

Now comes the fun part: creating an input file containing the various parameters which control a GYRE run. Using a text editor, create the file $GYRE_DIR/work/gyre.in with the following contents:

&constants
/

&model
	model_type = 'EVOL'  ! Obtain stellar structure from an evolutionary model
	file = 'spb.mesa'    ! File name of the evolutionary model
	file_format = 'MESA' ! File format of the evolutionary model
/

&mode
	l = 2                ! Harmonic degree
/

&osc
        outer_bound = 'VACUUM' ! Use a zero-pressure outer mechanical boundary condition
/

&num
	diff_scheme = 'COLLOC_GL4' ! 4th-order collocation scheme for difference equations
/

&scan
        grid_type = 'INVERSE' ! Scan for modes using a uniform-in-period grid; best for g modes
        freq_min = 0.5        ! Minimum frequency to scan from
	freq_max = 1.0        ! Maximum frequency to scan to
	n_freq = 250          ! Number of frequency points in scan
/

&grid
	alpha_osc = 10  ! Ensure at least 10 points per wavelength in propagation regions
	alpha_exp = 2   ! Ensure at least 2 points per scale length in evanescent regions
	n_inner = 5     ! Ensure at least 5 points between center and inner turning point
/


&ad_output
        summary_file = 'summary.txt'                            ! File name for summary file
	summary_file_format = 'TXT'                             ! Format of summary file
        summary_item_list = 'M_star,R_star,l,n_pg,omega,E_norm' ! Items to appear in summary file
        mode_template = 'mode.%J.txt'                		! File-name template for mode files
	mode_file_format = 'TXT'                   		! Format of mode files
        mode_item_list = 'l,n_pg,omega,x,xi_r,xi_h'   		! Items to appear in mode files
/

&nad_output
/

(The file doesn't have to be called gyre.in, but that's the convention followed in this documentation.)

This file is an example of a Fortran ''namelist'' file, containing multiple namelist groups. Each group begins with the line &name (where name is the name of the group); a list of name-value pairs follows, and the group ends with a slash /. Detailed information on the namelist groups expected in GYRE's input files can be found here; for now, let's just focus on some of the more-important aspects of the file above:

  • The constants namelist group is used to override constants such as the gravitational constant; here it's empty, indicating that default values should be used
  • The model namelist group instructs GYRE to read an evolutionary model, in MESA format, from the file spb.mesa
  • The mode namelist group instructs GYRE to consider quadrupole (ℓ=2) modes
  • The osc namelist group instructs GYRE to apply a zero-pressure outer mechanical boundary condition in the oscillation equations
  • The scan namelist group instructs GYRE to scan a region of dimensionless angular frequency space typically occupied by gravity modes
  • The grid namelist group instructs GYRE to perform calculations on a refinement of the model grid (see Understanding Grids for details on how this works)
  • The ad_output namelist group instructs GYRE to write out summary data to the file summary.txt, and individual mode data to files having the prefix mode-
  • The nad_output namelist group is empty, telling GYRE not to write out any non-adiabatic data (see Non-adiabatic Calculations for more info)

Running the Code

With the hard work done, it's now trivial to run the code:

cd $GYRE_DIR/work
../bin/gyre gyre.in

As it runs (on multiple cores, if you have a multi-core machine; see here for more details), the code will print lots of data to the screen. Let's break down this output, chunk by chunk. First, GYRE prints out its version number, tells us (in OpenMP threads) how many cores it is running on, and indicates which file it is reading parameters from (here, gyre.in):

gyre [5.0]
----------

OpenMP Threads   : 12
Input filename   : gyre.in

Next, GYRE loads the stellar model from the file spb.mesa. This model comprises 1814 points, and extends from the surface all the way to the center (which is why GYRE decides not to add a central point).

Model Init
----------

Reading from MESA file
   File name spb.mesa
   File version 1.00
   Read 1814 points
   No need to add central point

Next, GYRE prepares to searching for modes with harmonic degree l=2 and azimuthal order m=0 (not specified in gyre.in, but assumed by default), by setting up a frequency scan and a spatial (x) grid. The frequency scan specifies the frequencies at which GYRE samples the discriminant function D(omega), whose roots are the stellar eigenfrequencies (see Townsend & Teitler 2013); while the x grid specifies the spatial distribution of points used to discretize the oscillation equations.

Mode Search
-----------

Mode parameters
   l : 2
   m : 0

Building frequency scan
   added scan interval :  0.5000E+00 ->  0.1000E+01 (250 points, INVERSE)

Building x grid
   Found inner turning points, x range 0.1041 -> 0.1048
   Adding 0 inner point(s)
   Adding 21 global point(s) in iteration 1
   Adding 0 global point(s) in iteration 2
   Final grid has 1 segment(s) and 1835 point(s):
      Segment 1 : x range 0.0000 -> 1.0000 (1 -> 1835)

Next, GYRE attempts to bracket roots of the discriminant function by looking for changes in its sign:

Starting search (adiabatic)

Root bracketing
  Time elapsed :     3.314 s

Finally, for each sign change found GYRE uses a root solver to converge to the eigenfrequency. Each row of output here corresponds to an eigenfrequency:

Root Solving
   l    m    n_pg    n_p    n_g       Re(omega)       Im(omega)        chi n_iter
   2    0     -16      0     16  0.51863442E+00  0.00000000E+00 0.2679E-13      6
   2    0     -15      0     15  0.55636128E+00  0.00000000E+00 0.8611E-13      5
   2    0     -14      0     14  0.59457157E+00  0.00000000E+00 0.6864E-14      7
   2    0     -13      0     13  0.62301181E+00  0.00000000E+00 0.1796E-13      8
   2    0     -12      0     12  0.67563541E+00  0.00000000E+00 0.8523E-13      6
   2    0     -11      0     11  0.74334524E+00  0.00000000E+00 0.4626E-13      5
   2    0     -10      0     10  0.79690725E+00  0.00000000E+00 0.3680E-12      6
   2    0      -9      0      9  0.87153108E+00  0.00000000E+00 0.8141E-13      7
   2    0      -8      0      8  0.99747127E+00  0.00000000E+00 0.2521E-13      7
  Time elapsed :      0.580 s

The columns appearng are as follows:

  • l : harmonic degree
  • m : azimuthal order
  • n_pg : radial order (in the Eckart-Osaki-Scuflaire-Takata scheme)
  • n_p : acoustic-wave winding number
  • n_g : gravity-wave winding number
  • Re(omega) : dimensionless eigenfrequency (real part)
  • Im(omega) : dimensionless eigenfrequency (imaginary part; zero here because we've performed adiabatic calculations)
  • chi : convergence parameter
  • n_iter : number of iterations required for convergence

These values are printed to screen primarily to give an idea of GYRE's progress; more-detailed information about the modes found is given in the output files discussed below. Some things to watch out for:

  • The convergence parameter chi, defined as the ratio of discriminant values before and after the root finding, should small (on the order of 1E-9 to 1E-13). If it is significantly larger than this, the mode may not be properly converged; and if it is significantly smaller than this, there may be numerical issues with the discretization scheme.
  • The number of iterations n_iter should be moderate; values above 20 or so indicate that GYRE is having problems converging.
  • The mode index n_pg should be monotonic-increasing. Departures from this behavior can happen for a number of reasons:
    • Missing values can indicate that GYRE has skipped a mode in frequency space; the fix is to use a finer omega grid.
    • Missing values together with duplicate and/or non-monotonic values can indicate that GYRE isn't resolving the spatial structure of eigenfunctions; the fix is to use a finer x grid.
    • Missing values together with duplicate and/or non-monotonic values can also incdicate problems with the input stellar model --- for instance, incorrect values for the Brunt-Vaisala frequency across density discontinuities; the fix is to stop expecting GYRE to give sensible output when fed crap stellar models.

Interpreting Output Files

Overall properties of all modes found (eigenfrequencies, inertias, etc.) are collected together in the file summary.txt. For each mode GYRE also writes a file with the name mode.NNNNN.txt, containing data (eigenfrequency, eigenfunctions, etc.) specific to the mode. Here, NNNNN denotes a 5-digit index which increments (starting at 00001) for each mode found. Note that this index bears no relation to the radial order n_pg; it merely serves as a unique label for the modes.

Both the sumamry file and the mode files are text-based (it's possible to write HDF5-format files instead; see the Output Files page for details). Running

head summary.txt

will print out the first 10 lines of the summary file, which should look something like this:

                                                                                                                                                                                                                                                                
                        1                        2
                   M_star                   R_star
  0.9945999999999999E+034  0.3016908790335515E+012
                        1                        2                        3                        4                        5
                        l                     n_pg                Re(omega)                Im(omega)                   E_norm
                        2                      -16  0.5186344189658061E+000  0.0000000000000000E+000  0.1083833699332999E-002
                        2                      -15  0.5563612831705177E+000  0.0000000000000000E+000  0.1378396850031880E-002
                        2                      -14  0.5945715662438735E+000  0.0000000000000000E+000  0.3226917642759320E-002
                        2                      -13  0.6230118072964337E+000  0.0000000000000000E+000  0.3598212959967747E-002

The first three lines give column numbers, labels, and values for the scalar data — here, the stellar mass M_star and radius R_star, expressed in cgs units. The next two lines give column numbers and labels for the per-mode data (E_norm is the normalized mode inertia, and the other columns are the same as described above for the screen output); the subsequent lines then give the corresponding values (one line per mode). The mode files have a similar layout, with scalar data followed by array data representing the eigenfunctions (one line per radial grid point).

The choice of which data appear in output files isn't hardwired, but rather determined by the summary_item_list and mode_item_list parameters of the &ad_output and &nad_output namelist groups. Changing these parameters allows you to tailor the files to contain exactly the data you need. For a full list of possible items, consult the Output Files page.

Updated