python

Radiative transfer and ionisation code

Python is a Monte-Carlo radiative transfer code designed to simulate the spectrum of biconical (or spherical) winds in disk systems. It was origianally written by Long and Knigge (2002) and was intended for simulating the spectra of winds in cataclysmic variables. Since then, it has also been used to simulate the spectra of systems ranging from young stellar objects to AGN.

The name Python is today unfortunate, and changing the name is an ongoing debate within the development team. The program is written in C and can be compiled on systems runining various flavors of linux, including OSX on Macs.

The code is is available on github

Documentation

Various documentation exists:

  • A Quick Guide describing how to install and run Python (in a fairly mechanistic fashion).

For more information on how this page was generated and how to create documentation for python, look at the page for documentation on the documentation.

Authors

The authors of the python code and their institutions are:

Knox Long
Space Telescope Science Institute, 3700 San Martin Drive, Baltimore, MD 21218, USA Eureka Scientific, Inc., 2452 Delmer St., Suite 100, Oakland, CA 94602-3017, USA
Christian Knigge
Department of Physics and Astronomy, University of Southampton, Southampton, SO17 1BJ, UK
Stuart Sim
School of Mathematics and Physics, Queen’s University Belfast, University Road, Belfast, BT7 1NN, UK
Nick Higginbottom
Department of Physics and Astronomy, University of Southampton, Southampton, SO17 1BJ, UK
James Matthews
University of Oxford, Astrophysics, Keble Road, Oxford, OX1 3RH, UK
Sam Mangham
Department of Physics and Astronomy, University of Southampton, Southampton, SO17 1BJ, UK
Edward Parkinson
Department of Physics and Astronomy, University of Southampton, Southampton, SO17 1BJ, UK
Mandy Hewitt
School of Mathematics and Physics, Queen’s University Belfast, University Road, Belfast, BT7 1NN, UK

Getting Started

What machines will python run on? We have run python various versions of linux and on Mac. It is compiled using mpicc, with an option to compile with gcc. It uses the Gnu Scientific Libraries (gsl)

Installation

Python and the various routines associated are set up in a self-contained directory structure. The basic directory structure and the data files that one needs to run Python need to be retrieved and compiled.

If you want to obtain a stable (!) release, go to the Releases page.

If you want to download the latest dev version, you can zip up the git repository by clicking on the zip icon to the right of the GitHub page. Alternatively, clone it directly as

$ git clone https://github.com/agnwinds/python.git

If you anticipate contributing to development we suggest Forking the repository and submitting pull requests with any proposed changes.

Once you have the files, you need to cd to the new directory and set your environment variables

$ export PYTHON = /path/to/python/
$ cd $PYTHON
$ ./configure
$ make install
$ make clean

One can run a more rigorous clean of GSL with make distclean, or remove the compiled GSL libraries altogether with make rm_lib.

note that export syntax is for bash- for csh use

$ setenv PYTHON /path/to/python/

Atomic data is stored in our data repository with it’s own releases page- one should unzip these files and place them in a $PYTHON/data folder.

A development user may want to work on atomic data as part of their work, and pull in changes as they are made, in which case we recommend cloning the data repository:

$ cd $PYTHON; git clone https://github.com/agnwinds/data data

Running python

To run python you need to add the following to your $PATH variable:

$PYTHON/bin

You can then setup your symbolic links by running

$ Setup_Py_Dir

and run the code by typing, e.g.

$ py root.pf
Directory structure

The python directory structure is fairly simple:

source
Location of source code
bin
Location of executables
docs
Location of documentation, including sphinx docs, doxygen, parameters and documentation for the python programs in py_progs.
data
Location for all datafiles. Files that are mainly for reference should be gzipped to save space. Such files are not recreated in
bin
The location of the executables. (It is a good idea to put this directory in your path)
software
This directory contains libraries which are used in in python that must be recompiled when creating an installation on a new machine, primarily Bill Pence’s cfitsio package and the GNU scientific library gsl
py_progs
python programs for helping analyse the code. We recommend adding this directory to your PATH and PYTHON_PATH environment variables.
examples
A directory with a few examples of python runs. (Note that the input files will have changed and so one may not be able to run these examples without some changes in the input files.)
Please help by reporting bugs in installation

This can be done by submitting a bug under the Issues page

Running Python

The normal way to run Python is simply to enter

py xxx

where xxx is the root name of a parameter file. (The full name xxx.pf can also be entered).

However Python features a number of command line options which can be used to modify it’s operation. These include the following:

-h Causes Python to print out a brief help message and quit. The help message principally describes the command line options
-i (or –dry-run)
Causes Python to read and verify the inputs, writing a clean version of the input file xxx.pf to the output file xxx.out.pf, and then stop. This option is useful for setting up a proper .pf file. (Often one will want to copy xxx.out.pf back to xxx.pf before proceeding.
-t time_max Limits a run of python to approximately time_max in sec. This switch is used in situations where one would like to check whether the routine is operating properly be continuing, or where one needs to checkpoint the program after a certain period of time (due for example to time limits placed on jobs in a Beowulf cluster). The time is checked at the end of ionization and spectral cycles, immediately after saving the binary files that describe a model, and so one needs to leave a cushion between time_max and the maximum time one wants the program to run
-r Restarts a run that has been interrupted or halted, by reading a the xxx.windsave and xxx.specsave file (if it exists).
-v n Changes the amount of information printed to the screen by Python during a run. The default is 4. Larger numbers increase this. Smaller numbers decrease it. The log files are not affected.
--rseed Causes Python to use a random number seed that is time-based, rather than fixed. In most cases, a fixed seed is preferred so that problems can be replicated, but if is repeating the same calculation multiple times, then one may want a random seed.
--version Causes Python to print out the version number and commit hash (and whether uncommitted files exist, and then stop.
-p n_steps

Changes the number of photons generated during ionization cycles so that the number increases logarithmically to the maximum value. The number n_steps is optional, and specifies the number of decades over which the increase takes place.

Todo

NEED TO VERIFY THIS

Special switches

Python has a number of other switches that are not intended for the general user, but which may be useful in certain special cases. These include:

-d Enables a variety of specialized diagnostic inputs which have been implemented to help with solving various problems, and were regarded (by someone) as useful enough to maintain in the program. The user is then queried regarding which of these diagnostics to enable for a specific run. These diagnostic queries all start with @ (and can co-exist in the .pf file, with normal commands.
-e n Where n is a number, changes the number of errors of a specific type that are allowed to occur before the program gives up. For a variety of reasons, errors are expected during Python runs. Most of these errors are harmless in the sense that they occur rarely. But if an error occurs too often, something is seriously and so Python halts at that point. The default is \(10^{5}\) (per thread).
-e write n
Changes the number of times an error message of a specific type is written to a diagnostic file. When errors occur, a line describing the error is written to the diagnostic file the first n times the error occurs. After that statistics are maintained as to the number of times the error occurred, but it is not printed to the diagnostic file. The default is 100 (per thread)

Inputs

Todo

Fill in

Overview

Python uses a keyword based parameter file the specify a model. A portion of a parameter file (which must have the extension .pf) is as follows:

Wind.radiation(yes,no)                          yes
Wind.number_of_components                  1
Wind.type(SV,star,hydro,corona,kwd,homologous,yso,shell,imported)                   sv
Wind.coord_system(spherical,cylindrical,polar,cyl_var)          cylindrical
Wind.dim.in.x_or_r.direction               30
Wind.dim.in.z_or_theta.direction           30

Each line begins with a keyword followed optionally by a comment in parentheses, and then a value, e.g

  • Keyword: Wind.type
  • Comment: SV,star,hydro,corona,kwd,homologous,shell,imported
  • Value: SV

The comment generally specifies a set of valid choices or the units in which information is expected.

When a series of choices is presented, one does not need to enter the complete word, just enough to provide unique match to the choice.

One does not need to create a parameter file before running Python. Instead, assuming one is not working from a template parameter file, one simply invokes Python.

py my_new_model

or

py -n my_new_model

Python then queries the user for answers to a series of question, creating in the process a pf file, my_new_model.pf, that can be edited and used in future runs.

An example of a line presented to the user in interactive mode is:

Disk.mdot(msol/yr) (1e-08) :

There the number in the second set of parenthesis is a suggested value of the parameter. The user types in a new value and a carriage return, or, if the the suggested value seems appropriate, responds with a carriage return, in which case the suggested value will be used.

The -n switch above indicates that Python should accumulate all of the necessary inputs, write out the parameter file, and exit, which is useful if one is not completely sure what one wants.

System Description

The first set of parameters which Python needs are information about the overall system

System_type(star,cv,bh,agn,previous)                   bh

### Parameters for the Central Object
Central_object.mass(msol)                        10
Central_object.radius(cm)               8.85667e+06
Binary.mass_sec(msol)                           15
Binary.period(hr)                               72

### Parameters for the Disk (if there is one)
Disk.type(none,flat,vertically.extended)                 flat
Disk.radiation(yes,no)                          yes
Disk.rad_type_to_make_wind(bb,models)                   bb
Disk.temperature.profile(standard,readin)             standard
Disk.mdot(msol/yr)                            1e-6
Disk.radmax(cm)                                1e13

### Parameters for Boundary Layer or the compact object in an X-ray Binary or AGN
BH.radiation(yes,no)                            yes
BH.rad_type_to_make_wind(bb,models,power,cloudy,brems)                power
Boundary_layer.lum(ergs/s)              4.72063e+39
Boundary_layer.power_law_index                 -1.5

System_type is starting point, a basic classification of the type of object one is trying to model. This is used to guide further questions about the object and to set defaults.

Most of the other parameters are fairly self-explanatory, and are documented fully in the various Parameters entries.

Wind Models

Python allows for various types of models, which are defined by the following parameters

### Parameters describing the various winds or coronae in the system
Wind.radiation(yes,no)                          yes
Wind.number_of_components                  1
Wind.type(SV,star,hydro,corona,kwd,homologous,shell,imported)                   sv
Wind.coord_system(spherical,cylindrical,polar,cyl_var)          cylindrical
Wind.dim.in.x_or_r.direction               30
Wind.dim.in.z_or_theta.direction           30

Wind.radiation (WHICH PROBABLY WILL BE MOVED) allows for wind not only to scatter and absorb photons, but also to emit them by various processes, bound-bound, free-free, and recombination. It is the default for simple radiative transfer.

Wind.number_of_components is usually 1, but can be greater if one wishes to construct a wind from a combination of several wind models, for example a fast flow near the poles of a system, and a slow for near the disk. If the number of components exceeds 1, then the remaining questions relating to the wind will be posed multiple times.

The wind models incorporated into Python currently are:

SV
The Shlosman and Vitello parameterization of a bi-conical flow
Stellar_wind
A fairly standard parameterization of a spherical outflow for a hot star
Hydro
A special purpose mode used by us for importing models from Zeus and Pluto
Corona
A simple model for a corona above the disk
KWD
The Knigge Woods and Drew parameterization of a bi-conical flow
Homologous
A homologous expansion law useful for simulating SNe
Shell
A model of a thin shell useful for diagnostic studies
Imported
A general purpose mode for importing a wind from an ascii file

Todo

Update paths as they move

Importing Models

Python can read 1D or 2.5D grids of density and velocity, instead of setting up the model from an analytic prescription. Caution should be exercised with this mode, as it is still in a development phase, and the mode requires the user to ensure that things like mass and angular momentum conservation are enforced.

This mode is activated via wind type option “imported”, which triggers an extra question, e.g.

Wind.type(SV,star,hydro,corona,kwd,homologous,shell,imported)             imported
Wind.coord_system(spherical,cylindrical,polar,cyl_var)          cylindrical
Wind.model2import                    cv.import.txt

An example in cylindrical geometry, cv_import.pf, is given with a supplementary grid file in examples/beta/. The format expected in the grid input file for such a cylindrical model is as follows, although the column headers lines are actually not read.

i  j  inwind    x      z     v_x   v_y   v_z     rho   t_r
-- -- ------  -----  -----  ----- ----- -----   ----- -----
0  0    -1    1.4e9  3.5e9   0.0   0.0   6e5     0.0   0.0
0  1     0    1.4e9  3.5e10  1e5   0.0   2e6     1e9   0.0

where all physical units are CGS. i and j refer to the rows and columns of the wind cells respectively, while inwind tells the code whether the cell is in the wind (0), or out of the wind (-1). If a partially in wind flag is provided (1), the code defaults to treating this cell as not in the wind. This could in principle be adapted, but means that for the moment this mode is most useful when using models with sufficiently high resolution or covering factors that partially in wind cells are unimportant.

The other input files have slightly different formats. The best way to see the format is use the process described at the end of the page.

Creating your own model

In order to create your own model, there are a few important things to consider:

  • all units should be CGS (except for indices and flags, which are integers)
  • x and z for cylindrical (or r and theta for spherical polar) coordinates are supplied at the edges, rather than centres, of cells. Thus, a given cell is described by the location of it’s bottom left hand corner in (x,z) space.
  • Ghost cells must be included. This means that additional rows and columns of cells must be included at the edges of the grid, and they must be excluded from the wind so that their velocities and densities are set to zero, but have a velocity that python can interpolate with.
  • i and j correspond to rows and columns respectively, so that the first row of cells at the disk plane has i = 0.
  • rho the density of the cell in cgs units
  • The t_r column is currently not used, although it could be in future

Although cv_import.pf is designed to closely match the cv_standard.pf model, it does not match the model perfectly as the imported model does not deal with ‘partially in wind’ cells. As such, we generally recommend imported models are used for either wind models that entirely fill the grid or that have sufficiently high resolution that the partial filled cells are relatively unimportant.

Generating example inputs for testing and familiarizing oneself with Python’s import capability

If one is trying to use the import capability of Python for the first time, it will be useful to familiarize oneself with the process, and the file format for a particular coordinate system, by running first running Python on a model that is something similar to model to be imported, but which takes advantage of one of the kinematic models available with the code.

For example, suppose you have a hydrodynamical simulation of an AGN wind which is in polar coordinates and you want to use Python to calculate the spectrum. Then you might create a model of an AGN with a similar coordinate system using, say, a Knigge Wood & Drew wind (and similar atomic data). For specificity, suppose this model has the root name “test”

Once you have run the model, you can create an import file file by first running the routine windsave2table, or more specifically:

windsave2table test

This produces a large number of ascii tables, which are described elsewhere

In the py_progs directory, you will find 3 scripts, import_1d.py, import_cyl.py and import_rtheta.py, which will convert one of the output files test.0.master.txt to an import file, test.import.txt, that can be used with the import mode of Python. The 3 different routines are for 1d spherical coordinates, and polar (r-theta) coordinates respectively.

Assuming the py_progs directory is in your PATH, and given that our example is for cylindrical coordinates, one would run:

import_cyl.py test

At that point, you can test this import file, by modifying the first .pf file to import mode (imported). Running Python on this file, will result in your being asked the name of the import file, and give you a “baseline” to import the hydrodynamical simulation to work.

Note that one should not assume that spectra produced by the original run of Python and the run of the imported model will be identical. There are several reasons for this:

First, in creating the original model, Python accounts for the possibility that some cells are partially in the wind. This is not possible in the imported models. Only cells that are complete in the wind are counted.

Second, within Python, positions and velocities are assumed defined at the corners of cells, whereas densities are assumed to be cell centered. If one provides a table where all of the quantities are at the same exact position (namely density is at the same position as x), there will be a slight discrepancy between the way in model as calculated internally and as represented within Python.

Parameters

Todo

Fill in

System_type

The parameter is provides the program with a broad overview of the type of system that will be simulated, and is used by Python to initialize certain variable, and to control what variables are asked for later.

Type
Enumerator
Values
star
System in which the central object is a star
cv
System with a secondary star, which can occult the central object and disk depending on phase
bh
System with a black hole binary
agn
AGN
previous
In this case, one is starting from a previous run with python, and one want to either continue the run or change some parameters associated with radiation sources
File
python.c
Child(ren)
Central object

Todo

Fill in

Binary

Todo

Fill in

Binary.mass_sec

In binary systems the mass of the secondary. This is used along with the period to establish the Roche lobes, so that one can see the effects of eclipses on the system

Type
Double
Unit
M☉/year
Values
Greater than 0
File
setup_star_bh.c
Parent(s)
Binary.period

The perids of a binary system. Along with a mass, the binary period is used to define the Roche lobe of the system, which in turn can be used to see the effect of eclipses on the spectrum. Defining the system as a secondary also initializes the outer radius of the disk.

Type
Double
Unit
Hours
Values
Greater than 0
File
setup_star_bh.c
Parent(s)
Boundary_layer

Todo

Fill in

Boundary_layer.luminosity

The luminosity of the boundary layer.

Type
Double
Unit
ergs/s
Values
Greater than 0
File
setup_star_bh.c
Parent(s)
Boundary_layer.power_law_cutoff

This is a low frequency cutoff for an AGN-style power law spectrum of a form $L_nu=Knu^alpha$, as applied to the boundary layer of a star. It prevents the power-law being applied to low frequencies and giving an odd SED.

Type
Double
Unit
Hz
Values
Greater than 0
File
setup_star_bh.c
Parent(s)
Boundary_layer.power_law_index

The exponent 𝛼 in a power law SED applied to an AGN-style power law source for a non-AGN system. central source of the form $L_nu=Knu^alpha$.

Type
Double
Values
Any - but sign is not assumed, so for negative index use a negative value
File
setup_star_bh.c
Parent(s)
Boundary_layer.rad_type_in_final_spectrum

Determines the luminosity and SED of the boundary layer. The code can cause a source to radiate differently in the ionisation and spectral cycles. This variable allows a boundary layer source to radiate differently from Boundary_layer.rad_type_to_make_wind during the cycles used to calculate the output spectra. This can be

Type
Enumerator
Values
bb
Black-body radiation. The boundary layer radiates as a black-body source with surface luminosity set by its effective temperature (Boundary_layer.temp) and resulting in a total luminosity proportional to its surface area.
models
Radiate according to a model. Python can support tabulated models that output with a binned luminosity distribution depending on system properties like temperature and gravity. See Input_spectra.model_file. The total luminosity will be set by Boundary_layer.luminosity.
uniform
Available for System_type of star or cv. Photons are generated with a random, uniformly-distributed wavelength between Spectrum.wavemin and Spectrum.wavemax. Can in some cases substitute for a Kurcuz spectrum. This mode is only available when generating final spectra.
File
python.c
Parent(s)
Child(ren)
Boundary_layer.rad_type_to_make_wind

Determines the luminosity and SED of the boundary layer. The code can cause a source to radiate differently in the ionisation and spectral cycles. This variable allows a boundary layer source to radiate differently from Boundary_layer.rad_type_in_final_spectrum during the cycles used to calculate the wind ionisation state and temperature.

Type
Enumerator
Values
bb
Black-body radiation. The boundary layer radiates as a black-body source with surface luminosity set by its effective temperature (Boundary_layer.temp) and resulting in a total luminosity proportional to its surface area.
models
Radiate according to a model. Python can support tabulated models that output with a binned luminosity distribution depending on system properties like temperature and gravity. See Input_spectra.model_file. The total luminosity will be set by Boundary_layer.luminosity.
power
Radiate following a power-law model as $L_nu=Knu^alpha$. The total luminosity will be set by Boundary_layer.luminosity.
File
setup_star_bh.c
Parent(s)
Child(ren)
Boundary_layer.radiation

Says whether the boundary layer will radiate.

Type
Boolean (yes/no)
File
setup_star_bh.c
Parent(s)
Child(ren)
Boundary_layer.temp

The temperature of the boundary layer when radiating as a black body.

Type
Double
Unit
Kelvin
Values
Greater than 0
File
setup.c
Parent(s)
Central_object

Todo

Fill in

Central_object.blackbody_temp

If the AGN/BH is radiating as a black body, what temperature should it radiate at?

Type
Double
Unit
Kelvin
Values
Greater than 0
File
setup_star_bh.c
Parent(s)
Central_object.bremsstrahlung_alpha

The frequency exponent 𝛼 in bremstrahlung SED of the form $L_nu=nu^{alpha}e^{-hnu/kT}$

Type
Double
Values
Any - sign is not assumed so use negative if you want negative
File
setup_star_bh.c
Parent(s)
Central_object.bremsstrahlung_temp

The temperature T in bremstrahlung SED of the form $L_nu=nu^{alpha}e^{-hnu/kT}$

Type
Double
Unit
K
Values
Greater than 0
File
setup_star_bh.c
Parent(s)
Central_object.cloudy.high_energy_break

This is a command to define a cloudy type broken power law SED - mainly used for testing the code against cloudy. This SED has hardwired frequency exponents of 2.5 below the low energy break and -2.0 above the high energy break. This parameter defines the energy of the high energy break.

Type
Double
Unit
eV
Values
Greater than Central_object.cloudy.low_energy_break
File
setup_star_bh.c
Parent(s)
Central_object.cloudy.low_energy_break

This is a command to define a cloudy type broken power law SED - mainly used for testing the code against cloudy. This SED has hardwired frequency exponents of 2.5 below the low energy break and -2.0 above the high energy break. This parameter defines the energy of the low energy break.

Type
Double
Unit
eV
Values
Greater than 0
File
setup_star_bh.c
Parent(s)
Central_object.geometry_for_source
If the central source in an AGN/BH system is radiating, what geometry should it radiate from?
This is applicable even for black-body sources, where the luminosity depends on Central_object.radius.
Type
Enumerator
Values
lamp_post
The central source radiates from two point sources located on the system axis above and below the disk plane. Emission is completely isotropic.
sphere
The central source radiates from a spherical surface with radius Central_object.radius. Emission is cosine-weighted in the direction of the sphere normal at the point of emission. Photons that would be spawned in an extended disk (if Disk.type is vertically.extended) are re-generated.
File
setup_star_bh.c
Parent(s)
Child(ren)
Central_object.lamp_post_height

The distance above and below the disk plane of the two point sources used in the lamp-post model.

Type
Double
Unit
Central_object.radius
Values
Greater than 0
File
setup_star_bh.c
Parent(s)
Central_object.luminosity

The luminosity of a non-blackbody AGN central source. This is defined as the luminosity from 2-10keV.

Type
Double
Unit
ergs/s
Values
Greater than 0.
File
setup_star_bh.c
Parent(s)
Central_object.mass

Mass of the central object. This is very important, affecting wind speeds, gravitational heating and such.

Type
Double
Unit
M☉
Values
Greater than 0
File
setup_star_bh.c
Central_object.power_law_cutoff

Adds a low-frequency cutoff to the power law spectrum. Whilst this is required for power-law emission modes, it’s set globally and also used in cloudy broken-power-law emission modes!

Type
Double
Unit
Hz
Values
Greater than or equal to 0
File
setup_star_bh.c
Parent(s)
Central_object.power_law_index

The exponent 𝛼 in a power law SED applied to a power law source of the form $L_nu=Knu^alpha$.

Type
Double
Values
Greater than 0
File
setup_star_bh.c
Parent(s)
Central_object.rad_type_in_final_spectrum

Determines the SED of the central object in the spectral cycles. The luminosity is set by the options for the ionisation cycles, however.

Type
Enumerator
Values
bb
Available for System_type of star or cv. Black-body radiation. The boundary layer radiates as a black-body source with surface luminosity set by its effective temperature (Central_object.temp) and resulting in a total luminosity proportional to its surface area.
models
Available for System_type of star or cv. Radiate according to a model. Python can support tabulated models that output with a binned luminosity distribution depending on system properties like temperature and gravity. See Input_spectra.model_file. The total luminosity will be set by Central_object.luminosity.
uniform
Available for System_type of star or cv. Photons are generated with a random, uniformly-distributed wavelength between Spectrum.wavemin and Spectrum.wavemax. Can in some cases substitute for a Kurcuz spectrum. This mode is only available when generating final spectra.
brems
Available for System_type of agn or bh. Central object radiates with SED of a brehmsstralung spectrum as $L_nu=nu^{alpha}e^{-hnu/kT}$. This was originally developed to allow comparison to spectra generated according to Blondin heating and cooling rates.
cloudy
Available for System_type of agn or bh. Central object radiates with a ‘broken’ power law, intended largely for testing purposes against Cloudy. The SED form is $L_nu=Knu^alpha$, but beyond the provided high and low energy breakpoints the luminosity falls off sharply.
power
Available for System_type of agn or bh. Radiate following a power-law model as $L_nu=Knu^alpha$. The total luminosity will be set by Boundary_layer.luminosity.
File
python.c
Parent(s)
Child(ren)
Central_object.rad_type_to_make_wind

Multi-line description, must keep indentation.

Type
Enumerator
Values
bb
Black-body radiation. The boundary layer radiates as a black-body source with surface luminosity set by its effective temperature (Central_object.temp) and resulting in a total luminosity proportional to its surface area.
models
Radiate according to a model. Python can support tabulated models that output with a binned luminosity distribution depending on system properties like temperature and gravity. See Input_spectra.model_file. The total luminosity will be set by Central_object.luminosity.
brems
Available for System_type of agn or bh. Central object radiates with SED of a brehmsstralung spectrum as $L_nu=nu^{alpha}e^{-hnu/kT}$. This was originally developed to allow comparison to spectra generated according to Blondin heating and cooling rates.
cloudy
Available for System_type of agn or bh. Central object radiates with a ‘broken’ power law, intended largely for testing purposes against Cloudy. The SED form is $L_nu=Knu^alpha$, but beyond the provided high and low energy breakpoints the luminosity falls off sharply.
power
Available for System_type of agn or bh. Radiate following a power-law model as $L_nu=Knu^alpha$. The total luminosity will be set by Boundary_layer.luminosity.
File
setup_star_bh.c
Parent(s)
Child(ren)
Central_object.radiation

A booliean variable stating whether of not the central object should radiate. This will enable different follow-up questions depending on the system type.

Type
Boolean (yes/no)
File
setup_star_bh.c
Child(ren)
Central_object.radius

Radius of the central object in the system, e.g the white dwarf or black hole

Type
Double
Unit
cm
Values
Greater than 0
File
setup_star_bh.c
Central_object.temp

Temperature of the central star. Physically, this is used in blackbody radiation, shock heating and disk heating in YSO models. It is also used to help determine the frequency bands in which photons are emitted.

Type
Double
Unit
Kelvin
Values
Greater than zero
File
setup_star_bh.c
Parent(s)
Disk

Todo

Fill in

Disk.T_profile_file

When the user chooses to read in the temperature profile as a function of radius, the user is asked the name of the file that contains the desired profile.

Type
String
File
setup_disk.c
Parent(s)
Disk.mdot

The mass transfer rate in the disk when considering a standard Shakura-disk.

Type
Double
Unit
M☉/year
File
setup_disk.c
Parent(s)
Disk.rad_type_in_final_spectrum

Multi-line description, must keep indentation.

Type
Enumerator
Values
bb
Multi-line description, must keep indentation.
models
Multi-line description, must keep indentation.
uniform
Multi-line description, must keep indentation.
File
python.c
Parent(s)
Child(ren)
Disk.rad_type_to_make_wind

Multi-line description, must keep indentation.

Type
Enumerator
Values
bb
Multi-line description, must keep indentation.
models
Multi-line description, must keep indentation.
File
setup_disk.c
Parent(s)
Child(ren)
Disk.radiation

Multi-line description, must keep indentation.

Type
Boolean(yes/no)
File
setup_disk.c
Parent(s)
Child(ren)
Disk.radmax

The outer edge of the disk. Photons inside this radius are absorbed or re-radiated. Photons which are outside this radius pass through the disk plane.

Type
Double
Unit
cm
Values
Greater than 0
File
setup_disk.c
Parent(s)
Disk.temperature.profile

The choice of disk temperature profile

Type
Enumerator
Values
standard
A Shakura - Sunyaev disk, with a hard inner boundar
readin
Read the profile in from a file; the user will be queried for the name of the file
yso
YSO???
analytic
DEPRECATED??? A profile designed for the situation where the disk is being illuminated by star
File
setup_disk.c
Parent(s)
Child(ren)
Disk.type

Parameter defining whether there is a disk in the system

Type
Enumerator
Values
none
No disk
flat
Standard flat disk
vertically.extended
Vertically extended disk
File
setup_disk.c
Child(ren)
Disk.z0

Fractional height at maximum radius. The physical height at the outer disk will be this * Disk.radmax.

Type
Double
Values
Greater than 0
File
setup_disk.c
Parent(s)
Disk.z1

For a vertically extended the disk, the height of the disk is set to be Disk.z0 * Disk.radmax * (r/Disk.radmax)**Disk.z1 where Disk.z1 is the power law index

Type
Double
Values
Greater than 0
File
setup_disk.c
Parent(s)
Wind

Todo

Fill in

Corona

Todo

Fill in

Corona.base_den

The coronal model is defined in terms of a base density and a scale height

Type
Double
Unit
number/cm**3
Values
Greater than 0
File
corona.c
Parent(s)
Corona.radmax

The corona is a box-shaped region which sits immediately above the disk. radmax defines the outer edge of the box.

Type
Double
Unit
cm
Values
Greater than Central_object.radius
File
corona.c
Parent(s)
Corona.radmin

The corona is a box-shaped region which sits immediately above the disk. radmin defines the inner edge of the box.

Type
Double
Unit
cm
Values
Greater than Central_object.radius
File
corona.c
Parent(s)
Corona.scale_height

The coronal model is defined in terms of a base density and a scale height

Type
Double
Unit
cm
Values
Greater than 0
File
corona.c
Parent(s)
Corona.vel_frac

For the coronal model, the azimuthal velocity is given by the velocity of the underlying disk. One can also give the corona a radial velocity, which is a fraction of the disk velocity. (As coded, if this number is positive, the velicty is the r direction is toward the central object).

Type
Double
Unit
Disk velocity
Values
Any, 0 implies no radial velocity.
File
corona.c
Parent(s)
Corona.zmax

The corona is a box-shaped region which sits immediately above the disk. zmax defines the height of the box.

Type
Double
Unit
cm
Values
Greater than that the radius of the central object
File
corona.c
Parent(s)
Homologous

Todo

Fill in

Homologous.boundary_mdot

The mass loss rate at the base of the wind in a homlogous flow model, a flow in which the velocity is proporional to the radius. In general, mdot will decline with radius, depending on the exponent of the power law that describes the trend in density.

Type
Double
Unit
M☉/yr
Values
Greater than 0
File
homologous.c
Parent(s)
Homologous.density_exponent

The power law exponent which defines the decline in density of a homologous flow as a function of radious.

Type
Double
Values
Greater than 0 for a density that declines with radius
File
homologous.c
Parent(s)
Homologous.radmax

Maximum extent of the homologous wind.

Type
Double
Unit
cm
Values
Greater than Homologous.radmin
File
homologous.c
Parent(s)
Homologous.radmin

The starting point of for madel of a homologous flow, a model in which the velocity at any radius is proportional to the radius

Type
Double
Unit
cm
Values
Greater than or equal to Central_object.radius
File
homologous.c
Parent(s)
Homologous.vbase

Velocity at the base of the wind

Type
Double
Unit
cm
Values
Greater than 0
File
homologous.c
Parent(s)
Hydro

Todo

Fill in

Hydro.file

Relative path to a hydrodynamic snapshot file to be imported.

Type
String
File
hydro_import.c
Parent(s)
Hydro.thetamax

The maximum theta value to be read in from a hydrodynamic snapshot. This is typically used to excise a dense disk from the midplane of such a snapshot. Use a negative value to tell the code to use all the data.

Type
Double
Unit
Degrees
Values
-1 use all data
X
use up to that angle (typically less than 90)
File
hydro_import.c
Parent(s)
KWD

Todo

Fill in

KWD.acceleration_exponent

Sets the length scale over which the accleration to v_inf is accomplished. It is the value of the exponent beta for the Caster & Lamers equation of a stellar wind, v(r) = v_0 + (v_inf - v_0) * (1 - R_s/r) ** beta.

Type
Double
Values
Greater than 0
File
knigge.c
Parent(s)
KWD.acceleration_length

The size of the acceleration length scale for a disk wind described by the KWD model.

Type
Double
Unit
cm
Values
Greater than 0
File
knigge.c
Parent(s)
KWD.d

The ratio d/d_min is used to describe the degree of geometric collimation of the disk wind in the KWD model. However, d (the distance to the focal point in central object radii) is used as this provides a more natural parameter.

Type
Double
Unit
Central_object.radius
Values
Greater than 0
File
knigge.c
Parent(s)
KWD.mdot_r_exponent

The exponent for the mass loss rate as defined in the KWD model, m_dot(r) = F(r) ** alpha = T(r) ** (4 * alpha). F is the local luminous flux and T is the local temperature at a radius R. A value of 0 sets a uniform mass loss rate.

Type
Double
Values
Greater than or equal to 0
File
knigge.c
Parent(s)
KWD.rmax

The radius at which the disk wind terminates, in units of central object radii. This has to be greater than rmin.

Type
Double
Unit
Central_object.radius
Values
Greater than KWD.rmin
File
knigge.c
Parent(s)
KWD.rmin

The radius at which the disk wind begins, in units of central object radii.

Type
Double
Unit
Central_object.radius
Values
Greater than 1
File
knigge.c
Parent(s)
KWD.v_infinity

The velocity at large distances of a steller wind described by the KWD model, in units of escape velocity. Described in terms of Castor & Lamers equation, v(r) = v_0 + (v_inf - v_0) * (1 - R_s/r) ** beta.

Type
Double
Unit
Escape velocity
Values
Greater than 0
File
knigge.c
Parent(s)
KWD.v_zero

Multiple of the local sound speed at the base of the wind, this results in the initial velocity of the wind being able to be greater or less than the local sound speed.

Type
Double
Unit
Sound speed at wind base
Values
Greater than 0
File
knigge.c
Parent(s)
SV

Todo

Fill in

SV.acceleration_exponent

Power-law acceleration exponent (i.e. alpha) of a line driven wind in a Shlosman & Vitello (SV) CV disk wind model. Sets the length scale over which the accleration to v_inf is accomplished. This value is a constant; when equal to 1 the results resemble those of a linear velocity law. Typically for an SV type wind this power law exponent is 1.5. See equation (2) Shlosman & Vitello 1993, ApJ 409, 372.

Type
Double
Values
Greater than 0
File
sv.c
Parent(s)
SV.acceleration_length

The size of the acceleration length scale for a disk wind described by the Shlosman Vitelo model. See equation (2) Shlosman & Vitelo ApJ (1993),409,372

Type
Double
Unit
cm
Values
Greater than 0
File
sv.c
Parent(s)
SV.diskmax

The outermost radius from which the wind rises in a Shlossman-Vitello type disk wind. This radius is measured along the radial disk (r) direction i.e. zero describes the centre of the central object (white dwarf) See figure 1 of Shlosman & Vitello 1993, ApJ 409,372.

Type
Double
Unit
cm
Values
Greater than or equal to SV.diskmin (inner radius disk wind)
File
sv.c
Parent(s)
SV.diskmin

The innermost radius from which the wind rises in a Shlossman-Vitello type disk wind. This radius is measured along the radial disk (r) direction i.e. zero describes the centre of the central object (white dwarf) See figure 1 of Shlosman & Vitello 1993, ApJ 409,372.

Type
Double
Unit
cm
Values
Greater than or equal to Central_object.radius
File
sv.c
Parent(s)
SV.mdot_r_exponent

The exponent for the mass loss rate as defined in the Shlosman Vitelo model, See lambda in equation (4) Shlosman & Vitelo,ApJ,1993,409,372.

Type
Double
Values
Greater than or equal to 0. 0 sets a uniform mass loss rate.
File
sv.c
Parent(s)
SV.thetamax

The angle at which the wind rises from the outermost launching radius in a Shlossman-Vitello type disk wind. This angle is measured with respect to the vertical (z) direction i.e. zero describes a vertical wind. See figure 1 of Shlossman & Vitello 1993, ApJ 409,372.

Type
Double
Unit
Degrees
Values
Greater than sv.thetamin
File
sv.c
Parent(s)
SV.thetamin

The angle at which the wind rises from the innermost launching radius in a Shlossman-Vitello type disk wind. This angle is measured with respect to the vertical (z) direction. I.e. zero descirbes a vertical wind. See figure 1 of Shlossman & Vitello 1993, ApJ, 409, 372.

Type
Double
Unit
Degrees
Values
Greater than 0
File
sv.c
Parent(s)
SV.v_infinity

Asymptotic (i.e. final) velocity of a line driven wind in a Shlosman & Vitello CV disk wind model. Assumed to scale with the local velocity at the base of the streamline. See equation (2) Shlosman & Vitello 1993, ApJ 409, 372.

Type
Double
Unit
Escape velocity
Values
Greater than 0
File
sv.c
Parent(s)
SV.v_zero

The velocity at the wind base.

Type
Double
Unit
[‘Speed of sound in the wind’, ‘cm/s’]
Values
Greater than 0
File
sv.c
Parent(s)
SV.v_zero_mode

Multi-line description, must keep indentation.

Type
Enumerator
Values
fixed
Multi-line description, must keep indentation.
sound_speed
Multi-line description, must keep indentation.
File
sv.c
Parent(s)
Child(ren)
Shell

Todo

Fill in

Shell.wind.acceleration_exponent

Exponent beta for the Caster and Lamers description of a stellar wind v(r)=v_o + (v_inf - v_o) (1+R_s/r)**beta for a shell wind.

Type
Double
Values
Greater than or equal to 0
File
shell_wind.c
Parent(s)
Shell.wind.radmax

Multi-line description, must keep indentation.

Type
Double
Unit
cm
Values
Greater than Shell.wind.radmin
File
shell_wind.c
Parent(s)
Shell.wind.radmin

The innermost edge of a diagnostic type of wind made up of a single (ideally thin) shell.

Type
Double
Unit
cm
Values
Greater than 0
File
shell_wind.c
Parent(s)
Shell.wind.v_at_rmax

The velocity of a shell wind at the outer edge of the shell - the variation of the velocity in the shell is set by the velocity law exponent. It allows a gradient to be enforced.

Type
Double
Unit
cm/s
Values
Greater than or equal to 0
File
shell_wind.c
Parent(s)
Shell.wind_mdot

The mass loss through a diagnostic shell type wind. One normally sets this experimentally in order to get a required hydrogen density in the shell

Type
Double
Unit
M☉/year
Values
Greater than 0
File
shell_wind.c
Parent(s)
Shell.wind_v_at_rmin

The velocity of a shell wind at the inner edge of the shell - the variation of the velocity in the shell is set by the velocity law exponent. It allows a gradient to be enforced.

Type
Double
Unit
cm/s
Values
Greater than or equal to 0
File
shell_wind.c
Parent(s)
Stellar_wind

Todo

Fill in

Stellar_wind.acceleration_exponent

Exponent beta for the Caster and Lamers description of a stellar wind v(r)=v_o + (v_inf - v_o) (1+R_s/r)**beta

Type
Double
Values
Greater than or equal to 0
File
stellar_wind.c
Parent(s)
Stellar_wind.mdot

Mass loss rate for a wind modelled in terms of the Caster and Lamemers prescription for a stellar wind.

Type
Double
Unit
M☉/year
Values
Greater than 0
File
stellar_wind.c
Parent(s)
Stellar_wind.radmax

Multi-line description, must keep indentation.

Type
Double
Unit
cm
Values
Greater than or equal to Stellar_wind.radmin
File
stellar_wind.c
Parent(s)
Stellar_wind.radmin

Inner edge in cm for a stellar wind, normally the radius of the star.

Type
Double
Unit
cm
Values
Greater than or equal to Central_object.radius
File
stellar_wind.c
Parent(s)
Stellar_wind.v_infinity

The velocity at large distance of a stellar wind described in terms of the Casters and Larmers equation v(r)=v_o + (v_inf - v_o) (1+R_s/r)**beta

Type
Double
Unit
cm/s
Values
Greater than 0
File
stellar_wind.c
Parent(s)
Stellar_wind.vbase

Multi-line description, must keep indentation.

Type
Double
Unit
cm/s
Values
Condition e.g. greater than 0 or list e.g. [1, 2, 5]
File
stellar_wind.c
Parent(s)
Wind

Todo

Fill in

Wind.coord_system

The coordinate system used for a describing a component of the wind.

Type
Enumerator
Values
spherical
Spherical
cylindrical
Cylindrical
polar
Spherical polar
cyl_var
Cylindrical varying z
File
setup_domains.c
Parent(s)
Wind.dim.in.x_or_r.direction

Winds are calulated on spherical, cylindrical, or polar grids. This input variable gives the size of the grid in the x or r direction. Because some grid cells are used as a buffer, the actual wind cells are contained in a slightly smaller grid than the number given.

Note that in some situations there may be more than one wind component, known technically as a domain. In that case the user will be queried for this value mulitple times, one for each domain

Type
Integer
Values
Greater than or equal to 4, to allow for boundaries.
File
setup_domains.c
Parent(s)
Wind.dim.in.z_or_theta.direction

Winds are calulated on spherical, cylindrical, or polar grids. This input variable gives the size of the grid in the z or theta direction. Because some grid cells are used as a buffer, the actual wind cells are contained in a slightly smaller grid than the number given.

Note that in some situations there may be more than one wind component, known technically as a domain. In that case the user will be queried for this value mulitple times, one for each domain

Type
Integer
Values
Greater than 0
File
setup_domains.c
Parent(s)
Wind.filling_factor

The volume filling factor of the outflow. The implementation of clumping (microclumping) is described in Matthews et al. (2016), 2016MNRAS.458..293M. Asked once per domain.

Type
Double
Values
0 < f <= 1, where 1 is a fully smooth wind.
File
setup_domains.c
Parent(s)
Wind.fixed_concentrations_file

The filename for the fixed ion concentrations if you have set Wind_ionization to 2 (fixed). This file has format [atomic_number ionizationstage ion fraction].

Type
String
File
setup.c
Parent(s)
Wind.ionization

Multi-line description, must keep indentation.

Type
Enumerator
Values
LTE_te
Multi-line description, must keep indentation.
LTE_tr
Multi-line description, must keep indentation.
ML93
Multi-line description, must keep indentation.
fixed
Multi-line description, must keep indentation.
matrix_bb
Multi-line description, must keep indentation.
matrix_pow
Multi-line description, must keep indentation.
on.the.spot
Multi-line description, must keep indentation.
File
setup.c
Child(ren)
Wind.mdot

Multi-line description, must keep indentation.

Type
Double
Unit
M☉/year
Values
Greater than 0
File
[‘knigge.c’, ‘sv.c’]
Parent(s)
Wind.model2import

The name of a file to containing a generic model to read in to python from an ascii file. (Note that this is not the same as reading in a model generated by python, but is intended to allow one to read in a generic model in a variety of formats with only a limited amount of information required).

Type
String
File
import.c
Parent(s)
Wind.number_of_components

While most simple description of a wind consist of a single region of space, Python can calculate radiative transfer through more complicated structres, where one region of space is described with one prescription and another region of space with a second prescription. For example, one might want to place a disk atmoosphere between the disk and a wind. This parameter describes the number of components (aka domains) of the wind.

Type
Integer
Values
Greater than 0
File
python.c
Parent(s)
Child(ren)
Wind.old_windfile

The rootname of a previously saved model in a calculation one wishes to continue (with the possiblity of making changes to some of the details of the radiation sources, or to extract spectra from different inclinations)

Type
String
File
python.c
Parent(s)
Wind.radiation

Whether or not the wind should radiate.

Type
Boolean (yes/no)
File
python.c
Wind.radmax

Multi-line description, must keep indentation.

Type
Double
Unit
cm
Values
Greater than Central_object.radius and any minimum wind radii in the system.
File
setup_domains.c
Parent(s)
Wind.t.init

Starting temperature of the wind.

Type
Double
Unit
Kelvin
Values
Greater than 0
File
setup_domains.c
Parent(s)
Wind.type

Multi-line description, must keep indentation.

Type
Enumerator
Values
SV
Multi-line description, must keep indentation.
corona
Multi-line description, must keep indentation.
homologous
Multi-line description, must keep indentation.
hydro
Multi-line description, must keep indentation.
imported
Multi-line description, must keep indentation.
kwd
Multi-line description, must keep indentation.
shell
Multi-line description, must keep indentation.
star
Multi-line description, must keep indentation.
yso
Multi-line description, must keep indentation.
File
setup_domains.c
Parent(s)
Child(ren)
Radiative Transfer & Ionisation

Todo

Fill in

Atomic_data

Python uses an atomic data file, as found in the agnwinds/data repository. This is the relative path to the Atomic Data header file on disk. See Atomic Data

Type
String
File
setup_line_transfer.c
Parent(s)
Ionization_cycles

The number of ionization cycles to execute - these are cycles to determine the ionization and thermal state of the wind

Type
Integer
Values
Greater than 0
File
setup.c
Line_transfer

The way in which line transfer and scattering is dealt with in the code. Governs whether we adopt any approximations for radiative transfer, whether to use the indivisible packet and macro-atom machinery, and whether to use isotropic or anisotropic scattering.

Thermal trapping mode is recommended for non-macro atom runs, while thermal trapping macro-atom mode is recommended for macro-atom runs.

Type
Enumerator
Values
pure_abs

Pure absorption

The pure absortion approximation.

pure_scat

Pure scattering

The pure scattering approximation.

sing_scat

Single scattering

The single scattering approximation.

escape_prob

Escape probability

Resonance scattering and electron scattering is dealt with isotropically. free-free, compton and bound-free opacities attenuate the weight of the photon wind emission produces additional photons, which have their directions chosen isotropically. The amount of radiation produced is attenuated by the escape probability.

thermal_trapping

Escape probability + anisotropic scattering

We use the ‘thermal trapping method’ to choose an anistropic direction when an r-packet deactivation or scatter occurs.

macro_atoms

Macro-atoms

use macro-atom line transfer. Packets are indivisible and thus all opacities are dealt with by activate a macro-atom, scattering, or creating a k-packet. the new direction following electron scattering or deactivation of a macro atom is chosen isotropically.

macro_atoms_thermal_trapping

Macro-atoms + anisotropic scattering

as macro_atoms, but we use the ‘thermal trapping method’ to choose an anistropic direction when an r-packet deactivation or scatter occurs.

File
setup_line_transfer.c
Child(ren)
Photons_per_cycle

Multi-line description, must keep indentation.

Type
Double
Values
Greater than 0
File
setup.c
Surface.reflection.or.absorption

When photons hit the disk, there are several options

Type
Enumerator
Values
reflect
The photons are scattered back into the wind
absorb
The photons are simply lost from the calculation
thermalized.rerad
The photons are absorbed, in the next ionization cycle energy lost is treated as extra heat, and the effective temperature of the ring in the disk will be increased accordingly
File
setup.c
Wind_heating

Todo

Fill in

Wind_heating.extra_luminosity

This is a very special option put in place for modelling FU Ori stars, and should be used with extreme caution. Determines the shock factor.

Type
Double
Values
Condition e.g. greater than 0 or list e.g. [1, 2, 5]
File
setup.c
Parent(s)
Wind_heating.extra_processes

Multi-line description, must keep indentation.

Type
Enumerator
Values
adiabatic
Multi-line description, must keep indentation.
both
Multi-line description, must keep indentation.
none
Multi-line description, must keep indentation.
nonthermal
Multi-line description, must keep indentation.
File
setup.c
Child(ren)
Wind_heating.kpacket_frac

Multi-line description, must keep indentation.

Type
Double
Unit
None
Values
Condition e.g. greater than 0 or list e.g. [1, 2, 5]
File
setup.c
Parent(s)
Spectrum

Todo

Fill in

Spectrum.angle

The inclination angle with respect to the polar axis for obtaining a spectrum. This question will be repeated once for each desired incliniation

Type
Double
Unit
Degrees
Values
0 to 90 degrees, where 0 is normal to the disk and 90 is on the disk plane
File
setup.c
Parent(s)
Spectrum.live_or_die

Normally in creating detailed spectrum Python “extracts” photons in a certain direction reweighting them to account for the fact that they have been extracted in a certain direction. It is possible to just count the photons that are emitted in a single angle range. The two methods should yield the same or very similar results but the extraction method is much more efficient and live or die is basically a diagnostic mode.

Type
Enumerator
Values
live.or.die
Count only those photons that escape within a small angle range towards the observer
extract
Extract a component of all photons that scatter towards the observer
File
setup.c
Parent(s)
Spectrum.no_observers

The number of different inclination angles for which spectra will be extracted.

Type
Integer
Values
Greater than 0
File
setup.c
Parent(s)
Child(ren)
Spectrum.orbit_phase

For binary systems, the orbital phase at which the spectrum is to be extracted (so the effects of an eclipse can be taken into account in creating the spectrum). Phase 0 corresponds to inferior conjunciton, that is with the secondary in front (or depending on inclination angle, partially in front of) the primary

Type
Double
Values
Between 0 and 1
File
setup.c
Parent(s)
Spectrum.select_azimuth

Advance command which along with several other parameters specifies a spherical region of space in cylindrical coordinates. This parameter desribes the azimuth of the region. When this general option is used, a detailed spectrum is constructed just from photons that originate or scatter int he region

Type
Double
Unit
Degrees
Values
Between 0, and 360 or -180 to 180
File
setup.c
Parent(s)
Spectrum.select_location

One of several related parameters that permit one to apply additional conditions on the location of photons extracted in the detailed spectrum. The location refers here to the either where the photons was created or where it last scattered

Type
Enumerator
Values
all
Select photons regardless of where they are generated
below_disk
Select only photons generated from below (-z) the disk
above_disk
Select only photons orginating above the disk
spherical_region
Select photons by defining a spherical region
File
setup.c
Parent(s)
Child(ren)
Spectrum.select_photons_by_position

Advanced command associated with adding conditions for the detailed spectra that are extracted. This command simply asks whether one would like to select photons by position. If so one will be asked to define a spheical region in interms of its cylindrical coordinates.

Type
Boolean (yes/no)
File
setup.c
Parent(s)
Child(ren)
Spectrum.select_r

Part of a set of parameters which define a spherical region of space from which photons are to be extracted. select_r defines the radius of the spherical region

Type
Double
Unit
cm
Values
Greater than 0
File
setup.c
Parent(s)
Spectrum.select_rho

Advanced command which defines a spherical region of space from which photons are to be extracted in constructing a detailed spectrum. The region is defined by a cylindrical distance, and z height and an aximuth, and a radius r. This parameter defines the rho coordiante of the region.

Type
Double
Unit
cm
Values
Condition e.g. greater than 0 or list e.g. [1, 2, 5]
File
setup.c
Parent(s)
Spectrum.select_scatters

Advaned command that allows one to extract photons that have undergone a certain number of scatters. If n > MAXSCAT, that is to say a very large number then all scatters are slected. If lies between 0 and MAXSCAT then photons will be extracted only at the point a photon has undergone this number of scatters. If n is < 0 then photons with n or greater scattters will be extracted.

Type
Integer
Values
Greater than 0
File
setup.c
Parent(s)
Spectrum.select_specific_no_of_scatters_in_spectra

Advanced command which allows one to place additional constraints on the detailed spectra which are extract. This includes selectiong photons from above or below the disk, only photons which have scttered, etc.

Type
Boolean (yes/no)
File
setup.c
Parent(s)
Child(ren)
Spectrum.select_z

Advanced command which defines a spherical region of space from which photons are to be extracted in constructing a detailed spectrum. The region is defined by a cylindrical distance, and z height and an aximuth, and a radius r. This parameter defines the z coordiante of the region.

Type
Double
Unit
cm
Values
Within the z range of the model
File
setup.c
Parent(s)
Spectrum.type

The type of spectra that are produced in the final spectra. The current choices are flambda, fnu, or basic, where basic implies simply summing up the energy packets that escape within a particularly wavelength/ frequency bin.

Type
Enumerator
Values
flambda
λF(λ)
fnu
νF(ν)
basic
F(λ)
File
setup.c
Parent(s)
Spectrum.wavemax

The maximum wavelength of the detailed spectra that are to be produced

Type
Double
Unit
Angstroms
Values
Spectrum.wavemin
Greater than
File
setup.c
Parent(s)
Spectrum.wavemin

The minimum wavelength of the final spectra in Angstroms

Type
Double
Unit
Angstroms
Values
Greater than 0
File
setup.c
Parent(s)
Other

Todo

Fill in

Diag

Todo

Fill in

Diag.adjust_grid

Choose whether or not you would like to adjust the scale length for the logarithmic grid. Advanced command.

Type
Boolean (yes/no)
File
setup_domains.c
Parent(s)
Child(ren)
Diag.extra

Decide whether or not to use extra diagnostics in advanced mode. If set, this triggers a many extra questions that allow one to investigate things such as photon cell statistics, the velocity gradients in cells and the resonant scatters in the wind

Type
Boolean (yes/no)
File
python.c
Child(ren)
Diag.fractional_distance_photon_may_travel

The distance photon may travel in a cell is limited to prevent a photon from moving such a long path that the velocity may change non-linearly. This problem arises primarily when the photon is travelling azimuthally in the grid. This changes the default for the fraction of the maximum distance in a cell.

Type
Double
Values
0 to 1
File
diag.c
Parent(s)
Diag.keep_ioncycle_windsaves

Decide whether or not to keep a copy of the windsave file after each ionization cycle in order to track the changes as the code converges. Produces files of format python01.wind_save and so on (02,03…) for subsequent cycles.

Type
Boolean(yes/no)
File
diag.c
Parent(s)
Diag.keep_photoabs_in_final_spectra

This advanced options allows you to include or exclude photoabsorpiotn in calculating the final spectra. (but ksl does not know what the default is)

Type
Boolean (yes/no)
File
diag.c
Parent(s)
Diag.lowest_ion_density_for_photoabs

For efficiency reasons, Python does not try to calculate photoabsorption for an ion with an extremly low density. This advance parameter changes this density limit

Type
Double
Unit
n/cm**3
Values
Greater than 0
File
diag.c
Parent(s)
Diag.make_ioncycle_tables

Multi-line description, must keep indentation.

Type
Boolean (yes/no)
File
diag.c
Parent(s)
Diag.print_dvds_info

Print out information about the velocity gradients in the cells to a file root.dvds.diag.

Type
Boolean (yes/no)
File
diag.c
Parent(s)
Diag.save_cell_statistics

Choose whether to save the statistics for a selection of save_cell_statistics. If yes, it looks for a file called “diag_cells.dat” which contains the cells to track, and saves the photon details (weights, frequencies) for those that interact in the cell. Useful for checking the detailed MC radiation field in a cell.

Type
Boolean (yes/no)
File
diag.c
Parent(s)
Diag.save_extract_photons

Multi-line description, must keep indentation.

Type
Boolean (yes/no)
File
diag.c
Parent(s)
Diag.save_photons

Multi-line description, must keep indentation.

Type
Boolean (yes/no)
File
diag.c
Parent(s)
Diag.track_resonant_scatters

Multi-line description, must keep indentation.

Type
Boolean (yes/no)
File
diag.c
Parent(s)
Diag.use_standard_care_factors

Advanced command which allows one to change various other defaults associated with radiative transfer, inclusing the fractional distance in a cell that a photon can travel

Type
Boolean (yes/no)
File
diag.c
Child(ren)
Diag.write_atomicdata

Choose whether to write the atomic data that is being used to an output file.

Type
Boolean (yes/no)
File
setup_domains.c
Photon_sampling

Todo

Fill in

Photon_sampling.approach

Choice of whether and how to use stratified sampling in creating photons during the ionization stage of the calculation.

Type
Enumerator
Values
T_star
Sets a single band based on the temperature given
cv
Traditional cv setup
yso
YSO setup
AGN
Test for balance matching the bands we have been using for AGN runs
min_max_freq
Mode 1 sets a single wide band defined by f1 and f2
user_bands
User-defined bands
cloudy_test
Set up to compare with cloudy power law table command note that this also sets up the weight and photon index for the PL, to ensure a continuous distribution
wide
Test for balance to have a really wide frequency range
logarithmic
Generalized method to set up logarithmic bands
File
bands.c
Child(ren)
Photon_sampling.band_boundary

When the user specifies what bands are used for stratfied sampling, this parameter specifies the boundaries between energy bands in which a minimum fraction of photons will be generated. The number of times this parameter is request depends upon the number of energies bands being used.

Type
Double
Unit
eV
Values
Greater than 0, monotonically increasing
File
bands.c
Parent(s)
Photon_sampling.band_min_frac

When specifying manually the bands used for generating photons during the ionization phase, this parameter specifies the The minimum fraction of photons to be generated in this energy band. The number of times this parameter will be reqested depends upon the number of bands. The summ of the fractions need not sum to 1, in which case the remaining photons will be distributed according to the luminosity in the energy bands

Type
Double
Values
Greater than 0 and should sum to less than 1 over all bands
File
bands.c
Parent(s)
Photon_sampling.high_energy_limit

Stratified sampling is used during ionization cycles to generate photons. This parameter specifies the high energy limit for the frequencies of photons to be generated.

Type
Double
Unit
eV
Values
Greater than Photon_sampling.low_energy_limit
File
bands.c
Parent(s)
Photon_sampling.low_energy_limit

During the ionization phase, stratified sampling is used to provide good coverage of the full ionizing spectrum. This parameter sets the lowest envergy (frequency) of for phtoons to be generated whne the user wants to customize the bands.

Type
Double
Unit
eV
Values
Greater than 0
File
bands.c
Parent(s)
Photon_sampling.nbands

Python uses stratified samplign to generate photons during the ionization phase. This parameter allows the user to define the number of bands for stratified sampling, if s/he wants to customize the bands used for the generation of photons

Type
Integer
Values
Greater than 0
File
bands.c
Parent(s)
Child(ren)
Reverb

Todo

Fill in

Reverb.angle_bins

Used when generating 3d .vtk output files for visualisation. Sets the number of angle bins used in the output. Aesthetic only; bigger makes prettier meshes with larger filesizes.

Type
Integer
Values
Greater than 0
File
setup_reverb.c
Parent(s)
Reverb.disk_type

Setting for how photons generated in the disk are treated when generating path distributions for wind cells.

Type
Enumerator
Values
correlated
This mode assumes that disk emission is correlated with the central source. Photons generated in the disk start with a delay equal to the direct distance to the central source. This assumes that the ionisation state and luminosity of the disk surface layer is mostly determined by unscattered photons from the central source.
uncorrelated
This mode generates photons with a delay of 0 wherever in the disk they come from. This mode is of slightly questionable use and should be ignored in preference to 0 or 2. It will, in practise, generally work out similar to type 0 as most of the UV photons are generated close-in to the CO.
ignore

This mode assumes that disk photons do not correlate with the central source (i.e. disk surface ionisation state and emissivity is driven not by irradiation from the CO but by the mass inflow). This means that whilst they contribute to heating the wind, they do not strongly contribute to the lags for a given line. Photons generated by the disk do not contribute to the path distributions in the wind in this mode.

By removing the (generally) short-delay disk photons from the wind path distributions, this will slightly bias them towards the longer delays associated with wind self-heating/excitation.

File
setup_reverb.c
Parent(s)
Reverb.dump_cell

Position for a cell, listed as a pair of R:Z coordinates. Will accept any position that falls within a grid, will error out on ones that don’t. This can be slightly awkward and you may want to run a quick test then use py_wind to idenfity where wind locations are.

Type
Float:Float
Unit
cm:cm
Values
>0:>0
File
setup_reverb.c
Parent(s)
Reverb.dump_cells

Number of cells to dump. When dumping the path distribution info for a range of cells, this specifies the number of lines of Reverb.dump_cell that will be provided.

Type
Integer
Values
Greater than or equal to 0
File
setup_reverb.c
Parent(s)
Child(ren)
Reverb.filter_line

Line number of one line to include in the output .delay_dump file. This is the python internal line number. It can be found using either the macro-atom mode (which prints out the line number once it’s found one) or by doing an exploratory run with Reverb.filter_lines = -1, then looking through the delay dump file for photons of the right wavelength to see what their line is. This should almost certainly be changed to be specified using a species and wavelength!

Type
Integer
Values
Any valid line index
File
setup_reverb.c
Parent(s)
Reverb.filter_lines

Whether or not to filter any lines out of the output file. This is used to keep output file sizes down, and avoid them overwhelming the user.

Type
Int
Values
0

No filtering

Include all photons that contribute to the spectra in the output file. Not recommended as it leads to gargantuan file sizes.

-1

Filter continuum

Include all photons whose last interaction was scatter or emission in a line. Recommended setting for exploratory runs where you’d like to identify which lines are the easiest to process.
N

Filter lines

Include N Reverb.filter_line entries, each specifying one line to keep in the output file. If Reverb.matom_lines is >0, all macro-atom lines of interest are automatically included in the filter list.

File
setup_reverb.c
Parent(s)
Child(ren)
Reverb.matom_line

Specifies a line associated with a given macro-atom transition. The species and transition involved are specified. The internal line associated with this transition will be printed to standard-out for use when processing outputs. A line is specified as Element:Ion:Upper level:Lower level.

Type
Int:Int:Int:Int
Values
>0:>0:>1:>0
File
setup_reverb.c
Parent(s)
Reverb.matom_lines

Number of macro-atom lines to track paths for individually. This many reverb.matom_line entries are required, and the line associated with each has the path of photons deexciting into it recorded in its own array. Note: This doesn’t give rise to any noticable differences to the pure wind mode in most simulations.

Type
Integer
Values
Greater than or equal to 0
File
setup_reverb.c
Parent(s)
Child(ren)
Reverb.path_bins

Number of bins for photon paths. Reverb modes that record the distribution of path lengths in every wind cell bin them in this number of bins. Bins are logarithmically spaced between the minimum scale in the system (the smallest ‘minimum radius’ in any domain) and the 10 * the maximum scale in the system (10 * the ‘maximum radius’ in any domain). Default value is 1000, going much higher does not lead to qualitative differences in TF, going lower makes the bin boundaries show up in the TF.

Type
Integer
Values
Greater than 0
File
setup_reverb.c
Parent(s)
Reverb.type

Whether to perform reverberation mapping. Reverberation mapping tracks the path of photons emitted in the simulation as they travel through the geometry, assuming that any delays from recombination etc. are negligible and all delays are due to light travel time. For each final spectrum, all contributing photons are output to a ‘.delay_dump’ file that can then be processed using our ‘tfpy’ Python (no relation) library.

Type
Enumerator
Values
none
Off
photon
Each photon is assigned an initial path based on its distance from the central source (assuming emission in the disk and wind is correlated with emission from the CO).
wind
CO photons are assigned paths as in Photon mode, disk photons are assigned paths as set by the reverb.disk_type parameter. Photons generated in the wind are assigned a path based on the distribution of paths of photons that have contributed to continuum absorption in that cell.
matom

This works as wind mode, but for a number of specified macro-atom lines paths are tracked for those photons who cause a deexcitation into a given line. When a photon is emitted in one of those lines, the path is drawn from that specific distribution. This distribution is build up not just from the last cycle of the simulation, but from all cycles after the wind achieves >90% convergence. This is necessary as some lines are poorly-sampled.

This mode gives pretty much identical results to wind mode, but at least we made it to check rather than just assuming it would be fine.

This requires that Line_transfer is either macro_atoms or macro_atoms_thermal_trapping

File
setup_reverb.c
Child(ren)
Reverb.visualisation

Which type of visualisation to output, if any. Reverb modes that keep arrays of photon paths per cell can output them either as averages in a 3d model, or as a selection of flat text files with full bin-by-bin breakdowns. Useful for diagnostics.

Type
Enumerator
Values
none
No visualisation.
vtk
Mesh visualisation. Outputs mean incident path per cell, photon count per cell, and mean observed delay to ‘.vtk’ format, readable using a range of programs including (my preferred option) VisIt, available at https://visit.llnl.gov/.
dump
Outputs distributions of paths for continuum heating and each line to a range of ‘dump cells’ specified by X & Z position.
both
Outputs both vtk and dump.
File
setup_reverb.c
Parent(s)
Child(ren)
geo

Todo

Fill in

geo.xlog_scale

Choose the logarithmic scale length for the grid in the x-direction.

Type
Double
Unit
cm
File
setup_domains.c
Parent(s)
geo.zlog_scale

Choose the logarithmic scale length for the grid in the z-direction.

Type
Double
Unit
cm
File
setup_domains.c
Parent(s)
Input_spectra.model_file

In addition to being able to generate several types of spectra, such as blackbodies and power laws, Python can read in a series of spectra which are tabulated and are describable in terms of (usually) temperature and gravity). This parameter defines the name of the file which gives the location of the individual spectra and the temperate and gravity associated with each spectrum. (One may wish to use the same files for several radiation sources, viz the disk and the star) Python actually only reads in the data the first time.

Type
String
File
setup.c
Parent(s)

Outputs

Todo

Fill in

Diagnostic files

Python logs a considerable amount of information as it runs. Some of this information is printed to the screen but a much more voluminous version of progress of the program is placed in a sub- directory, named diag_whatever, where whatever is the root name of the model being run.

In this directory one will find log files, e.g. whatever_0.diag, whatever_1.diag, … where the in a multiprocessor run, the number refers to the log of a specific thread.

Inspecting these logs is important for understanding whether a Python run is successful, and for understanding why if failed if that is the case.

Model

As Python is run, it repeatedly writes out two binary files that contain essentially all information about the wind as calculated in the ionization phase of the program, along with status of the program at the last point where the file was written. These files along with the parameter file are sufficient to restart the program, if for example, one wants to check point the program after a certain time, and restart where one left off, or to add spectral cycles to get better spectra.

.wind_save
A binary file that contains essentially all information about the wind including ion densities, temperatures, and velocities in each cell, along with status of the program at the last point where the file was written.
.spec_save
A binary file that contains all of the information about the spectra that have created. This file is not of interest to users directly. It is used when restarting

Two routines exist as part of the Python distribution allow the user to gain insight into the actual model

windsave2table

Executed from the command line with windsave2table rootname.

Produces a set of standard set ascii tables that that show for each grid cell quantities such as wind velocity, \(n_e\), temperatures, and densities of prominent ions.

py_wind

Executed from the command line with py_wind rootname

Allows the user to query for information about the model interactively. The results can be written to ascii files for future reference

Spectra Files

Python is intended to produce simulated spectra. These spectra are all ascii tables intended to be accessible with software packages such as astropy.

All of the ascii begin with commented headers that contain all of the parameters of associated with a run, along with the date of the run and the specific version of Python used to make the run. In principle, if one still has access to any of the spectra, one can reproduce the entire run again.

Broad band spectra are created from the last ionization cycle. Detailed are calculated from all of the spectral cycles.

For a model with root name cv, the following broadband spectra will be created:

  • cv.spec_tot - various spectra
  • cv.log_spec_tot
  • cv.spec_tot_wind
  • cv.log_spec_tot_wind
File types
.spec_tot

An ascii file that contains various spectra from the ionization-calculation phase of the program on a linear frequency scale. The first few lines of the file (omitting the header) are as follows:

# Freq.        Lambda  Emitted  Star+BL  Disk     Wind     HitSurf  Scattered
3.023938e+14 9913.975 1.07e+33 2.03e+31 1.05e+33 1.05e+30 4.11e+31        0
3.049952e+14 9829.418  1.1e+33 2.24e+31 1.07e+33 3.97e+30 4.42e+31        0
3.075965e+14 9746.292 1.09e+33  2.1e+31 1.07e+33 1.22e+30 3.63e+31        0
3.101978e+14 9664.559 1.11e+33 1.97e+31 1.09e+33 1.33e+30 4.34e+31        0
3.127991e+14 9584.186 1.08e+33 2.03e+31 1.06e+33 1.27e+30 4.75e+31        0

The first two columns are fairly obvious. Lambda is in Angstroms. The remainder indicate the luminosity of the system in specific bands. Emitted is the total emergent spectrum, Star+BL is the emergent spectrum from photons bundles originating on the Star or BL, Disk and Wind are the same for photons originating in the disk and wind respectively. HitSurf represents photons that did not escape the system but ran into a boundary, and Scattered are photons that somewhere along their path out of the system were actually scattered.

.log_spec_tot
An ascii file which contains the same information as .spec_tot, but with a logarithmically space frequency intervals. This gives better sampling of the SED in a lot of cases and is much better for plotting things such as the input spectrum.
.spec_tot_wind
Identical to .spec_tot but just including photons that were generated in the wind or scattered by the wind
.log_spec_tot_wind
A logarithmic version of .spec_tot_wind
.spec

an ascii file that contains the final detailed spectra for the wavelengths of interest at a distance of 100 pc.

Photons bundles are generated in cycles in python and the .spec file is actually written out at the end of each cycle as the program is running in the spectrum-generation phase of the program. So one can inspect the spectrum as it is building up.

The beginning of the file (omitting the header) is as follows:

Freq.        Lambda    Created  Emitted   CenSrc     Disk     Wind  HitSurf Scattered A10P0.50 A28P0.50 A45P0.50 A62P0.50 A80P0.50
1.620713e+15 1849.757  3.8401e-12 3.6348e-12 9.1429e-14 3.5434e-12          0 8.8693e-14 1.7753e-13 9.2741e-12 7.6342e-12 6.3434e-12 2.3932e-12  9.382e-13
1.620925e+15 1849.514  4.8471e-12 4.7931e-12 2.7382e-13 4.4306e-12 8.8704e-14 1.8213e-13 2.4885e-13 1.0177e-11 7.7666e-12 3.2906e-12 3.4296e-12 1.3389e-12
1.621138e+15 1849.272  5.3058e-12  5.182e-12 9.1477e-14 4.9992e-12 9.1404e-14  2.674e-13 3.5847e-13 1.2354e-11 6.9236e-12 5.9863e-12 3.3748e-12 1.7905e-12
1.621351e+15 1849.029  3.9346e-12 3.9028e-12          0 3.8127e-12 9.0124e-14 8.9142e-14 2.6728e-13 1.1158e-11 6.4932e-12 5.1452e-12 3.9074e-12 8.1597e-13

where the first line indicates the version of python used to generate the spectrum, the second gives a brief description of each column, and the remainder of the file is the spectrum. The most important columns are 1 and 2, which are respectively the frequency and wavelength and the columns that begin with, which give the spectrum that would be observed from the object at various inclination angles and orbital phases (The number of columns is therefore variable depending on how many inclinations one asked for in the input files). The other columns Emitted, Star+BL, Disk, Wind, HitSurf and Scattered have approximately the same meaning as in the .spec_tot file.

The remaining files parallel those generated for the broadband spectra.

Spectra Files

Python is intended to produce simulated spectra. These spectra are all ascii tables intended to be accessible with software packages such as astropy.

All of the ascii begin with commented headers that contain all of the parameters of associated with a run, along with the date of the run and the specific version of Python used to make the run. In principle, if one still has access to any of the spectra, one can reproduce the entire run again.

Broad band spectra are created from the last ionization cycle. Detailed are calculated from all of the spectral cycles.

For a model with root name cv, the following broadband spectra will be created:

  • cv.spec_tot - various spectra
  • cv.log_spec_tot
  • cv.spec_tot_wind
  • cv.log_spec_tot_wind
File types
.spec_tot

An ascii file that contains various spectra from the ionization-calculation phase of the program on a linear frequency scale. The first few lines of the file (omitting the header) are as follows:

# Freq.        Lambda  Emitted  Star+BL  Disk     Wind     HitSurf  Scattered
3.023938e+14 9913.975 1.07e+33 2.03e+31 1.05e+33 1.05e+30 4.11e+31        0
3.049952e+14 9829.418  1.1e+33 2.24e+31 1.07e+33 3.97e+30 4.42e+31        0
3.075965e+14 9746.292 1.09e+33  2.1e+31 1.07e+33 1.22e+30 3.63e+31        0
3.101978e+14 9664.559 1.11e+33 1.97e+31 1.09e+33 1.33e+30 4.34e+31        0
3.127991e+14 9584.186 1.08e+33 2.03e+31 1.06e+33 1.27e+30 4.75e+31        0

The first two columns are fairly obvious. Lambda is in Angstroms. The remainder indicate the luminosity of the system in specific bands. Emitted is the total emergent spectrum, Star+BL is the emergent spectrum from photons bundles originating on the Star or BL, Disk and Wind are the same for photons originating in the disk and wind respectively. HitSurf represents photons that did not escape the system but ran into a boundary, and Scattered are photons that somewhere along their path out of the system were actually scattered.

.log_spec_tot
An ascii file which contains the same information as .spec_tot, but with a logarithmically space frequency intervals. This gives better sampling of the SED in a lot of cases and is much better for plotting things such as the input spectrum.
.spec_tot_wind
Identical to .spec_tot but just including photons that were generated in the wind or scattered by the wind
.log_spec_tot_wind
A logarithmic version of .spec_tot_wind
.spec

an ascii file that contains the final detailed spectra for the wavelengths of interest at a distance of 100 pc.

Photons bundles are generated in cycles in python and the .spec file is actually written out at the end of each cycle as the program is running in the spectrum-generation phase of the program. So one can inspect the spectrum as it is building up.

The beginning of the file (omitting the header) is as follows:

Freq.        Lambda    Created  Emitted   CenSrc     Disk     Wind  HitSurf Scattered A10P0.50 A28P0.50 A45P0.50 A62P0.50 A80P0.50
1.620713e+15 1849.757  3.8401e-12 3.6348e-12 9.1429e-14 3.5434e-12          0 8.8693e-14 1.7753e-13 9.2741e-12 7.6342e-12 6.3434e-12 2.3932e-12  9.382e-13
1.620925e+15 1849.514  4.8471e-12 4.7931e-12 2.7382e-13 4.4306e-12 8.8704e-14 1.8213e-13 2.4885e-13 1.0177e-11 7.7666e-12 3.2906e-12 3.4296e-12 1.3389e-12
1.621138e+15 1849.272  5.3058e-12  5.182e-12 9.1477e-14 4.9992e-12 9.1404e-14  2.674e-13 3.5847e-13 1.2354e-11 6.9236e-12 5.9863e-12 3.3748e-12 1.7905e-12
1.621351e+15 1849.029  3.9346e-12 3.9028e-12          0 3.8127e-12 9.0124e-14 8.9142e-14 2.6728e-13 1.1158e-11 6.4932e-12 5.1452e-12 3.9074e-12 8.1597e-13

where the first line indicates the version of python used to generate the spectrum, the second gives a brief description of each column, and the remainder of the file is the spectrum. The most important columns are 1 and 2, which are respectively the frequency and wavelength and the columns that begin with, which give the spectrum that would be observed from the object at various inclination angles and orbital phases (The number of columns is therefore variable depending on how many inclinations one asked for in the input files). The other columns Emitted, Star+BL, Disk, Wind, HitSurf and Scattered have approximately the same meaning as in the .spec_tot file.

The remaining files parallel those generated for the broadband spectra.

Issues with Generating Spectra

With the current machinery to create spectra, it is possible to come across the situation where models with large optical depth or wind velocities will generate spectra with different flux normalisation depending on the wavelength range.

This problem was originally encountered whilst modelling Tidal Disruption Events. Two spectra for the same model were generated over two wavelength ranges; a restricted (1100 - 2600 A) and a broader (500 - 5000 A) range. The problem encountered was that the broad range spectrum had more flux than the spectrum with the restricted range. The figure below shows the same model, but over two wavelength ranges - as well as two spectra where the maximum number of scatters a photon can undergo is changed,

  • tde_flux_small_range: The restricted wavelength range
  • tde_flux_large_range: The broad wavelength range
  • tde_flux_small_range_maxscat: The restricted wavelength range with a value MAXSCAT = 50
  • tde_flux_no_maxscat: The restricted wavelength range with no MAXSCAT limit
_images/spectrum_generation_large_optical_depth.png

Example spectra showing differing flux totals

The problem here is not caused by a bug with the code, but is a consequence of the large wind velocities and optical depths of the model. We currently believe that there are two reasons why the flux differs between these two wavelength ranges.

Doppler Shifting out of the Spectrum Wavelength Range

At the edges of the restricted spectrum above, the flux is reduced. This is due to photon frequencies being shifted outside of the wavelength range of the spectrum. If a significant number of photons are removed from the spectrum in this way, then the following Error is printed,

spectrum_create: Fraction of photons lost: 0.10 wi/ freq. low, 0.19 w/freq hi

This tells one the fraction of the photon sample which does not contribute towards the spectrum due to to the photon frequencies being larger or smaller than the defined spectrum range, due to Doppler shifting. In models with large wind velocities (0.2 - 0.5 c) and a small spectral range, the fraction of photons lost is large and the flux at the edge of generated spectra is reduced - as can be seen above in the above figure. However, when the wind has a more moderate velocity, the number of photons lost due to being shifted out of the range is much lower and does not produce a noticeable effect on the flux normalisation of the spectra.

Removing Photons due to Too Many Scatters

As well as edge effects, flux can be lost due to photons being removed from the photon sample due to scattering too many times. In Python, when a photon has undergone MAXSCAT = 500 scatters, a photon is assumed to have become stuck in the wind and hence it is terminated and no longer tracked.

In models with large optical depths, the number of photons terminated in this way can become large. During spectrum generation, these photons will never fully escape the system but will only contribute partially to the spectrum due to extract - they will never contribute if Live or Die is used instead.

At current, there is no logic to detect this and hence no error is given. However, it is often insightful to read the output from the Photons contribution to the various spectra table, as shown below,

Photons contributing to the various spectra
Inwind   Scat    Esc     Star    >nscat    err    Absorb   Disk    sec    Adiab(matom)
     0       0    3455       0       0       0       0       0       0       0
     0       0    3455       0       0       0       0       0       0       0
     0       0     427       0       0       0       0       0       0       0
     0       0    1598       0       0       0       0       0       0       0
     0       0    1430       0       0       0       0       0       0       0
     0       0       0       0       0       0       0       0       0       0
     0       0    3313       0   17209       0     169       0       0       0
     0  282914  487831    1474       0       0     129  223605       0       0
     0  441756  336514    2223       0       0     135  215325       0       0
     0  609395  180082    5677       0       0      94  200705       0       0
     0  750672   61308   12010       0       0      59  171904       0       0
     0  838923   26143   29057       0       0      55  101775       0       0

In the above table, one can see that 17,209 photons which scattered more than MAXSCAT times contributed to the the scattered spectrum, suggesting that a large number of photons were terminated due to too many scatters.

Note

The photon numbers presented in this table are only for the master MPI process. Hence, if running in multiprocessor mode, the number here will never equal the total number of photons in the simulation, but only the number of photons in the current process.

GitHub Issue

The original GitHub issue discussing this problem can be found here; #471.

Atomic Data

The purpose of this document is to clarify where the atomic data used in Python comes from and how to translate the raw data into data which is read. All data used by Python is accessed by routines in get_atomic_data.c and stored in structures.

Bound-bound electron collision strengths

Source

We use the Chianti atomic database, specifically the *.scups files. These “contain the effective electron collision strengths scaled according to the rules formulated by Burgess & Tully 1992, A&A, 254, 436 The values in the file are functional fits to \(\Upsilon(T)\)Omega` in our calculations for collisional de-excitation rate coefficient

\(q_{21}=\Omega\frac{8.629\times10^{-6}}{g_u\sqrt{T}}\)

In the g-bar formulation

\(\Omega=4.77\times10^{16}g_l\overline{g}\frac{f_{abs}}{\nu}\)

These values of \(\Upsilon\) simply replace \(\Omega\).

Translation to Python format

It is necessary to link each line in our line list with the relevant electron collision strength. This is achieved using the python script “coll_stren_lookup.py” which first reads in the “lines_linked_ver_2.py” line list, then attempts to work out which lines are which by comparing the energy and the oscillator strength of the line. If these match to within a factor of 10% then the code logs this as a possible match. If better matches come along, then the code adopts those instead.

Each matched line get a line in the data file which is basically all of the line data for the matched line. This is to give Python the best chance of linking it up with the line internally.

Data format

Each record has three lines. The first line has a keyword CSTREN and this contains all the line data for the line to which this collision strength refers as the first 10 fields. These fields are identical to the 10 fields that appear in a standard line file. The next 8 fields are

  • 1-2 Upper and lower level of transition - Chianti nomenclature
  • 3 Energy of transition - Rydberg
  • 4 Oscillator strength x lower level multiplicity (GF)
  • 5 High temperature limit value
  • 6 Number of scaled temperatures - 5 or 9
  • 7 Transition type cite{1992A&A…254..436B} nomenclature
  • 8 Scaling parameter (C) (Burgess & Tully 1992) nomenclature

The next two lines, with labels SCT and SCTUPS are the 5 or 9 point spline fits to \(\Upsilon\) vs T in reduced units y,x.

There are four different types of transitions, each with their own scaling between temperature of the transition and $Upsilon$

For example, for type 1 (the most common)

\(x=1-\frac{lnC}{ln\left(\frac{kT}{E_ij}+C\right)\)

and

\(y(x)=\frac{\Upsilon}{ln\left(\frac{kT}{E_{ij}}+e\right)\)

So, to get \(\Upsilon\) for a given T, one converts T to x via the correct equation,then linearly interpolate between values of y(x), then convert back to $Upsilon$

Python structure

The data is stored in Python in the Coll_stren structure which has memebers

  • int n - internal index
  • int lower,upper - the Chianti levels, not currently used
  • double energy - the energy of the transition
  • double gf - the effective oscillator strength - just oscillator strength x multiplicity
  • double hi_t_lim - the high temperature limit of y
  • double n_points -The number of points in the spline fit
  • int type - The type of fit, this defines how one computes the scaled temperature and scaled coll strength
  • float scaling_param - The scaling parameter C used in the Burgess and Tully calculations
  • double sct[N_COLL_STREN_PTS] -The scaled temperature points in the fit
  • double scups[N_COLL_STREN_PTS]- The sclaed coll sttengths in ythe fit

There is also a member in the line structure (coll_index) which points to the relevant record

comments

This data has been generated to match the Verner line list. If we want to use the Kurukz line list, then a new set of data should be made

Bound Bound (Kurucz)

This is the data for computing bound bound, or line interactions in simple atoms. There are two main sources of dtaa currently used in Python. This page discusses the “Kurucz” data

Source

The Kurucz data used to create simple lines was taken from the Kurucz website. The file gfall.dat was used, though a few extra lines known to have been missing were likely added.

Translation to Python format

There are several steps to creating the data used in Python from that in gfall.dat, that are carried out by py_read_kurucz and py_link. The first routine reads the gfall.dat file and creates two output files, a file containing the lines and the associated such as the effective oscillatory strength and a file which contains information about the ion levels. py_read_kurucz chooses only a portion of the Kurucz lines, namely those associated with ions with ionization potentials in a certain range and lines with gf factors exceeding a certain value. The second program py_link attempts to create a model ion with links between the levels and the ions. Both of these routines are driven by .pf files, similar to what are used in python. Examples of the .pf files are in the directory py_kurucz

Data format
Label Element Ion \(\lambda(\AA)\) f \(g_l\) \(g_u\) \(e_l\) \(e_u\) \(index_l\) \(index_u\)
Line 1 1 926.2 0.003 2 4 0.0 13.4 0 9
Line 1 1 930.7 0.005 2 4 0.0 13.3 0 8
Line 1 1 937.8 0.008 2 4 0.0 13.2 0 7

where f is the oscillator strength of the transition, the g nuumbers are the multiplicities of the upper and lower levels, the e values are the energy levels relative to ground state of the levels in eV and the two indices relate to the levels structure.

Python structure

Where the data is stored internally in Python

Comments

The version of gfall.dat that is used in Python is out of date, though whether this affects any of the lines actually used in python is unclear. The website listed above has a link to newer versions of gfall.dat.

Bound Bound (Verner)

This is the data for computing bound bound, or line interactions in simple atoms. There are two main sources of dtaa currently used in Python. This page disucsses the “Verner” data

Source

Give original references

Translation to Python format

Similar to the Kurucz list there is a program py_read_verner that can produce a line list from the Verner data. Once done then py_link can be used to link the the levels with the lines.

Data format

The data format is identical to that of produced for the Kurucz data

Python structure

Where the data is stored internally in Python

Comments

In practice we have not used these data for any Python publications. At some point early in the AGN project, NSH increased the number of lines, and generated lines_linked_ver_2.py and levels_ver_2.py. I think this was because there was a small bug which meant the oscillator strength cut that was stated was not that which was applied.

Bound Free Macro-atom (empty)

Source
Translation to Python format
Data format
Python structure
Comments

Bound Free (Topbase)

Source

Obtained from The Opacity project. See also Cunto et at 1993, A&A, 275, L5

Translation to Python format

ksl - It’s not clear that we are now making use of the topbase data in this way but my original attempt to incorporate topbase photoinization data into Python is contained in the directory topbase. Processing of these files was done by py_top_phot. My feeling is that we can replace these remarks with those that are more up to date, once Nick and James discuss this section, and delete any mention of my original attempt on this in the data-gen archive.

Data format

Explain the ascii format of the file which is read into Python

Python structure

Where the data is stored internally in Python

Comments

Extrapolation to higher energies

Around about Python 76 we discovered that some topbase cross-sections have maximum energy cutoffs, for no good reason. This causes unrealistic edges to appear in spectra. To counter this, JM wrote a script to extrapolate the cross-section to higher energies. This is done by calculating the gradient in log-space at the maximum energy and extrapolating to 100 keV. A number of cross-sections had unrealistic gradients at the original maximum energy, and were identified by eye and then forced to have a \(\nu^{-3}\) shape. This is the shape of a hydrogenic cross-section and whilst it is not accurate for non-hydrogenic ions, it is more realistic (and conservative) than some of the unphysically shallow gradients that were being found. This is also briefly described in section~3.7.2 of Matthews PhD thesis. The python scripts can be found in progs/extrapolate_xs/ with docstrings describing their use.

Bound Free (Verner)

This is data for bound free or photoionization data. There is information for both inner shell (auger) and outer shell PI.

Source

There are three sources for this data

Translation to Python format

Tabulation Process

The raw VFKY data comes in a series of fit parameters. We decided, circa Python 78, to tabulate this data instead. Partly, I think I because the on the fly method was time consuming (yes, computing all the pow() commands to commute the cross sections on the fly took a huge amount of time) and we decided that tabulating pre program was better than doing it in the program, so that everything was of the same format.

The script which does this is progs/tabulate_xs/photo_xs.py – it creates a file like photo_vfky_tabulated.data.

Inner and Outer Shells

For the ground states, we split the cross sections up into outer shell and inner shell cross sections. This allows us to calculate possible auger ionization as ions relax after an inner shell ionization. This is done using the python script “verner_2_py.py. This script takes the normal verner cross sections, which truncate at the first inner shell edge and firstly appends the outer shell data from VY to that to make a full outer shell cross section. These are written out into “vfky_outershell_tab.data” It then writes out the inner shell cross sections into “vfky_innershell_tab.data”. There is a lot of complicated machinery to try and work out the exact shell that is being ionized to allow these rates to be linked up to the relevant electron yield (and flourescent) records.

Data format

Explain the ascii format of the file which is read into Python

VFKY_outershell_tab.data

Label z state islp level threshold_energy n_points
PhotVfkyS 1 1 -999 -999 1.360e+01 100

This data is linked to the relevant ion via z and state, islp and level are not used. the last number n_points, says how many points are used for the fit, and the next n_points lines in the file, each preceeded by the label PhotVfky are pairs of energy (in eV) vs cross section which make up that fit.

VY_innershell_tab.data

label z state n_shell l_subshell threshold_energy n_points
InnerVYS 3 1 1 0 6.439e+01 100

This data is linked to the relevant ion via z and state. the n_shell and l_subshell numbers are used to cross reference to the electron yield records. As above, the last record shows how many points are in the fit, and the data pairs making up the fit follow with the keyword InnerVY.

Python structure

Where the data is stored internally in Python

Comments

The manner in which this data is read into Python is a bit labyrinthine at the moment. The intention is to use a combination of VFKY and VY for all ground states, an

Direct Ionization

This is the data to compute ionization rates from collisions between ions and hot electrons.

Source

The data comes directly from Dere 2006, A&A, 466, 771 . This paper gives direct ionization and excitation-autoionization rate coefficients for many ions as a function of temperature for Maxwellian electron distributions.

Translation to Python format

The data table is downloaded in its entirety from the data table associated with the paper. All that happens is that the table is saved to a text file, and the keyword DI_DERE is just prepended to each row.

Data format

Each line starts with the label DI_DERE and then follows

  • Nuclear Charge - z - used to identify the ion
  • Ion - state in our normal notation, so 1=neutral
  • Number of splines N- the number of spline points for the fit of rate coefficients vs scaled temperature
  • Scaled temperatures - there are N of these
  • Scaled Rate coefficients - N of these

The scaled temperatures are given by

\(x=1-\frac{\log{f}}{\log(t+f)}\)

where t=kT/I. I is the ionization potential, and f=2.0. The rate coefficient R(T) is recovered from the scaled rate coefficient in the table, $rho$ using

\(\rho=t^{1/2}I^{3/2}R(T)/E_{1}(1/t)\)

where \(E\_{1}\) is the first exponential integral. In python we use the gsl_sf_expint_E1 routine in gsl.

Python structure

This data is stored in the dere_di_rate structure with members

  • int nion- Internal cross reference to the ion that this refers to
  • int nspline - the number of spline points that the fit is evaluated over
  • double temps[DERE_DI_PARAMS]- temperatures at which the rate is tabulated
  • double rates[DERE_DI_PARAMS]- rates corresponding to those temperatures
  • double xi - the ionization energy of this ion
  • double min_temp -the minimum temperature to which this rate should be applied
Comments

This data is also in Chianti , although in a different form. So we could potentially use this data as part of a push to just use Chianti for all our data uses. An updated set of DI data is available here

Auger Electron Yields

This data is linked with the inner shell photoionization data. It gives probabilities for different numbers of electrons to be ejected following inner shell ionizations.

Source

This data comes from Kaastra and Mewe 1993, A&A, 97, 443 . The data is downloaded from the vizier site linked and put into a file called “electron_yield.data”

Translation to Python format

The translation takes place using the python script “kaastra_2_py.py” which takes the saved raw data file “electron_yield.data” and compares it line by line to the inner shell cross section data in “vy_innershell_tab.data”(see above). The n shell and l subshell to which each record applies is coded in the KM data and needs to be decoded. This is what the script does, and all the script then does is output the yield data into a new file “kaastra_electron_yield.data” which contains the n and l cross reference.

Data format

This is the data format of the electron yield data

Label z state n l IP(eV) <E>(eV) P(1e) P(2e)
Kelecyield 4 1 1 0 1.15e+02 9.280e+01 0 10000
Kelecyield 5 1 1 0 1.92e+02 1.639e+02 6 9994
Kelecyield 5 2 1 0 2.06e+02 1.499e+02 0 10000

The data is linked to the correct inner shell photoionization cross section (and hence rate) via z, state, n shell and l subshell. The IP is not used. <E> is the mean electron energy ejected, used to calculate the PI heating rate in radiation.c. The last ten columns in the file (2 shown in the table above) show the chance of various numbers of electrons being ejected in units of 1/10000.

Python structure

The data is stored in python in the inner_elec_yield structure which contains

  • int nion - Index to the ion which was the parent of the inner shell ionization
  • int z, istate - element and ionization state of parent ion
  • int n, l - Quantum numbers of shell
  • double prob[10] - probability for between 1 and 10 electrons being ejected
  • double I - Ionization energy
  • double Ea - Average electron energy
Comments

Elements and Ions

The first file that must be read into textsc{python} is the file that defines the elements and ions. The

Translation to python:

The original data and the translation can be found in py_verner. A simple awkscript converts the downloaded data to Python format.

Datafile - elem_ions_ver.py:

There are two sections to the file, first elements are defined:

Label z Symbol Abundance Atomic Weight
Element 1 H 12.00 1.007940
Element 2 He 10.99 4.002602

and then the ions.

Label Symbol z state g $xi$ max lev max nlte . config
IonV H 1 1 2 13.59900 1000 10 1s(2S_{1/2})
IonV H 1 2 1 1.0000e+20 0 0 Bare.
IonV He 2 1 1 24.58800 1000 10 1s^2(1S_0)$
IonV He 2 2 2 54.41800 1000 10 1s(2S_{1/2})
IonV He 2 3 1 1.0000e+20 0 0 Bare
Python structure:

This data is held in Python in various fields in structures elements and ions.

Comments:

Supernova models

Supernovae (SNe) do not have solar abundances. SS included an additional file, texttt{elem_ions_ver_sn.py} for use with SN models. This is accessed through the texttt{standard_sn_kurucz} masterfile and as far as I know is just added by hand to match expected Type Ia abundances and specifically the abundances used by Tardis.

ksl - The abundances used by Verner are not necessarily the best values today. This is one of the the items we should consider updating.

Free-Free Emission

Source

The free-free Gaunt factors are taken from Sutherland 1998, MNRAS, 300, 321. The data is available for download here where three files exist

  • gffew.dat : Free-Free Emission Gaunt factors as a function of scaled electron and photon energy.
  • gffgu.dat : Free-Free Emission Gaunt factors for Maxwellian electrons as a function of scaled temperature and scaled frequency.
  • gffint.dat : Total Free-Free Emission Gaunt factors for Maxwellian electrons.

The last file is the one we use to calculate free-free emission, since this in integrated gaunt factor over a population of electrons with a Boltzmann distribution of velocities. The other two files could be of use in the future should we wish to have gaunt factor corrections for the heating rates,in which case we should use the gffgu.dat data file. However generally speaking free-free heating is never important and there would be significant overhead in calculating a gaunt factor for each photon.

Translation to python

The file is simply modified by hand to put a label “FF_GAUNT” at the start of each data line and a hash at the start of each comment line.

Datafile - gffint.dat:

The format of the data file to be read into python is as follows:

Label \(\log(\gamma^2)\) \(<g_{ff}(\gamma^2)>\) \(s_1\) s_2$ s_3
FF_GAUNT -4.00 1.113E+00 1.348E-02 1.047E-02 -1.855E-03
FF_GAUNT -3.80 1.117E+00 1.744E-02 9.358E-03 5.565E-03
FF_GAUNT -3.60 1.121E+00 2.186E-02 1.270E-02 4.970E-03

where \(\gamma^2\) is the “scaled inverse temperature experienced by an ion” and the other four numbers allow the free-free Gaunt factor to be computed at any scaled inverse temperature

\(x=\frac{Z^2}{k_BT_e}\frac{ 2\pi^2e^4m_e}{h^2}\)

through spline interpolation between the two bracketing values of \(\log(\gamma^2)\)

\(<g_{ff}(x)>=<g_{ff}(\gamma^2)>+\Delta\left[s_1+\Delta\left[s_2+\Delta s_3\right]\right]\)

where

\(\Delta=\log(x)-\log(\gamma^2)\)

Python structure

This data is held internally in Python in the structure gaunt_total which has members

  • log_gsqrd
  • gff
  • s1, s2, s3
Comments

We currently just use the total free free emission gaunt factor as a function of temperature for a Maxwellian distribution of electrons. This is OK for cooling, however we should really use the frequency dependant gaunt factor for free free heating. If we ever have a model where free-free heating is dominant, this should be looked into.

Macro-atom Lines (empty)

Source
Translation to Python format
Data format
Python structure
Comments

Auger Photon Yields

When inner shell (or Auger) ionization takes place - there is a chance of photons being eected as the inner shells are re-filled. This data provies the information to compute the photons thus made. It is currently not used.

Source

This data comes from Kaastra and Mewe 1993, A&A, 97, 443 . The data is downloaded from the vizier site linked and put into a file called “fluorescent_yield.data”

Translation to Python format

The translation takes place using the python script “kaastra_2_py.py”. All identical to electron yield, but input file is “fluorescent_yield.data” and output is “kaastra_fluorescent_yield.data”

Data format

This is the data format of the electron yield data

Label z state n l photon_energy(eV) yield
Kphotyield 5 1 1 0 1.837e+02 6.000e-04
Kphotyield 5 1 1 0 1.690e+01 7.129e-01
Kphotyield 6 1 1 0 2.768e+02 2.600e-03

The data is linked to the correct inner shell photoionization cross section (and hence rate) via z, state, n shell and l subshell. The photon energy field is thew energy of the fluorescent photon in eV, and yield is the number of said photons emitted per ionization multiplied by \(10^4\).

Python structure

The data is stored in python in the inner_fluor_yield structure which contains

  • int nion - Index to the ion which was the parent of the inner shell ionization
  • int z, istate - element and ionization state of parent ion
  • int n, l - Quantum numbers of shell
  • double freq - the rest frequency of the photon emitted
  • double yield - number of photons per ionization x \(10^4\)
Comments

This data is not currently used

Photon Banding Strategies

Photon packets are emitted from a number of different radiation sources, such as the accretion disk or from the wind itself. When a photon is created, it is defined by its frequency \(\nu\) and weight \(w\). Photons are generated at the beginning of each cycle and can either be generated uniformly over the entire frequency range, or can be generated in pre-defined frequency bands, where certain frequency bands are biased to have more photons.

Uniform Sampling

In the most simple case, the frequency of a photon is sampled uniformly over the entire frequency range. The total weight of all photons is equal to the luminosity of the system and each photon has weight has equal weight given by,

\[w_{i} = \frac{\sum_{i = \text{sources}} ~ \int_{\nu_{\text{min}}}^{\nu_{\text{max}}} L_{\nu, i} d\nu}{N},\]

where \(N\) is the total number of photons, \(\nu_{\text{min}}\) and \(\nu_{\text{max}}\) define the frequency range and \(L_{\nu, i}\) is the luminosity for a radiation source. Note that a summation is used to find the luminosity for each radiation source and that \(w\) has units of \(\text{ergs s}^{-1}\).

Banded Sampling

In practice, uniform sampling is generally an inefficient approach to generating photon packets. For example, it is often desirable to produce a sufficient number of photons within a specific frequency range, i.e. around a photoionisation edge. However, if these frequencies lie on the Wien tail of a blackbody distribution, it is unlikely that a sufficient number of photons will be generated as most of the luminosity is generated at lower frequencies. It is possible to get around this problem by generating an increasingly large number of photons. But, this is computationally expensive and inefficient.

In order to cope with cases this, Python implements importance sampling which effectively increases the number of photons which are sampled from specific frequency bands considered important. Photons are now generated with the weight,

\[w_{j} = \frac{\sum_{j = \text{sources}} ~ \int_{\nu_{i}}^{\nu_{i + 1}} L_{\nu, j} ~ d\nu}{f_{i} N},\]

where, again, this is a summation over all radiation sources. \(N\) is the total number of photons, \(f_{i}\) is the fraction of photons emerging from frequency band \(i\), \(\nu_{i}\) and \(\nu_{i+1}\) are the lower and upper frequency boundaries for each frequency band and \(L_{\nu, j}\) is the luminosity of the radiation source. Hence, more photons from frequency bands with a larger fraction \(f_{i}\) will be generated. However, photons from important bands (where more photons are sampled from) will have reduced weight, whilst photons from the less important frequency bands will have increased weight.

This scheme has the benefit of allowing the generation of a lower number of photons whilst still sufficiently sampling important frequency ranges, decreasing the computational expense of a simulation. As such, this is the preferred sampling method.

Available Sampling Schemes

Python currently implements seven pre-defined frequency bands and and two flexible run time banding schemes. The parameter used to define the photon sampling scheme is,

Photon_sampling.approach(T_star,cv,yso,AGN,min_max_freq,user_bands,cloudy_test,wide,logarithmic)

Minimum and Maximum Wavelengths

At present, the largest wavelength a photon can be is hardwired to 20,000 Angstroms. The smallest wavelength a photon can take is defined by the temperature of hottest radiation source, but is at least 115 Angstroms - twice that of the Helium edge.

T_star

Create a single frequency band given a temperature T, which is the temperature of the hottest radiation source in the model. All photons will then be generated from this single frequency band.

CV

Pre-defined bands which have been tuned for use with CV systems, where a hot accretion disk (~100,000 K) is assumed to exist. In this scheme, there are four bands where the majority of photons are generated with a wavelength of 912 Angstroms or less.

YSO

Pre-defined bands which have been tuned for use with YSO systems. In this scheme, there are four bands.

AGN

Pre-defined which have been tuned for use with AGN system. In this scheme, there are ten bands, with a minimum frequency of \(1 \times 10^{14}\) Hz and a maximum frequency of \(1 \times 10^{20}\) Hz.

min_max_freq

Create a single band using the minimum and maximum wavelength as described by the minimum and maximum wavelengths calculated for the current model.

user_bands

This allows a user to create their own frequency bands, defined by photon energies measured in electron volts. The first band has the lowest photon energy and each subsequent band must have a larger energy than the previous band. Each band also requires a minimum fraction of photons to be sampled from this band, where the sum of the fractions for each band must be equal to or less than one.

Maximum Number of Bands

Currently, a maximum of 20 frequency bands can be defined. If a user attemps to specify more than than 20 bands, Python will create an error message and fallback to using 20 bands.

cloudy_test

This set of bands were created for use in testing against the photoionisation and spectral synthesis code Cloudy.

wide

Pre-defined bands which have very wide frequency range. The purpose of this band is for testing, hence is best to avoid using this band for a working model.

logarithmic

This is the same as user_bands, however the frequency bands are now defined in log space. This allows one to better sample a frequency range which spans many orders of magnitude.

Maximum Number of Bands

Currently, a maximum of 20 frequency bands can be defined. If a user attemps to specify more than than 20 bands, Python will create an error message and fallback to using 20 bands.

Minimum Fraction

For logarithmic user defined bands, the fraction of each band is set to 1 / nbands.

Convergence

Ionization cycles in Python are intended to establish the correct ion densities and temperature of the various cells in the wind. The degree to which this happens for a given number of ionization cycles depends on how far the initial guess of electrons temperatures in various portions of the wind and the number of photons generated during each photoionization cycles. Furthermore, the accuracy of the final model depends on the number of photons that pass through each cell. As a result, the accuracy with which ion abundances and temperature are determined will differ on a cell by cell basis. In a typical model with a biconical outflow, a small cells at the outer edge of the accretion disk will record fewer photon passages than one in the middle of the grid that is exposed to large numbers of photons from the disk.

To estimate whether a a model calculation has reached a steady state, Python carries out three tests, one comparing the difference in the radiation temperature between the current and the preceding ionization cycle, one comparing the electron temperature in the same manner and once comparing heating and cooling rates in the current cycle. More specifically, if a cell satisfies the following 3 tests,

\[\left | \frac{t_e^n-t_e^{n-1}}{t_e^n+t_e^{n-1}} \right | < \epsilon\]
\[\left | \frac{t_r^n-t_r^{n-1}}{t_r^n+t_r^{n-1}} \right | < \epsilon\]
\[\left | \frac{Heat_{tot}^n-Cooling_{tot}^{n}}{Heat_{tot}^n+Cooling_{tot}^{n}} \right | <\epsilon\]

where \(\epsilon = 0.05\), then cell is said to have converged.

The routine run_check.py produces two plots related to convergence, one showing the fraction of cells that have passed each of the tests above as a function of cycle, and the other showing the number of failed checks for each cell in the wind.

Evaluation

Determining whether Python has run successfully from a a scientific point of view depends very specifically on one’s goals. Did the spectra turn out to be what one expected? Here by evaluation we mean, did my run complete without significant errors and did the ionization structure converge to a “steady state” given the number of ionization cycles, the number of photons, and the frequency distributions of the photons we chose.

Convergence

A very basic question about a particular run is, has it reached a “steady state” and if it is in a steady state are cells stable in the sense that fluctuations are small.

Todo

Define what converged means here, and provide examples

Note that it is not always important that all cells be converged. The Monte Carlo process preferentially picks out the cells which affect the emergent radiation field. Portions of the grid which do not get many photons are typically the ones that are “not converged”, but since they don’t contribute much to the emergent radiation, one does not need them to be converged (except if one wants to make nice plots of the temperature as a function of position in the wind or of the density of a particular ion species). On the other hand, if one is using Python in conjunction with a hydrodynamical code one wants all the cells to be converged.

Errors

Python is designed to continue to run unless something catastrophic happens. As it runs, it logs error messages that can be found in the .diag files. These messages are a combinations of warnings, and/or unusual occurrences, that if they start occurring often suggest a real problem.

These error messages are all of the form:

Error: wind2d: Cell   0 ( 0, 0) in domain 0 has 1 corners in wind, but zero volume

that is they begin with the word Error. followed by the subroutine in the code where the error occurred followed by what is hopefully a helpful. If one is concerned about a particular message, one can hopefully determine what is happening by looking for that message in the log files.

Python keeps a count of the number of times a particular message has occurred and at the end of the program, and the very end of the diag files contain a listing of how many times a particular error has occurred.

Error summary: End of program, Thread 2 only
Recurrences --  Description
     7 -- getatomic_data: line input incomplete: %s
   128 -- get_atomicdata: Could not interpret line %d in file %s: %s
     1 -- error_count: This error will no longer be logged: %s
     1 -- get_wind_params: zdom[ndom].rmax 0 for wind type %d
     1 -- wind2d: Cell %3d (%2d,%2d) in domain %d has %d corners in wind, but zero volume
     1 -- check_grid: velocity changes by >1,000 km/s in %i cells
     1 -- check_grid: some cells have large changes. Consider modifying zlog_scale or grid dims

As indicated here, these are the errors for only thread 2 of a program. In order to get a summary of all the threads, there is a script py_error.py that be run as py_error.py rootname from the main run directory. Note that in many cases, the summary will be the number times an error occurred in one thread times the number of threads, but not always.

One should be aware of these errors, and watch out for situations where the number of errors of a particular type is much larger than usual.

The homolgous wind model

In the homolgous model the wind/outflow is assumed to have spherical symmetry and to have a particularly simply velocity law: specifically, homolgous expansion

\[v \propto r\]

This sort of velocity law has the advantage of being very simple to work with, and is generally a good approximation for supernovae.

In pracise, the outflow (wind) is assumed to extend from an inner radius \(r_{min}\) to an outer radius \(r_{rmax}\). The physical idea here is not necessarily that the wind stops at the maximum radius, but rather that it is sufficiently dilute that spectrum formation beyond this point becomes unimportant.

In our implementation, the specifics of the velocity law are determined by giving the outflow speed at \(r_{min}\) via a parameter \(v_{base}\). It then follows that the velocity at all other points in the wind is

\[v = v_{base} \frac{r}{r_{min}}\]

The density in the wind is determined by setting a mass flux at the inner boundary ( \(boundary\_mdot\), in solar masses per year). The variation of the density at larger radii is the controlled by an exponent (\(density\_exponent\)) such that

\[\rho \propto r^{- density\_exponent}\]

The Knigge, Wood, and Drew prescription of a bi-conical wind

In the KWD paramterization for a bi-conical flow the wind is envisioned to have poloidal streamlines that all point to a position a distance d below the disk, as is shown below:

Todo

This figure needs modification to account for the fact that we allow rmin and rmax to be specified.

images/kwd.png

As descried by KWD95, streamlines emerge thoughout the entire disk, with the innermost streamline just grazing the surface of the central object and the outermost streamline emerging from the outer radius of the disk. In the current version of Python, while this is the default choice, the wind region can be restricted to streamlines that arise from between \(r_{min}\) and \(r_{rmax}\). For fixed values of \(r_{min}\) and \(r_{rmax}\), the wind will tend to be more collimated the larger the value of d.

In the KWD parameritization, the mass loss per unit area per unit area of the disk is given by

\[\frac{\delta \dot{m}}{\delta A} \propto T(R)^{4\alpha}\]

where T(R) is the temperature of the disk at radius R. With this parameterization, the mass loss rate per unit area is constant for \(\alpha=0\) and is porportional to the luminoous flux is \(\alpha=1\).

The KWD95 model incorporates a velocity law reminiscent of a stellar wing, viz.

\[v_l=(nc_{s}) + (v_{\infty} - nc_{s})\left(1- \frac{R_{v}}{l+R_{v}} \right)^{\beta}\]

where \(nc_s\) is the speed at the base of flow,a multiple of the sound speed, \(R_v\) is the scale length, \(beta\) is the exponent that determines the number of scale lengths over which the wind acclerates, and \(v_{\infty}\) is defined as a multiple of the escape velocity at the footpoint of each stream line.

For the model, the sound speed of the disk is defined to be

\[c_s(R) = 10 \sqrt{\frac{T_{eff}(R)}{10^4 K}} km s^{-1}\]

Meta-documentation

How to document Python

This documentation is written in ReStructured Text, and parsed by Sphinx. A general guide to ReStructured Text can be found here. We’re trying to maintain a roughly consistent format for the documentation.

Installing the documentation tools

This guide is produced using Sphinx. Sphinx is written in python and available from the pip package manager. We require the Python 3 version of Sphinx. Install it, and the other modules required, as:

cd docs/sphinx
pip3 install -r requirements.txt
Building the documentation

Once Sphinx is installed, you can make the documentation using a Makefile as:

make html

You can tell if the documentation was built successfully by looking at the output of make html. You should see:

build succeeded.

The HTML pages are in html.

If the build was successful then the documentation can be viewed by opening docs/sphinx/html/index.html. Many errors will not stop the build process. Details on the build errors can be found in the section on Common errors & warnings.

You can make minor changes to the documentation and recompile using make html again. If you add new pages or move existing ones, the table of contents will need to be regenerated. Do this via:

make clean
make html

General documentation

Conventions

Each file should have a title, with subsections within it created by underlining the titles as:

Title
#####

Section
=======

Subsection
----------

Subsubsection
^^^^^^^^^^^^^

When referring to a parameter, link to the documentation page as:

The number of domains can be set by :ref:`Wind.number_of_components`.

When referring to files, code (e.g. shell script) or values used by the code, render it as monospaced text as:

Run the program using ``py``.
Set the parameter :ref:`Wind.type` to ``SV``.
Outputs can be found in ``filename.rst``.

When referring to a library or program name, render it in bold as:

Though this program is called **Python**, it is written in **C**, using the **GSL** library.

Content of interest to developers but not users should be broken out into a callout as:

.. admonition :: Developer Note

    This value is only stored to single-precision

Developer Note

This is a developer note

Documentation that needs expanding should be indicated with a to-do callout as:

.. todo :: Expand this section

Todo

This is a to-do note

Content relating to a specific GitHub issue/pull request can be linked directly to it as #1/#56:

This arose due to issue :issue:`1`, which was fixed by :user:`kslong` using :pr:`56`.

When writing a table, use the full form where possible as:

+----+----+
|Name|X   |
+----+----+
Name X

Parameter documentation

Formatting

Parameters are documented in a consistent way. They have a set of properties. Not every parameter will have all properties but you should fill them all in where possible. A full example outline is:

Title
=====
Description.
Use :ref:`Parameter.name` to link to other parameters, or other pages within the documentation.


Type
  Enumerator

Values

  option
    Description
    Multi-line if desired

  other
    More description

    Child(ren)
    * :ref:`Corona.radmin`

  yet_another
    More description

    Child(ren)
    * :ref:`KWD.rmin`
    * :ref:`KWD.rmax`

File
  `filename.c <https://github.com/agnwinds/python/blob/master/source/filename.c>`_

Parent(s)
  * :ref:`System_type`: `agn`, `binary`

The sections we expect are entered as a definition list. A definition list consists of titles followed by a definition block indented by 2 characters. The headings, in the order we expect, are:

Name
The parameter name, as used by Python input files.
Description

A description of the parameter and its function. This can include links to other pages and parameters, using the format

Use :ref:`Parameter.name` to link to other parameters, or other pages within the documentation.
Type
This is whether the parameter is an integer, float, or enumerator (a list of choices).
Unit
This is the unit. It can be something like cm, m or even derived from other parameters (e.g. Central_object.radius).
Values

If the parameter is an integer or float, this should describe the range of values it can take. For example, Greater than 0 or 0-1.

If the variable type is Enumerator, then instead it should include a nested definition list of the possible choices. Where each choice implies a different set of possible children (e.g. Wind.type) then each choice should have its own Children definition list, e.g.

SV
  * :ref:`SV.thetamin`
  * :ref:`SV.thetamax`
File
The file the parameter is found in. This is a link to the file on the master branch.
Child(ren)
If the parameter implies any others. For example, Spectrum.no_observers has child parameters Spectrum.angle.
Parent(s)
If the parameter depends on another. For example, KWD.rmax is only required for a specific choice of Wind.type.
Locations

Parameters are stored in `docs/sphinx/source/inputs/parameters/.

If multiple parameters share a root (i.e. SV.radmin, SV.radmax), then they should be stored within a directory with the same root name as the parameters (i.e. SV/SV.radmin.rst, SV/SV.radmax.rst). In the level above that directory, there should be a .rst file with the same name that serves to link those files to the table of contents, as:

SV
==

Some description of the parameter group.

.. toctree::
   :glob:

   SV/*

Storing all the parameters in one folder would result in it being unreadably busy. Instead, we sift the parameters into groups. Where multiple different parameters or parameter folders fall into the same rough category (e.g. central object parameters, wind types and the like) we create subfolders to group them into. The order that these appear in the sidebar can be set if you enter the filenames explicitly in the docs/sphinx/source/input/parameters.rst file.

Common errors & warnings

Undefined Label
/path/to/file.rst:line_number:
WARNING: undefined label: label_name (if the link has no caption the label must precede a section header)

This warning occurs when the :ref:'location' format is used to link to a section that does not exist. Check the spelling

Duplicate Label
/path/to/file.rst:line_number:
WARNING: duplicate label label_name, other instance in /path/to/other/file.rst

This warning occurs when two sections have the same name. The autosectionlabel addon automatically creates a label for each title and section component. This is generally not a problem unless you really need to

Inline emphasis
/path/to/file:line_number:
WARNING: Inline emphasis start-string without end-string.

This warning occurs when a line contains an un-escaped * character, as * is used to denote emphasis. Either escape it with \ (i.e. \*) or wrap it in a :code: tag.

Bullet list ends without a blank line
/path/to/file.rst:line_number:
WARNING: Bullet list ends without a blank line; unexpected unindent.

This warning occurs when a bullet-list doesn’t have a blank line between it and the next bit of text. It commonly happens when you accidentally forget to space a bullet and the text following it, e.g.

* blah1
* blah2
*blah3
Inline substitution_reference
/path/to/file:line_number:
WARNING: Inline substitution_reference start-string without end-string.

This warning occurs when you have a word immediately followed by a pipe that is not part of a table, e.g. something|. It tends to occur during typos in table creation e.g.

+---+---+
| a||b  |
+---+---+

Python Scripts

There are several Python (the scripting language) scripts written to prepare input for and analyse the output of python (the C code).

These are not yet documented in a way Sphinx can read for inclusion into this documentation. Currently, to generate the documentation, you must navigate to docs/pydocs and run write_docs.py.

The resulting output file can be found here.

Quick Guide to Python

This quide is intended to allow users to install Python, to run Python as a computer program and then to check whether the run has completed as expected.

It does not describe except in passing any information about the physics of Python, the details of a particular wind model,or criteria for evaluating whether the inputs correspond to a plausible model of an astrophysical system.

Radiative Transfer Modes

Python has a number of different radiative transfer modes, which affect the treatment of lines and scattering, and also whether we use indivisible packet constraints or allow photon weights to be attenuated by continuum absorption. These modes are selected with the parameter Line_transfer. The different modes are briefly described on that parameter page. This page is designed to give an overview of the assumptions and concepts behind, as well as the basic operation of, the different techniques. The aim is that, in partnership with the parameter page and the atomic data documentation, all the information regarding the different radiative transfer modes should be present.

For introductions and references regarding Monte Carlo radiative transfer techniques generally, we recommend reading Noebauer & Sim 2019. For specifics regarding Python, we recommend reading Long & Knigge 2002 as well as PhD theses by Higginbottom and Matthews.

Sobolev Approximation

Python always uses the Sobolev approximation to treat line transfer. In this approximation, it is assumed that the thermal line width is small compared to the velocity gradient. The Sobolev approximation is described extensively in astrophysics literature, and so we do not describe it in detail here. We refer the users to section 8.2 of Noebauer & Sim 2019 and references there in for a discussion of the Sobolev escape probabilities approach.

Weight Reduction v Indivisible Packets

Python was originally written in such a way that photon packet weights were not indivisible and allowed to be attenuated. This is the way the code is described in the original Long & Knigge 2002 paper. In the standard, weight reduction mode, photon weights are attenuated by continuum opacities (free-free, bound-free). Conservation of energy is then hopefully achieved by calculating the emission from the wind .

In indivisible packet mode, there is a fundamental shift in philosophy. All energy packets are strictly indivisible and conserve energy whenever they undergo radiative processes (the only exception is adiabatic cooling). Thus, even bound-free absorption is dealt with at a single interaction point.

Indivisible packet mode is activated by setting the Line_transfer parameter to either macro_atoms or macro_atoms_thermal_trapping. The terminology adopted here is slightly confusing, since the line transfer mode does not explicitly include a macro-atom treatment of atomic species (see next subsection).

Developer Note

The radiative transfer mode is stored using the code variable geo.rt_mode.

Macro-atoms and 2-level-atoms

The macro-atom scheme was devised by Leon Lucy in the early 2000s (Lucy 2002, Lucy 2003). It involves the reformulation of the process of radiation transport through a plasma in radiative equilibrium into a traffic-flow problem. Lucy showed that, when in radiative equilibrium, the energy flows through a system depend only on the transition probabilities and atomic physics associated with the levels the energy flow interacts with. By quantising this energy flow into radiant (r-) and kinetic (k-) packets, we can simulate the energy transport through a plasma discretised into volume elements (macro-atoms), whose associated transition probabilities govern the interaction of radiant and kinetic energy with the ionization and excitation energy associated with the ions of the plasma.

Todo

add refs, describe properly.

Developer Note

Macro-atoms are identified using their atomic data, in particular by providing data with identifiers LevMacro, LinMacro, PhotMacro.

Simple-atoms still interact with r- and k-packets, but do not possess internal transition probabilities. As a result, they are analogous to the two-level atom treatment, as any excitation is immediately followed by a deactivation into an r- or k-packet. The choice of radiative or kinetic deactivation is made according to the relative rates in the two-level atom formalism.

Isotropic v Anisotropic Line Scattering

Python always treats electron scattering as an isotropic process, and continuum emission processes are also treated as isotropic, except for Compton scattering. For Compton scattering, the direction and energy change is calculated self-consistently according to the energy change formula \(E/E'=1+(h \nu/mc^2)(1+\cos\theta)\). We first draw a random cross section that our photon packet will see. This cross section represents an energy change and hence a direction. The distribution of angles is taken care of by using a differential cross section vs energy change function.

Caution

Compton scattering is currently not accounted for when using indivisible packet mode.

Line emission and scattering is isotropic unless one of the thermal_trapping line transfer modes is selected. In the thermal trapping mode, any line interaction or emission results in an anisotropic direction being generated. This direction is generated by a rejection method which samples the Sobolev escape probability in each direction from the line interaction region. Unless you specifically want to consider isotropic line emission, we recommend always using the anisotropic thermal trapping mode.

Todo

move the below to where we describe photon sources and generation?

In the case of isotropic emission, the direction of a photon packet is chosen so that the probability of emission in each bin of solid angle is the same. It follows that

\[p(\Omega)d\Omega \propto \cos \theta \sin \theta d\theta d\phi\]

where the angles are in polar coordinates and relative to the local outward normal. For a spherical emitting source, such as a star, one must first generate a location on the star’s surface and then calculate the photon direction relative to the normal at the point. For emission from optically thick surfaces the above equation can be modified to include linear limb darkening, \(\eta(\theta)\), such that

\[p(\theta, \phi) d\theta d\phi = \eta(\theta) \cos \theta \sin \theta d\theta d\phi.\]

The Eddington approximation is usually adopted in the code, so that $eta(theta)$ is given by

\[\eta(\theta) = a (1 - \frac{3}{2} \cos \theta).\]

The constant \(a\) is normalised such that the total probability sums to 1. Whenever a radiation packet undergoes an electron scatter, the new direction is chosen to be isotropic. However, when the photon is a line photon, the new direction is chosen according to a line trapping model, which samples a probability distribution according to the Sobolev escape probability in different directions.

Doppler Shifts and The Comoving Frame

When calculating opacities, the photon frequency must be shifted from the rest frame of the photon into the rest frame of the plasma. This shift depends on the before and after directions of the photon. Let us denote these two directions with unit vectors \(\vec{n}_i\) and \(\vec{n}_f\), respectively, and consider a situation when a photon scatters off an electron in a region of the wind moving at velocity \(\vec{v}\). The final frequency of the photon with initial frequency is

\[\nu_f = \nu_i ~\frac{1 - (\vec{v} \cdot \vec{n}_i) / c}{1 - (\vec{v} \cdot \vec{n}_f) / c}.\]

In the case of a resonance scatter with line transition u to j, the new frequency is

\[\nu_f = \frac{\nu_{uj}}{1 - (\vec{v} \cdot \vec{n}_f) / c}.\]

The above formulae are the non-relativistic case, which is currently used in the code. However, this should in general be improved to use the special relativistic formula. This would produce more accurate Doppler shifts for the fastest regions of an outflow, as the current treatment introduces errors of order 5 Angstroms at the blue edges of the highest velocity absorption lines in quasar and CV wind models.

When real photons resonantly (or electron) scatter off real plasma in a flow, they conserve energy and frequency in the co-moving frame of the plasma. In the case of an outflow, doing the frame transformation from system->flow->system over the course of an interaction results in a redshifting of a photon, and as a result an energy loss - in other words, the photon does work on the flow even though energy is conserved in the co-moving frame. Indivisible packet schemes (such as macro-atoms) often enforce strict energy conservation in the frame of a given cell (physically, but see also Lucy 2002). This means that, when keeping track of packets in the observer frame, one needs to correct the energies (not just the frequencies) using a Doppler shift. Python does not currently conserve energy in the co-moving frame.

Todo

test whether this is an issue.

The Shlosman & Vitello prescription of a bi-conical wind

In the SV93 prescription, the wind emerges between \(r_{min}\) and \(r_{rmax}\) along streamlines whose orientation with respect to the system are described an angle

\[\theta = \theta_{min} + (\theta_{max} - \theta_{min}) x^{\gamma}\]

where

\[x=\frac{r_o - r_{min}}{r_{max}-r_{min}}\]

and \(r_o\) refers to the footpoint of a streamline.

_images/sv.png

The geometry of a Shlosman & Vitello wind

The poloidal velocity along the streamlines is defined to be

\[v_l - v_o + (v_{\infty}(r_o)-v_o) \frac {(l/R_v)^{\alpha}}{(l/R_v)^{\alpha}+1}\]

The scale length \(R_v\) and the exponent \(\alpha\) control the acceleration of the wind between a velocity \(v_o\), at the base of the wind and the terminal velocity \(v_{\infty}(r_o)\). The initial velocity \(v_o\) can be set to either a constant, normally 6 km/s, or a multiple of the sound-speed at the streamline base. The terminal velocity of each streamline varies depending on the location of the streamline in the inner and outer disk, being characterized as a fixed multiple of the escape velocity at the footpoint of the streamline. Thus the poloidal velocity is greatest for stream lines that originate from the inner regions of the disk, since the gravitational potential that must be overcome is greatest there.

The mass loss per unit surface area \(\delta \dot{m}/\delta A\) of the disk is controlled by a parameter \(\lambda\) such that

\[\frac{\delta\dot{m}}{\delta A} \propto \dot{m}_{wind} r_o^{\lambda} cos(\theta(r_o))\]

With this prescription, the overall mass loss rate declines with radius if \(\lambda\) is somewhat less than -2.

To use the SV93 prescription, therefore, one must provide the basic parameters of the system, the mass of the WD, the accretion rate, the inner and outer radius of the disk, and in addition, for the wind \(\dot{m}_{wind}\), \(r_{min}\), \(r_{max}\), \(\theta_{min}\), \(\theta_{max}\), \(\gamma\), \(R_{\nu}\), \(\alpha\), \(\lambda\), and the multiple of the escape velocity to be used for \(v_{\infty}\).