>











Nightfall User Manual







by R. Wichmann
(rwichman@lsw.uni-heidelberg.de)


Contents

1 Introduction

NIGHTFALL is an interactive application that introduces into the fascinating realm of eclipsing binary stars. Apart from their light variations that make them interesting objects for observations, eclipsing binaries are of fundamental importance for astrophysics, e.g. for measuring the mass of stars. NIGHTFALL is capable of producing:

It is, however, not able to fry your breakfast egg on your harddisk.

Eclipsing binary stars are most often very close systems. In such systems, owing to tidal forces, the shapes of both stars can be highly nonspherical, up to the possible formation of an 'overcontact' system, where both stars form a single, dumbbell-shaped object.

NIGHTFALL is a mildly ultramundane program of baroque complexity (I like Verdi and Händel on lazy sunday mornings - friday evenings are better with Iron Maiden and a good whisky). NIGHTFALL is based on a physical model that takes into account the nonspherical shape of stars in close binary systems, as well as mutual irradiance of both stars, and a number of additional physical effects.

NIGHTFALL can handle a large range of configurations, including overcontact systems, eccentric (non-circular) orbits, surface spots and asynchroneous rotations (stars rotating slower or faster than the orbital period), and the possible presence of a third star in the system.

NIGHTFALL supports OpenGL (or MesaGL) for an animated display of binary systems. Individual frames of the animation can be saved as JPEGs to create am MPEG movie.

NIGHTFALL can be build in a parallelized version, if your compiler supports the OpenMP standard (e.g. sufficiently recent versions of gcc do).

NIGHTFALL supports the GNOME desktop (if installed), but does not require it.
Also, NIGHTFALL supports internationalization. Currently, besides the default language (english), only german is supported. The language is selected by the environment variable LANG (must be set before starting the program, in sh, bash: LANG=de; export LANG
in csh, tcsh: setenv LANG de).

1.1 Remarks

NIGHTFALL is not part of my research work - rather it is the result of a recreational activity aimed at distracting myself from daily research work.

I have tested NIGHTFALL with published data on several binary stars. Data for these systems are included in the source code distribution. As no two light curve programs use exactly the same algorithms, results are never identical; however, results from NIGHTFALL appear to be within the range of similar light curve programs used in the literature. I hope that I have found most flaws in the logic of the code, at least in the 'scientific' part. If you want to use NIGHTFALL for a publication, you may want to evaluate its performance by yourself. NIGHTFALL comes WITHOUT ANY WARRANTY (see also the copyright statement for more information).

1.2 Bugs

Several, probably. If you find a bug and can eliminate it, send me a diff. If you find a bug and can't cope with it, send me a report, and wait for the next version. If you would like a feature, tell me.

1.3 Copyright

NIGHTFALL is copyright (c) 1998-2008 Rainer Wichmann, (c) 2001-2002 Markus Kuster, (c) 2001-2002 Patrick Risse

NIGHTFALL is free software; you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation; either version 2 of the License, or (at your option) any later version.

NIGHTFALL is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.

You should have received a copy of the GNU General Public License along with this program; if not, write to the Free Software Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.

This document is considered part of the NIGHTFALL program.

2 Requirements and Installation

2.1 Requirements

NIGHTFALL has been developed with the aim of being able to compile and run the program on a typical Linux system without any need of downloading/installing additional libraries or programs. Although I have not tested this, NIGHTFALL should compile on other Unix systems as well.

For command line use/file output, all you need is a C compiler (e.g. gcc). For interactive use and plotting, NIGHTFALL requires additional libraries/programs. Note: for compiling a program (e.g. for compiling NIGHTFALL), you need so-called 'development' packages of the required libraries. These 'development' packages are usually not installed by default.

Note: With GNUPLOT, plots are not that nice sometimes, but animated mode runs much smoother with GNUPLOT than with PGPLOT. GNUPLOT does not support incremental plotting and display of images, thus a few options are only available with PGPLOT.

2.2 Installation

After downloading the source code, do:

$ gunzip -c nightfall-(version number).tar.gz $\vert$ tar -xvf -
$ cd nightfall-(version number)
$ ./DoInstall.sh

The 'DoInstall.sh' script will query you for some information (e.g. where to install the program), and then (optionally) build, test, and install the application. If you want to do it by hand, instead of './DoInstall.sh' run the following sequence of commands:

$ ./configure
$ make
$ make install

If you want to use PGPLOT, 'configure' might find it by itself; otherwise you might need the 'configure' option
'-with-pgplot-include=PFX' (where PFX should be the directory where the PGPLOT header file cpgplot.h is installed), and '-with-pgplot-lib=PFX' (where PFX should be the directory where the PGPLOT library files libpgplot.a, libcpgplot.a are installed).

If you want to build a parallelized version, you need the 'configure' option '-enable-openmp'.

2.2.1 Install locations

By default, the binary is installed to /usr/local/bin, the data/help/configuration files to
/usr/local/share/nightfall. If you don't like this, use the 'configure' option '-prefix=/where/to/install', where '/where/to/install' will replace '/usr/local'.

2.3 Compile/configure problems

If you want to use a different compiler, you can supply the option CC=myCCompiler to configure (and also F77=myFortranCompiler, if you compile with PGPLOT).

2.4 Compile-time configuration

This section is targeted at researchers wanting to adjust NIGHTFALL to their needs. Several aspects/limits of NIGHTFALL can be configured at compile time. Mostly, these are related to fixed array sizes. For most users, the defaults should be sufficient.

To configure, edit the file 'Light.h' and change the default value(s) of the '#define ITEM' statement(s) in the first section of this file (where 'ITEM' stands for any of the bold words on the list below). Configuration options include (but are not limited to):

3 Usage

3.1 Simplified 'educational' mode

As of version 1.64, NIGHTFALL supports a simplified 'educational' mode that is targeted at undergraduate lab settings. In this mode, NIGHTFALL always starts in interactive mode, with a much simpler GUI stripped of almost all more or less advanced options. In additions, the sizes of the stars have to be specified in units of solar radii instead of fractions of the critical Roche volume (which allows to use NIGHTFALL without having to explain the concept of the Roche volume).

To start NIGHTFALL in this mode, use the command line parameter '-edu', like:

$ nightfall -edu

3.2 General

If used in interactive mode, only the command line parameter '-U' is required:

$ nightfall -U

If used in non-interactive mode, unless a configuration file is read in at startup (see info on configuration files), NIGHTFALL requires at least the following six numerical arguments on the command line (in that order):

These six numerical arguments are always required, if NIGHTFALL is used in command-line (i.e. non-interactive) mode without reading in a configuration file (see below)

$ nightfall -U -C ty_boo.cfg

will read parameters from a configuration file and start NIGHTFALL in interactive mode. The configuration file is a simple text file that can be edited by hand. In interactive mode, you can also write out the current parameters to a configuration file.

$ nightfall (without arguments) will produce a full list of options (many).

By default, NIGHTFALL will do nothing more than run in non-interactive mode, compute the lightcurve, write it to an output file 'NightfallCurve.dat', and exit silently. If you want more (nifty plots, etc.), read on.

3.3 Plotting, Graphics

requires that NIGHTFALL has been compiled with support for PGPLOT or Gnuplot.

In interactive mode, there are several menu options available to produce plots. In non-interactive (command-line) mode, you can choose from the following:

Detailed output is always written to a file 'NightfallCurve.dat'.

3.4 Interactive usage - the graphical user interface (GUI)

Use command line option -U to choose this option. A GUI will (hopefully) come up.

In the menu bar, there is a File menu for reading data/configuration files. Please note that you need to reset the data memory [Reset Memory] before reading in data for another binary, otherwise you will end up in a mess ...
The Output menu allows to choose between several facilities for viewing the output of a computation and visualizing the binary system. Some of them allow interactive change of parameters like viewing angle etc. Please read the online help for more details.

In the toolbar below, there are buttons for computing a lightcurve [Compute], toggling animated view of the binary [Animate] while computing, writing out the current binary configuration [Write], and getting help on the currently active notebook page [Context].

The rest of the GUI is layed out as a 'notebook' with several different pages for setting options:

On the bottom of the GUI, there is a status/progress bar.

3.5 OpenGL

Note that OpenGL can lead to program crashes with some graphic cards.

If compiled with OpenGL support, the animated display will by default run in OpenGL mode. You can switch to GNUPLOT/PGPLOT mode in the PlotOptions page in the GUI. In OpenGL mode, you have the following options:

With textures enabled, the default is to use a pixmap as texture. All other textures display the variation of physical quantities (caused by the non-spherical shape and/or limb darkening) across the stellar surface. Available choices are: surface temperature, surface gravity, flux towards the observer, and the velocity of surface elements towards the observer relative to the star's centre of mass.

3.6 Parallelization

Internally, if compiled with -enable-openmp, nightfall will run the computation of mutual reflection in parallel, using as many threads as there are cores on the host unless the environment variable OMP_NUM_THREADS is set.

In addition, it is possible to distribute the workload of an optimization run across several hosts. This requires that you can login via ssh to those hosts, without using a password (e.g. by using a RSA/DSA key and ssh-agent).

In order to make use of the latter option, you need to have a file parallel.cfg in your current working directory (the one from which you are starting NIGHTFALL, with one line per host. The format of valid lines is:
HOST hostname maximum_number_of_instances

I.e. the line starts with the keyword HOST, followed by the hostname, followed by a number which is the maximum number of NIGHTFALL instances to run on that host. E.g.:

HOST host1 5
HOST host2 5
HOST host3 2

Furthermore,

So basically, you need to make sure that /path/to/nightfall [options] will work if executed in your home directory when you ssh into the remote machine.

Note that for the simplex fitting routine, the maximum number of instances spawned is equal to the number of free parameters (this would happen when the simplex shrinks towards the best solution and all vertices except the best are recomputed). For the simulated annealing global optimization, the first sweep will spawn at most 100 instances, while all subsequent sweeps will spawn 5 instances.

3.7 Problems

3.8 Configuration File

NIGHTFALL comes with a set of configuration files, each of which will set the system parameters for a particular (real) binary system, and automatically load some data files with observed data for this particular binary star (see next section).

For the file format, see the commented example file 'ty_boo.cfg' in the source code distribution. To read in such a file, on the command line use:

-C path/to/config/file

If you don't give the full path, NIGHTFALL will search (in this order of priority) the present working directory, an eventual subdirectory './cfg', and the default data directory set at compile time. If the configuration file is in one of these, only the name, not the full path is required.

3.9 Data files

NIGHTFALL comes with a set of sample observational data for several different eclipsing binaries. For each included binary system, the data comprise lightcurves in several filters, as well as radial velocity curves. Also, for each of these systems there is a configuration file (see details on configuration files). Loading this configuration file will set the appropriate parameters for that particular system, and also read in the data. You can then compute the lightcurve and visualize the binary, with the observed data overlayed on the lightcurve.

E.g. you might call NIGHTFALL as 'nightfall -U -C data/ty_boo.cfg', then press 'ANIMATE' to switch on the real-time animation, and press 'COMPUTE' to compute the lightcurve. With 'PlotCurve' you can then get a plot of the final lightcurve, with the data overlayed.

Of course, you can also experiment by yourself and try to fit the lightcurve by varying the system parameters, or use the automatic fitting option.

To read in a single data file, use

-I path/to/data/file.

If you don't give the full path, NIGHTFALL will search (in this order) the present working directory, an eventual subdirectory './data', and the default data directory set at compile time. If the data file is in one of these, only the name, not the full path is required. For details on the format of the files, see here.

Data are available for the following systems:

3.10 Environment variables

The run-time behaviour of the program can be modified by the value of some environment variables. To set an environment variable to some value, depending on your shell the following command is required:

in csh, tcsh: setenv VARIABLE=value

in sh, bash: VARIABLE=value; export VARIABLE

4 Introduction to Binary Stars

This section (and the next) provide(s) an introduction to the problem (at a popular science level, hopefully), the options, and the algorithm(s) used. More technical aspects are in italics. If you find the less technical part too arcane, feel free to supply (constructive) suggestions.

4.1 The Roche Geometry

Imagine two lakes, seperated by a ridge. There are about three possible configurations:

Replace 'water' by 'gas', 'lake' by 'star', and you have the possible configurations of a close binary star system.

The stellar shapes in such a system are given by the sum of the gravitational forces, and the centrifugal force due to the orbital motion. Instead of using the forces, it is easier to use a potential (forces can then be expressed as the derivative of the potential, if needed). In the case of binary stars, this potential is called the Roche potential, named after the French mathematician Edouard Albert Roche (1820-1883). The stellar surface is then given by an equipotential surface (a surface, on which the potential is constant). Thus, the introduction of the potential makes the computation of any forces superfluous in this particular application.

The largest size a single star can have in a binary system is given by the Roche lobe - a teardrop-shaped equipotential surface, whose cusp touches the cusp of the Roche lobe of the other star at a point called Lagrange 1 (L1) (there are four more Lagrange points, which are of less interest here). L1 is located between both stars, on the line connecting their centres. At L1, the sum the forces is zero, thus, if a star fills its Roche lobe, at L1 matter can flow into the Roche lobe of the other star (provided it is not filled as well).

Thus, L1 would correspond to 'lowest point of the ridge' in the two-lake example above, and the Roche lobes would correspond to the two valleys below that point, that potentially can be filled by the two lakes. But just as the two lakes in the example above may be smaller than their maximum size (before flowing together), also the stars in a binary system may be smaller than their respective Roche lobes. A star that completely fills its Roche lobe will assume its teardrop-like shape. A star filling only a small fraction of its Roche lobe will be more spherical - distortion increases with the Roche lobe filling factor.

Note that the relative size of the Roche lobe of the two stars in the system depends on their mass ratio. The absolute size of the Roche lobe also scales with the separation of the two stars. Thus, with a fixed mass ratio and fixed absolute sizes of the stars (e.g. in kilometers), a star may fill its Roche lobe (and have a highly distorted figure), if the stars are rather close, but the same star might fill just a tiny fraction of its Roche lobe (and thus would be nearly spherical) if the binary separation would be rather large.

To define the sizes of stars in NIGHTFALL, you have to give the 'Roche lobe filling factor', which is defined in units of the Roche lobe (actually, its polar radius). NIGHTFALL uses a dimensionless potential, i.e. the distance between both stars is arbitrarily set to unity. Thus, for a fixed absolute size of the stars, decreasing the 'Roche lobe filling factor' would be equivalent to increasing the distance.

4.2 Shape of the lightcurve

The shape of the lightcurve depends mainly on three factors:

4.3 Suggested experiments

5 More details

5.1 Which star is which ?

Stars in binary systems are labelled 'Primary' and 'Secondary'. In NIGHTFALL, the star called 'Primary' is the star which is eclipsing first, i.e. the star that passes in front of the other one at orbital phase zero (orbital phase indicates the position of the stars in the orbit on a scale from zero to one). The secondary is the star that is eclipsed first. In animated view, at start the Primary is left, the Secondary right. Note that this labelling of 'Primary' and 'Secondary' is inverse to the usual convention ... mea culpa. Maybe I fix it sometime.

To exchange the stars (if needed), you can either (i) swap the eclipses in your data, or (ii) swap the stars themselves.
For (i), add half an orbital period to the phase zeropoint in your datafile (for a circular orbit), or the time lag of the second eclipse (for an eccentric orbit).
For (ii), swap temperatures and Roche lobe fill factors, and replace the mass ratio q by 1/q. Don't forget to swap spots, if you have. For an eccentric orbit, add or subtract 180 deg to/from the Periastron length.

For a circular orbit, NIGHTFALL starts at orbital phase -0.25 - which is identical to +0.75 -, and calculates up to +0.75. In a circular orbit, eclipses are at orbital phase 0.0 and 0.5.).

NIGHTFALL uses a dimensionless potential, i.e. the distance between both stars is arbitrarily set to unity. The relative size of both Roche lobes then depends only on the mass ratio of both stars. The size of both stars has to be given as a fraction of the Roche lobe (actually, its polar radius). This fraction might be larger than one, if you want to specify an overcontact system. To specify a semi-detached system, set the Roche lobe filling fraction to 1.0 for one star, less for the other.

5.2 How is eclipse testing done ?

NIGHTFALL divides each stellar surface into a grid of a few thousand elements. Eclipse testing is done by checking - for each surface element individually- whether the line of sight towards that surface element intersects the other star. For overcontact systems, a star might eclipse its own throat region (the region connecting both stars). This condition is tested as well.

Although simple tests based on orbital phase or intersection of sperical regions suffice in most cases, still sometimes a rigorous and expensive test is needed.

5.3 Temperature and Brightness

Stellar surface temperatures can range from a few 1000 Kelvin to several 10000 Kelvin. The respective brightness can be calculated from the so-called blackbody law (a blackbody is an idealized thermal radiation source). The blackbody law is applicable to thermal radiation, such as infrared radiation of your own body, or radiation by stars; it is not applicable to non-thermal radiation sources like lasers.

By default, NIGHTFALL uses the blackbody law, which is neither terribly good nor terribly bad. It is possible to use, instead, light fluxes from detailed numerical computations of stellar atmospheres. Contrary to blackbody fluxes, these model atmosphere fluxes (like real stellar fluxes) depend not only on temperature, but also (mildly) on surface gravity, and on the chemical composition of the stellar atmosphere. As actual atmosphere calculations would be prohibitively expensive, model atmosphere fluxes are taken from tables that cover only a limited range in temperature (3000 K to 35000 K). Tables are hardcoded, and only available for one surface gravity value (a compromise value that should be ok for most cases), and solar chemical composition.

To switch on model atmosphere fluxes, use

-M

5.4 Output Lightcurves

Lightcurves are output in eight commonly used broad-band filters (UBVRIJHK from near-UV to near-infrared), and four narrow-band filters (Stroemgren uvby). The shape of the lightcurve depends on the filter passband.

The human eye itself also is not equally sensitive to all wavelengths of light - it also is a kind of 'filter' for light. To create a lightcurve that looks like the human eye would see it, you have to choose an output filter that matches as close as possible the sensitivity of the human eye. For broadband filters, the V filter gives the best match to the human eye. In the Stroemgren filter system, you might want to choose Stroemgren v.

Lightcurves are output in magnitudes. This is a relative unit commonly used in astronomy, and defined as

\begin{displaymath}m_1 - m_2 = -2.5 \times \log{\frac{flux_1}{flux_2}}. \end{displaymath}

(For magnitude differences smaller than about 0.4, the difference in magnitude times 100 is nearly equal to the percentage difference in flux, i.e. 0.1 mag is about 10 per cent difference.) The brighter an object, the smaller its magnitude (Sun is -26.7 in the V filter, while the faintest stars visible by naked eye are about +6). NIGHTFALL uses the brightness at quadrature (phase -0.25 in a circular orbit, both stars fully visible) as normalization (i.e. as $flux_2$ in the equation above).

Output goes to a file 'NightfallCurve.dat'. To select plotting of output light curve, use

-G

To select the filter for the plotted lightcurve, use

-Bfilter,

where filter should be one of UBVRIJHKuvby. The default for plotting is the V filter, which is the best match to the human eye. If you read in data files with observed lightcurves (see here for more info ), the default will be the filter of the first lightcurve read in.

To obtain a hardcopy (i.e. PS file), use

-H

5.5 Gravity Brightening

The nonspherical shape of both stars causes a non-constant surface gravity. This, in term, causes brightness variations, with regions of higher surface gravity having higher brightness. This effect is called gravity brightening or gravity darkening (depending on which article you read). Gravity brightening can produce deep minima in the lightcurve, even if there is no eclipse at all !! NIGHTFALL always takes care of this effect - there is no option to switch it off.


6 Disk

In a binary system where one star fills its Roche lobe, there may be mass transfer via Roche lobe overflow. I.e. at the L1 point between the two stars, matter from the Roche-lobe filling star can flow into the Roche lobe of the other star. For physical reasons (conservation of angular momentum), this matter will not just fall onto the other star, but rather form a so-called accretion disk.

Within this disk, the gas orbits the star and loses angular momentum by friction, thus slowly moving towards the stellar surface. If the star has a strong magnetic field, this field may disrupt the disk close to the star. I.e. the disk will not reach all the way to the stellar surface, but rather there will be a gap between the inner edge and the star, and the gas will take the last part of its journey along the magnetic field lines.

The disk is characterized by inner and outer radius, temperature, height, and the disk model. Furthermore, optionally a hot spot may be included in the model.

The outer and inner edge of the disk are modelled as vertical ``walls''.

6.1 Disk model

For all models, the height of the disk at radius r is computed as

\begin{displaymath}a + b * r^c. \end{displaymath}

6.1.1 Simple disk

For a simple disk, $c = 1$, $a$ is set equal to half of the thickness parameter, and $b$ is set equal to the H/R parameter. Setting H/R to zero results in a flat disk. The temperature is constant and equal to the respective input parameter.

6.1.2 Isothermal disk

While the simple disk model is isothermal, it's run of the height with radius is not physically correct for an isothermal disk. In the isothermal disk, the exponent $c$ is set to $c = 1.5$, which is the correct value for an isothermal disk. Furthermore, the H/R parameter is ignored, the constant $a$ is set to $a = 0.0$, and the factor $b$ is chosen such that at the inner disk radius, the disk thickness equals the thickness input parameter.

6.1.3 Reprocessing disk

The reprocessing disk is the model for a disk that is reprocessing (absobing/re-radiating) the radiation from the central star. The exponent $c$ is set to $c = 9/8$, which is the correct value for a reprocessing disk. Furthermore, the H/R parameter is ignored, the constant $a$ is set to $a = 0.0$, and the factor $b$ is chosen such that at the inner disk radius, the disk thickness equals the thickness input parameter. Lastly, the temperature falls off with radius as $r^{-frac{3}{4}}$.

6.2 Disk parameters

6.3 Hot spot parameters

The disk may have a hot spot (where the accretion stream from the donor star hits the outer edge of the disk). The hot spot is assumed to have constant temperature and constant depth (as measured from the outer edge of the disk).

To switch off the hot spot, set the extent in longitude to zero.


7 Advanced Options

7.1 Fractional Visibility

As eclipse testing is only done for the centres of surface elements, the coarse surface grid introduces numerical artifacts that are readily visible sometimes. To fix this problem, it is possible to compute fractional visibilities for the surface elements (i.e. compute what fraction of a surface element is eclipsed). To switch on, use

-F

Note: this option is not available if NIGHTFALL is compiled in high precision mode (-enable-high-precision), since flux is not sufficiently accurate conserved.

7.2 Reflection

Due to the mutual irradiance of both stars, in addition to its own light, each star will also reflect light of its companion. This can be a very important effect, especially if the two stars are close (large Roche fill factor) or the temperature difference is large. NIGHTFALL offers two options:
(1) by default, the irradiating star is treated as a point source. This is ok for low Roche lobe filling factors, but not very accurate for stars filling a large fraction of the Roche lobe. It is rather fast, however.
(2) it is possible to compute the mutual irradiance of all pairs of surface elements, with up to nine iterations (usually, two to three are sufficient). This is an $N^2$ algorithm (N the number of surface elements), and thus computationally very expensive. To switch on, use

-Rn

where 'n' is the number of iterations.

7.3 Overcontact systems

For overcontact systems, the second Lagrange point (L2) comes into play. It is located behind the less massive star (as seen from the more massive), and has the same property as L1, i.e. the force vanishes there, and matter might flow out of the common envelope of an overcontact system (I suppose, by now you know where L3 is. L4 and L5 are in the orbital plane, left and right of the line connecting the centres of both stars. L3-5 are not terribly important for binary stars.)

Thus, the surface of the common envelope of an overcontact system has to be between the two equipotential surfaces given by the potentials of L1 and L2 (remember, the force is the derivative of the potential, thus if the net force is zero, the potential can still have a non-zero value - it just has to be 'flat' locally).

To have an overcontact system, set the Roche lobe filling factor of one star larger than 1.0 (as there is only one surface now, the smaller Roche lobe filling factor will be ignored). If you choose a too large value, NIGHTFALL will adjust it. (The largest possible value depends on the mass ratio. Anything larger than about 1.3 is rather unreasonable for any mass ratio.)

The combination of overcontact and non-circular orbit or asynchroneous rotation (see below) is not supported. It would be rather unphysical anyway, as the strong interaction (tidal forces and friction) would circularize and synchronize the system extremely rapidly.

7.4 Asynchroneous Rotation

Tidal forces in close binaries will tend to enforce synchroneous rotation (both stars rotating with the orbital period, thus showing each other always the same side) on a timescale usually shorter than stellar lifetimes.

There are, however, occasions when stars might rotate asynchroneously, e.g. young stars that are not yet rotating synchroneously, or massive stars that have short livetimes anyway, too short for synchronization to occur during their lifetime. Also, tidal forces fall off very rapidly with increasing distance, thus wide binaries are likely candidates for asynchroneous rotation (the tidal force is inverse proportional to the cube of the distance - unlike e.g. gravity, which is inverse proportional to the square of the distance and thus falls off much less rapidly).

A (nonstellar) example for synchronization is the Earth-Moon system, where Moon's rotation has been synchronized already, while Earth's rotation is known to slow down gradually. This would eventually lead to Earth's synchronization as well, but I suspect it might take longer than the lifetime of our solar system ...

To switch on asynchroneous rotation, use

-fP fratio or -fS fratio

[-fP for Primary, -fS for Secondary). 'fratio' is the ratio between stellar rotation period and orbital period.

Asynchroneous rotation modifies the Roche potential (see the info on 'Roche lobe'), and hence the equipotential surface which defines the stellar shape. One effect is a 'flattening' of the star for faster rotation. Another effect is that the 'critical lobe', i.e. the largest possible surface, is modified. For faster rotation, it becomes smaller than the Roche lobe (which is the 'critical lobe' for a synchroneously rotating star, as discussed here). In this situation, you might have two binary component, both filling their 'critical lobes', but well separated from each other.

For slower rotation, the 'critical lobe' can become larger than the Roche lobe... NIGHTFALL will check (and complain), if both stars intersect. (While an overcontact system with synchroneous rotation is no problem, asynchroneous rotation, with contact, will cause heavy friction, I would think, presumably leading to rapid synchronization. Thus a contact system with asynchroneous rotation probably is unrealistic.)

Surface spots (if there are any) will move with the stellar rotation rate (just as they should).

7.5 Eccentric Orbit

According to Keplers laws, the shape of the orbit is an ellipse. Often, it is close to a special case of an ellipse - a circle. However, sometimes binary orbits are markedly eccentric, i.e. non-circular. In an eccentric orbit, the distance and the orbital velocity is not constant. The stars will gain velocity, as they fall toward each other, and move fastest at Periastron (closest approach). They will lose velocity again as they move away from each other, and move slowest at Apastron (largest distance). Thus the time between first and second eclipse in general is different from the time between second and first. Also, the width of the eclipses may be different, as in general the stars will move with different velocities during the two eclipses. NIGHTFALL comes with sample data for the star 'V541 Cygni', an eccentric binary system where you can observe both effects.

Similar to asynchroneous rotation, the comment applies that tidal forces will act towards circularization of the orbit. Eccentric orbits are more likely in wide binaries than in close ones.

The problem is, that in an eccentric orbit the changing distance is equivalent to a changing Roche lobe filling factor, and thus a changing shape of the star (remember, the larger the Roche lobe filling factor, the larger also is the nonspherical distortion of a star's shape). While in a circular orbit the stellar surface is calculated only once, and then just rotated in space, in an eccentric orbit the stellar surface must be re-calculated at each step in orbital phase, thus causing substantial computational overhead.

To switch on this option, use

-e eccentricity periastron_length

where 'eccentricity' (range 0-1, 0 is circular orbit) is defined as

\begin{displaymath}e = \frac{r_2 - r_1}{r_2 + r_1}, \end{displaymath}

with $r_2$ is the largest and $r_1$ the smallest distance. Clearly, $e = 0.0$ if $r_2 = r_1$, which is the case for a circular orbit. $e = 1.0$ is a degenerate case (a parabola), which cannot be handled by the program. 'periastron_length' is the length (in degree) of the periastron, i.e. the point of closest approach in the orbit. To find out how it is counted, you may set the Roche lobe filling factor to 1.0, and e to a high value, like 0.5, Then, using the animation option (-A), you can identify the periastron easily, as the star will fill the Roche lobe at closest approach.

The input Roche lobe filling factor is assumed at Periastron. Due to the variable orbital velocity, 'synchroneous' rotation is not really synchroneous - rotation will lag behind the orbital motion for part of the orbit, and advance for the other part of the orbit. Surface spots (if there are any) will move accordingly.

7.6 Limb Darkening

The depth to which you can see into a star's atmosphere (where the visible light comes from) varies with the viewing angle. As the temperature (which determines the light flux) increases with depth, you can see hotter (= brighter) layers of the atmosphere towards the centre of the star's disk, where you can look deeper. Towards the limb, you see shallower, cooler, and thus less bright layers of the atmosphere. This results in the limb of a star being darker than the centre of its visible disk, an effect that can be seen readily on good (!) photographies of the Sun.

Limb darkening, as a function of the cosine of the viewing angle towards the stellar surface, is well aproximated by simple expressions. NIGHTFALL offers three different options, with expressions that are linear or include additional square or square root terms. The square root law is probably most accurate, but you will find that there is not much difference.

The default is the linear law. To change this, use

-Ln

with 'n' a number in the range 0-2 (0 = default).

To set the limb darkening coefficients explicitely (e.g. for fitting, if you use an external fitting routine), use

-l[P/S] lc_1_a,lc_1_b,..,lc_12_b

where 1..12 refers to the passbands in the order UBVRIJHK,Strömgren uvby.

7.7 Surface Spots

Cool stars like the sun (i.e. stars with surface temperatures of a few thousand Kelvin only) often have surface spots, which are regions of somewhat lower temperature on the surface. Among such stars, some are known to have surface spots ('starspots') much larger than those shown by the Sun. In extreme cases, spots may cover a few tenths of the stellar surface. Usually, these are cool spots (i.e. a few 100 K cooler than the surrounding area), caused by magnetic activity (like on the Sun).

To include spots, use

-sP longitude latitude radius dimfactor or

-sS longitude latitude radius dimfactor

(-sP for a spot on the Primary, -sS for a spot on the secondary). Spots are circular. The arguments are longitude and latitude of the spot's centre, the radius (all in degree) and the factor, by which the surface temperature is changed in the spot area. You can have multiple spots on each star. For overlapping spot areas, 'dimfactor' is averaged.

It is possible, but physically very unrealistic, to set 'dimfactor' to rather low or high values (0.5 - 2.0). Temperature deviations of more than about 1000 K may be not realistic. Hot spots are seen only in exceptional cases.

7.8 Radial Velocities

Unlike lightcurves, which can be expressed in a relative unit, radial velocity curves only make sense in absolute units (km/s in NIGHTFALL). Thus, they require absolute dimensions as input for the system. (NIGHTFALL supplies default values, however, if you don't want to bother about this.) Use

-tP period or

-tM mass or

-tD distance

('period' in days, 'mass' ( = total mass of both stars) in solar masses, 'distance' in solar radii.) You need to give two of these; the third can (and will) be calculated from Kepler's third law.

Radial velocities are computed as the sum of the orbital velocity of a point mass plus corrctions (flux-weighted contributions from each surface element). In animated view, you can see the resulting sum as well as the correction term (the latter multiplied by 2).

7.9 Line profiles

NIGHTFALL can calculate spectral line profiles at each phase step, which are output to files (one for each phase step). You can specify the rest wavelength 'lambda_zero' for the line. The luminosities (of individual surface elements) used to compute the profile are those for the passband with the closest effective wavelength. To switch on, use

-Plamda_zero

In interactive mode, you have the option to view the line profiles, and change the phase interactively.

Some numerical artifacts present, probably due to finite surface grid.

7.10 Third Light

The presence of an additional light source in the system (e.g. a third star) will decrease the contrast between eclipsed and non-eclipsed parts of the lightcurve. To include this effect, use

-3filter fraction

where 'filter' is one of the supported filters (UBVRIJHKuvby) and 'fraction' is the relative contribution of third light to the total system luminosity. I.e.:

\begin{displaymath}L1 + L2 + L3 = 1.0,\end{displaymath}

where $L1 + L2$ is the combined light flux from Primary and Secondary, and $L3$ is 'third light'.

$\underline{\smash{\hbox{This option is not tested yet.}}}$

7.11 Replacing the built-in filter set

It is possible to replace any or all of the built-in filters with (a) different one(s). The command line option to replace a filter is:

-r X:path

where X is one of the built-in filters (UBVRIJHKuvby), and path is the path to a file containing the data for the replacement filter. The format of the file is as follows (also see the example Filter_OM_UVM2.txt in the data/ subdirectory):

The filter definition file is a text file, containing certain keywords and data. Comments start with '#', empty lines are ignored. In flux/LDC tables, values are separated by whitespace. The keywords are:

NAME filter_name

WAVELENGTH filter_central_wavelength_in_mu

# (second moment of passbands)$^2$ / FWHM
# may be approximated following the prescription by
# Young A.T. 1992, Astronomy and Astrophysics 257, 366; Equation 12:
# m2 = FullWidth(0.05 of peak)/4.4
# $=>$ (m2)$^2$/FWHM
CORR_FACTOR numerical_value

# Limb Darkening Coefficients. Argument is the LD law the coefficients are for.
LDC_START linear OR sqrt OR quadratic OR four
# followed by a table of coefficients, sorted by $\log{g}$,
# and then $T_{eff}$ (ascending)
$T_{eff}$ $\log{g}$ $coeff_1$ $coeff_2$ ...

# Model flux tables (optional). Arguments are $\log{g}$, lower and upper
# limits for $T_{eff}$. Data sorted by $T_{eff}$ (ascending)
FLUX_START $\log{g}$ $T_{eff}$_lower $T_{eff}$_upper
$T_{eff}$ $flux(T_{eff})$

8 Fitting observed data

NIGHTFALL offers the possibility of determining a best-fit model for observed data. Several datasets can be fit simultaneously. Both a local and a global optimization algorithm are available. Fit results are written to a file 'NightfallFitted.dat'.

8.1 Important notes

8.2 Reading in the data

In order to determine a best-fit model to observed data, you first have to read them into memory. Use

-I path/to/data/file.

Only one file is read. To read more files, prepend each of them with a -I. You can read in (only) one datafile for each filter.

Each row in the file should consist of two or three numbers. The first is the date of the observation (as decimal number, in any unit you like), the second the measurement value (in magnitudes for brightness, in km/s for radial velocity), and the third (optional) the estimated error of the measurement.

Lines starting with a '#' are ignored, with the following exceptions (no blank after '#' !!):

#P period gives the orbital period of the system (same unit as dates). The program will use this value to fold the data into orbital phase. The default is 1.0 (thus assuming that your data are already folded in phase).

#C is a flag that tells the program to not fold the data into orbital phase. Instead, the program will synthesize a lightcurve covering the full timespan of the data.

#Z zeropoint gives the zeropoint for orbital phase, i.e. the time of Primary eclipse (same unit as dates). The default is 0.0 (again assuming that your data are already folded in phase).

#B filter gives the filter (UBVRIubvyJHK) in which your data have been observed. Default is V. For radial velocities, use '1' for Primary, '2' for Secondary.

#W error gives the average estimated error of measurements (if you do not have individual ones). Default is 0.01 (brightness) or 1.0 (radial velocities). You can mix individual and average error estimates (e.g. if you have individual error estimates only for some of your data).

#V system_velocity For radial velocity curves, use this parameter to set the system velocity. Radial velocities will be set to $data - system\_velocity$.

#N normalisation_phase Use this parameter to set the phase at which to normalize the light curve. The default is the starting point of the lightcurve.

#S shift This is a particular important parameter that requires some care. As noted in the info on brightness, NIGHTFALL will normalize its lightcurves to light at the normalisation_phase (see above), which by default is the starting point of the lightcurve (the leftmost in the plotted lightcurve). Thus, at this point, the brightness is zero magnitudes. For a meaningful fit, your observed lightcurve must be shifted up/down to have the value zero at this point as well (within the measurement errors). NIGHTFALL will try to do it automatically, but some correction may be required. Plot the lightcurve (data will be overplotted, if there are any), and check. Note: incorrect use of this parameter may make a fit look better, but the fitted parameters might be meaningless. Do not shift to anywhere else than the starting point of the lightcurve.

The source code distribution includes example files for the binaries TY Boo (an overcontact system) and some more. If in doubt, look into them.

8.3 Finding a local optimum

Determining a best-fit model means to optimize the parameters (like inclination, temperature, mass ratio, etc.) in such a way that the mismatch between model and observation is minimized. This mismatch is measured by a suitable merit function (sometimes also called cost function). NIGHTFALL uses the chi-square function as merit function.

General problems of determining best-fit parameters are:

(1) uniqueness: there may be several/many (nearly) equally good solutions. Do a few trials. Restart fitting with last best-fit as starting point.

(2) overfitting: the information content of a lightcurve is limited. Fitting too many parameters might produce good-looking, but meaningless results. Use as few free parameters as are sufficient for a good fit.

(3) local vs. global optimum: just that you find a valley, doesn't mean this is the deepest valley on Earth. The same applies to optimization problems. Any local optimization algorithm will only find a local optimum. If the problem is well-behaved, this will be the one and only, global, optimum. If the problem is badly-behaved, you might need a few trials with different starting points to find (hopefully) the global optimum. If the problem is even worse, you might need an awful lot of computing time to find the global optimum. Most examples in textbooks are nice. Most real-world optimization problems turn out to be bad or even ugly.

NIGHTFALL uses the so-called Simplex algorithm for local optimization. This is a direct search algorithm that is not terribly fast (and not terribly slow either), but very robust. To switch on optimization, use

-Xparameters tolerance

where 'parameters' is a string of characters indicating the parameters you want to fit (all others kept fixed), and 'tolerance' is the stopping criterion (something like 0.1 or less would be appropriate - 0.001 has a special meaning, see the section on global fitting). Use 'nightfall' without options to get the character codes for fit parameters.

If more than one data file is input (e.g. lightcurves in different filters), NIGHTFALL will fit all data simultaneously. This will probably work well only if different datasets are properly weighted - i.e. if the error estimates (or at least their ratios) are ok.

Have a cup of coffee ready. Use the -Db option (switch on 'Busy' in interactive mode) to see what is going on meanwhile.

IMPORTANT: Restart fit with last best-fit as starting point. Continue this until you are sure that the solution has converged and does not improve anymore, i.e. the value of 'SDV (Chi square)' does not significantly change anymore. Otherwise, your results may be completely meaningless. (swich on the -Db option to see 'SDV (Chi square)' for each iteration - see the info on debug options.

Output is always written to a file 'NightfallFitted.dat'.

8.4 Goodness-of-fit

To evaluate how good a fit is, NIGHTFALL offers the following options:

(1) if you plot the lightcurve, residuals will be plotted as well. Look at them to check whether there are systematic trends ( = bad fit).

(2) for a good fit, the residuals should scatter randomly around the model, with no systematic trends. This can be quantified by computing the 'runs statistic' (the number of runs = occurences of two or more consecutive residuals above or below the model curve). Obviously, a large number of runs would occur for a strictly alternating sequence, which is very unlikely for a random sequence. Likewise, only two runs would occur for the first half of the data below, the second above the model curve - also suspiciously non-random. NIGHTFALL will print out the actual number of runs, and the lower and upper limits for a 90 percent confidence interval.

(3) in theory, the goodness-of-fit can be evaluated from the $\chi^2$ (Chi Square) value of the fit (the function actually minimized in parameter fitting), which should be close to unity for a good fit. However, this only works if the error estimates for the measurements are realistic - neither to high nor to low. This is very rarely the case in astronomy ...

8.5 Finding a global optimum

To find a global optimum, basically a stochastic search strategy is required (or an exhaustive seach of the complete parameter space ...). There are different possibilities, ranging from complete random search to some 'intelligent' variation of random search. NIGHTFALL offers 'Simulated Annealing', which is a kind of mathematical implementation of the cooling of matter (leading to crystallization, i.e. an energy optimum, if cooling is slow enough). Switch on by setting the fit tolerance to 0.001 (in command-line mode).

Be prepared for a computing time on the order of a day or more (if you have no other CPU-expensive job running).

I am not sure whether the algorithm is correctly implemented. In the present implementation, cooling might be too fast, or might stop at a too high 'temperature'. However, my own experiments were rather satisfactory most of the time.

Apparently, 'Simulated Annealing' does not mathematically guarantee that the global optimum is indeed found, unless 'cooling' proceeds infinetely slow ...

8.6 Mapping the Chi-Square function

This option will create a two-dimensional map of the merit function (i.e. in the case of NIGHTFALL the Chi-Square function that measures the goodness of a fit) with respect to two parameters. Start values are the current values, step values can be entered. The gridsize is fixed at compile time (default 16 x 16, i.e. 256 lightcurves will be evaluated). To switch on this option, use

-Xparameters step1 step2

Parameters are coded like in the fitting option (see above), but only two parameters should be chosen.

9 Miscellaneous

9.1 Debug options

-Dcharacterstring

Most debug options (selected by 'characterstring') are of little use unless you know the code rather well. Some will produce excessive (and excessively messy) output.

Exceptions are:

9.2 Output of the surface map

It is possible to obtain the 2D surface map of the stars, as seen by an observer, i.e. the map displayed in animated mode. You have to set the environment variable NIGHTFALL_SMAP_PATH to the path of the map file (the program will append the index of the phase to this, i.e. you will get a seperate file for each step in orbital phase). If NIGHTFALL_SMAP_PATH is undefined, no map will be printed.

By default, the map will be output for the V band; you can change that with the environment variable NIGHTFALL_SMAP_BAND (0..11 for UBVRIJHKuvby).

The map includes for each visible surface element:
(1) the index of the star,
(2) the index of the surface element,
(3, 4) x, y coordinates in the viewing plane of the observer, with the star at (0, 0),
(5, 6) x, y coordinates as above, but with the centre of mass at (0, 0),
(7) $\cos{\gamma}$, the line-of-sight angle,
(8) temperature,
(9) dimensionless gravity,
(10) area,
(11) flux.
The flux is not normalized to the area; it is the flux that this surface element contributes to the total flux.

9.3 User-defined wavelenghts

It is possible to compute monocromatic fluxes at up to twelve different, user-defined wavelengths. These wavelengths (unit: micrometer) must be provided as a comma-seperated list in the environment variable NIGHTFALL_MONO_WAVE. They will replace the wavelenghts used in the blackbody approximation.

10 Technical details

This section is intended as a technical reference for experts who want to familiarize themself with the nasty details and the algorithms used.

10.1 Geometry

The geometric setup is based on a paper [11] by Djuraševic (1992a). In a cartesian coordinate system $(x, y, z)$, the stars are located at $(0, 0, 0)$ and $(1, 0, 0)$, and the $z$-axis is perpendicular to the orbital plane. A normalized, dimensionless Roche potential is used, which at a point $P(x,y,z)$ is given by the value $C$ as

\begin{displaymath}C = \frac{1}{r_1} + q (\frac{1}{r_2} - x) + \frac{q + 1}{2} ( x^2 + y^2 )f^2, \end{displaymath}

with $r_1 = x^2+y^2+z^2$, $r_2 = (x-1)^2+y^2+z^2$, mass ratio $q = m_2/m_1$, and nonsynchronism parameter $f = w/w_k$ (i.e. the ratio of the angular velocity to the Keplarian angular velocity).
For practical purposes, a spherical coordinate system $(r, \eta, \phi)$ is used, which is defined by

\begin{eqnarray*} x & = & r\ \cos{\eta},\\ y & = & r\ \sin{\eta}\ \cos{\phi},\\ z & = & r\ \sin{\eta}\ \sin{\phi}.\\ \end{eqnarray*}


The surface of the star is then divided into elementary areas by a grid in $(\eta, \phi)$, and for each surface element the gravity acceleration $g(r, \eta, \phi)$, the area $dS(r, \eta, \phi)$, and the normal vector $(l, m, n)$ is computed as outlined in [11].

Notes: (i) Equation (1-16) in [11] has a minor typo (wrong sign for $l$).
(ii) Djuraševic [11] apparently uses equidistant steps in $\phi$, thus leading to very unhomogeneous elementary areas. NIGHTFALL avoids this by adjusting the number of steps $N_{\phi}$ as

\begin{displaymath}N_{\phi} = 10 + N_{\eta} sin{\eta},\end{displaymath}

with $N_{\eta}$ a (compile-time) constant.

10.1.1 Accretion disk

Eclipse testing for the accretion disk is based on the method described in [11] by Djuraševic (1992b).

10.2 Reflection and gravity darkening/brightening

There are two options available for the treatment of reflection. Both are bolometric corrections, in the sense that the bolometric flux of the irradiating component is used to modify the temperature distribution on the irradiated component. See [26] for a discussion of this issue.

The 'simple' option for the treatment of reflection is described in [11]. The correction should be exact for spherical stars.

The 'detailed' reflection treatment loops over all pairs of surface elements $(dS_1,dS_2)$ and sums up, for each surface element $dS_1$, the irradiation by all visible surface elements $dS_2$ of the other star. Again, bolometric irradiation is computed. Of course, the true temperature of the irradiating surface elements (including reflection) is not known, thus it is necessary to iterate the algorithm. Convergence is typically reached after 2-3 iterations. The algorithm is described in [13].

For both treatments, by default for convective stars (below 7700 K) an albedo of 0.5, and for radiative stars an albedo of 1.0 is used. It is possible to choose different albedo values at program start using the command line option(s) -aP value for the primary and -aS value for the secondary.

Gravity darkening is computed by

\begin{displaymath}T(r, \eta, \phi) = T_{eff}(\frac{g(r, \eta, \phi)}{g_{eff}})^{\beta}. \end{displaymath}

For the gravity brightening exponent, the results from [6] (Fig. 1 in the paper) are used, which provide a smooth transition to the Von Zeipel (1924) exponent of 0.25 for radiative stars.

Notes: (i) The 'simple' reflection treatment in [11] has been supplemented by a penumbral correction for the partial visibily of the other star, if it is at the horizon. This penumbral correction assumes that the horizon is flat and that the other star is spherical, i.e. that the visible part is a segment of a circle.
(ii) For the 'detailed' reflection treatment, instead of the quadratic limb darkening law used by [13], a square root law ([10])

\begin{displaymath}I(\mu) = I(1)[1 - c(1-\mu) - d(1 - \sqrt{\mu})] \end{displaymath}

is used (where $\mu$ is the cosine of the angle subtended by the emergent radiation and the direction perpendicular to the stellar surface), and bolometric limb darkening coefficients are taken from [7]. The normalization factor for the square root law was calculated as

\begin{displaymath}\frac{1}{\pi\ (1 - c/3 - d/5)}.\end{displaymath}

(iii) The temperature limit of 7700K dividing stars with convective/radiative envelopes can be changed by the environment variable NIGHTFALL_RADIATIVE (within the range 4000 - 12000K). This changes the albedo only (0.5 below, 1.0 above).

10.3 Doppler boost

Since version 1.88, the computed flux includes a correction for the Doppler boost caused by the radial velocity of the components relative to the observer. This correction is computed following Eq. (3) in [25] van Kerkwijk et al. 2010, using the effective wavelength of the filter passbands for $\lambda$.

10.4 Spots

Spots are always circular, and characterized by four parameters: longitude $\lambda$, latitude $\beta$, radius $r$, and a 'dimming' factor $A_p = T_p/T$, i.e. the ratio of the (local) temperature $T_p$ with spot to the temperature $T$ without spot. A detailed discussion of the trigonometric expressions used to identify surface elements within the spot area can be found in [11], section I-3. In the 'detailed' reflection treatment, reflection is calculated with spots, i.e. the spots are applied first, then the reflection is calculated. For overlapping spots, in the overlap area the mean value of their $A_p$ is used.

10.5 Output flux

In the blackbody approximation, for each stellar component and each filter a (temperature-dependent) effective wavelength is computed, following Equation 3.30 in [4]:

\begin{displaymath}{\lambda}_{eff} = {\lambda}_{0} + \frac{ 5 {\mu}_2 ({\lambda}_{p}-{\lambda}_{0}) }{\sqrt{{\lambda}_{0}}}, \end{displaymath}

where ${\lambda}_{0}$ is the filter wavelength, and ${\lambda}_{p}$ the wavelength of the blackbody peak for the effective temperature of the respective component. The required second moments ${\mu}_2$ of the filter passbands are computed following the prescription by [27], Equation 12:

\begin{displaymath}{\mu}_2 \simeq \frac{FW(0.05)}{4.4}, \end{displaymath}

where $FW(0.05)$ is the full width at the 0.05 level.

For the 'model atmosphere' option, fluxes for temperatures below 9800K are from Hauschildt et al. (private comm., 2006) PHOENIX models, otherwise flux tables from Kuruzc models ([19]). are used. All models are for solar abundances. PHOENIX models below 2000K incorporate dust formation, with dust remaining in situ (no settling). The (originally monochromatic) fluxes have been integrated over the filter passbands (as given in [20]). Surface gravities $\log{g} = 3.5 - 5.0$ (in steps of 0.5) are available.

Three different limb darkening approximations are available, a linear, a quadratic, and a square root one, which are given by the following expressions, respectively:

\begin{eqnarray*} I(\mu) & = & I(1)[1 - u(1-\mu)],\\ I(\mu) & = & I(1)[1 - a(1-... ... and}\\ I(\mu) & = & I(1)[1 - c(1-\mu) - d(1 - \sqrt{\mu})],\\ \end{eqnarray*}


where $\mu$ is the cosine of the angle subtended by the emergent radiation and the direction perpendicular to the stellar surface. Limb darkening coefficients are from [7]. A discussion of the relative merits of these three approximations can be found in [10].

Notes: (i) When tabulated values are used, for temperatures out of range (below 2000K and above 31000K for $\log{g} = 3.5$, 39000K for $\log{g} = 4.0$, 49000K for $\log{g} = 4.5$, 50000K for $\log{g} = 5.0$) the lowest/highest value is used (i.e. no extrapolation is attempted). (ii) It is possible to set the limb darkening coefficients on the command line, with the option -l[P/S] lc_1_a,lc_1_b,..,lc_12_b, where 1..12 refers to the passbands in the order UBVRIJHK,Strömgren uvby. This can be used for fitting the coefficient when an external fit routine is used.

10.6 Eclipse testing

Eclipse testing follows the method proposed in [11]. First, the contact angle of the Roche lobes is used to exclude eclipses, if possible.
If an eclipse is possible, for each surface element $dS$ the line of sight (LOS) towards $dS$ is tested for intersection with the smallest sphere enclosing the eclipsing star (no intersection = not eclipsed). Then the LOS is tested for intersection with the largest sphere within the eclipsing star (intersection = eclipsed). This procedure takes advantage of the simplicity of testing the intersection of a line with a sphere.
If the LOS intersects the outer, but not the inner sphere, as a last (and computationally expensive) resort, the minimum of the Roche potential along the LOS is searched, and compared against the surface potential of the eclipsing star.

Notes: (i) In [11], the osculating cone (i.e. the tangent cone to the Roche surface at the inner Lagrangian point $L_1$) is used as the first eclipse criterion. According to [5], this can lead to serious errors for mass ratios very different from unity, as the larger Roche lobe is concave at the $L_1$. NIGHTFALL uses tabulated values from [5] (the cone angle $\phi_{\rm max}$ as given in their Table 2). (ii) In [11], the coordinate frame of the eclipsed star is used, for eclipse testing, which seems not to work for asynchroneous rotation. NIGHTFALL therefore uses the coordinate frame of the eclipsing star, following [2]. (iii) In [11], it is proposed to evaluate the potential along the LOS at six steps only. NIGHTFALL uses a more rigorous approach with a minimum finding routine (Brent's algorithm [3]).

10.7 Fractional visibility

The eclipse testing routine assigns to each surface element a visibility of 0 (eclipsed) or 1 (visible). However, the visibility is only evaluated for the centre of the element. This can lead to 'spikes' in the light curve, if large numbers of surface element centres become visible at once. Therefore, an option is provided to compute a fractional visibility for surface elements on the shadow limb.

The algorithm first searches all pairs of elements on the shadow limb (i.e. pairs with one element eclipsed, the adjacent element uneclipsed), and determines the potential minimum along the LOS towards them (which is typically already available from the eclipse test). Then, a linear approximation is made: for each pair, the distances $d_i$ to the shadow limb are taken as proportional to the differences between the minimum Roche potential $p_{i}^{LOS}$ along the LOS and the surface potential $p_{*}$ of the eclipsing star:

\begin{displaymath}f = \frac{d_1}{d_2} = \frac{p_{1}^{LOS} - p_{*}}{p_{2}^{LOS} - p_{*}}.\end{displaymath}

The surface element with the larger distance is assumed to be completely eclipsed or visible, and the other one of the pair is assigned a fractional visibility $(1/2 + \vert f\vert)$ based on the above approximation. This is not exact for individual surface elements, as the dividing line between them may not be parallel to the shadow limb. However, averaged over all surface elements on the limb, the error should be negligeable.

10.8 Eccentric orbits

In an eccentric orbit, the complete geometry and temperature distribution is re-calculated at each step ( = step in mean anomaly). As the Roche potential is used in a dimensionless form, the change in distance is equivalent to a change in the Roche volume filling factor ( = surface potential) with an unchanged unit distance.

This implies that the correct new surface potential must be found as a function of the stellar volume, which is scaled up/down with the distance. To avoid an iterative numerical integration of the stellar volume (i.e. varying the surface potential until the correct volume is found), which would be prohibitively expensive, the algorithm uses analytical approximations from [18] to derive the new surface potential.

The following procedure is used:

  1. Solve the Kepler equation to determine the distance of the stars (and their position in the orbit).
  2. Re-scale the distance to unity.
  3. Scale the stellar volume by the the cube of the distance scale factor.
  4. Find the 'mean radius' as the root of an analytical expression for the stellar volume (of 11th order in this 'mean radius').
  5. Compute the new surface potential as a function of the 'mean radius'.

Numerical details, including the definition of the 'mean radius', can be found in [18].

10.9 Optimization

For local optimization, the Simplex algorithm is used. This is a direct search algorithm that does not require derivatives. For $N$ free parameters, the simplex is a polyhedron with $(N + 1)$ vertices (or points) in the $N$-dimensional parameter space. At each step, the simplex moves through this parameter space according to some rules, basically moving away from its worst point. Details of the algorithm can be found in [17,22].

For global optimization, an implementation of the 'simulated annealing' method is provided. Basically, 'simulated annealing' does a stochastic search of the parameter space. Replacing the current best point with a better point (i.e. a downhill step) is always allowed, replacing it with a worse point (an uphill step) is allowed with some probability depending on the (steadily decreasing) 'temperature' of the system.
The implementation is based on the 'Very Fast Simulated Re-Annealing' algorithm [16], however, the 're-annealing' part is not included. For random number generation, the 'Mersenne Twister' random number generator (period length $2^{19937} - 1$) by Makoto Matsumoto and Takuji Nishimura [21] is used. It has a Mersenne prime period of $2^{19937} - 1$ (about $10^{6000}$) and is equi-distributed in 623 dimensions.

10.10 Runs Test

The nonparametric runs test is computed to provide an estimate of the goodness-of-fit that (unlike chi-square) is independent of the correctness of the measurement error estimates. In an ordered sequence of two different symbols (here: positive and negative fit residuals, ordered by phase), a run is a succession of one or more identical symbols, preceded and followed by the other symbol.

The expected number R of runs can be calculated analytically. If $m$ and $n$ are the number of negative and positive residuals, respectively, then the expectation value for the number of runs R, and its variance, are ([23]):


\begin{displaymath}\mathscr{E}_{\mathrm{(R)}} = 1 + \frac{2mn}{(m + n)} \end{displaymath}

and


\begin{displaymath}\mathrm{var(R)} = \frac{2mn(2mn - m - n)}{(m + n - 1)(m + n)^2} \end{displaymath}

Bibliography

1
Alencar S.H.P., Vaz L.P.R. (1997), A&A 326, 257

2
Antokhina E.A. (1996), ARep 40, 483

3
Brent R. (1973), Algorithms for minimization without derivatives, Prentice-Hall

4
Budding E. (1993), An Introduction to Astronomic Photometry, Cambridge University Press

5
Chanan G. A., Middleditch J., Nelson J.E. (1972), ApJ 208, 512

6
Claret A. (2000), A&A 359, 289

7
Claret A. (2000), A&A 363, 1081

8
Diaz-Cordovés J., Claret A., Giménez A. (1995), A&ASS 110, 329

9
Claret A., Diaz-Cordovés J., Giménez A. (1995), A&ASS 114, 247

10
Diaz-Cordovés J., Giménez A. (1992), A&A 259, 227

11
Djuraševic G. (1992a), Ap&SS 196, 241

12
Djuraševic G. (1992b), Ap&SS 196, 267

13
Hendry, P. D., Mochnacki, S. W. (1992), ApJ 388, 603

14
Hauschildt P.H., Allard F., Baron E. (1999), ApJ 512, 377

15
Hauschildt P.H., Allard F., Ferguson J., Baron E., Alexander D.R. (1999), ApJ 525, 871

16
Ingber L. (1989), J.Math.Comput.Modelling 12, 967

17
Kallrath J., Linnel A.P. (1987), ApJ 313, 346

18
Kopal Z. (1989), The Roche Problem, Kluwer Academic Publishers

19
Kurucz R. L. 1998, http://kurucz.harvard.edu/grids/gridp00/ip00k2.pck19, accessed 2006-04-12

20
Landolt-Boernstein (1982), Numerical Data and functional Relationships in Science and Technology. New Series, Berlin: Springer

21
Matsumoto M. & Nishimura T., ACM Transactions on Modeling and Computer Simulation, Vol. 8, No. 1, p. 3

22
Nelder J.A., Mead R. (1965), Comput. J. 7, 308

23
Rohatgi, Vijay K. (1984), Statistical Inference, John Wiley & Sons, Inc.

24
Van Hamme W. (1993), AJ 106, 2096

25
van Kerkwijk, M. H. et al. 2010, ApJ 715, 51

26
Wilson R.E. (1990), ApJ 356, 613

27
Young A.T. 1992, A&A 257, 366

11 Command line Options

   
'Educational':  
-edu switch on highly simplified 'educational' mode
Mandatory:  
(q) mass(Secondary)/mass(Primary)
(i) inclination angle (degree)
(rf1) Primary Roche fill factor
(rf2) Secondary Roche fill factor
(t1) Primary temperature
(t2) Secondary temperature
   
Interactive:  
-U Interactive mode
   
Graphic Output:  
-A Animated view
-V[v,i,c,a] Visualize geometry (default: v)
  v: view of stars
  i: image of potential
  c: contourmap of potential
  a: all of the above
-G[P,S,1,2] Graph of lightcurve (default: 1)
  P,S: close-up of Primary/Secondary eclipse
  1,2: display 1/2 orbital cycles
-B[U/B/V...] Bandpass to display in graph (default: V)
-H Hardcopy (postscript plot)
   
Files:  
-I datafile Read in a data file containing observed data
-C cfgfile Read in a configuration file
   
Advanced System Parameters:  
-f[P/S] F asynchroneous rotation ratio (Period/Period_Orbit)
-s[P/S] longitude latitude radius dimfactor Spot on Primary/Secondary
-e e w eccentric orbit, e = eccentricity, w = periastron length
-t[P/M/D] value period/total mass/separation in days/solar masses/solar radii
   
Debugging Options:  
-D[vwb] Debug [verbose,warn,busy]
   
Computation Options:  
-Plamda_zero Profile of absorption line at rest wavelength lamda_zero (nm)
-Nnn nn steps for lightcurve (default 80)
-M use Model atmosphere
-O[P/S] value $\log{g}$ (surface gravity) Primary/Secondary
-F compute Fractional visibility
-L[0-2] Limb darkening method (default: linear = 0)
  0: linear
  1: quadratic
  2: square root
-R[1-9] Reflection treatment
  0: Point source
  1-9: iterations for mutual reflection
-a[P/S] value Albedo of Primary/Secondary
  (default: 0.5 for T $<=$ 7700, 1.0 otherwise)
-3CM Third light, C: colour code, M: magnitude
   
Data Fitting:  
-X[..] Tolerance Fit the parameters coded in string [..]
  012345: q, i, rf1, rf2, t1, t2
  67: e, w (eccentric orbit)
  89: F(Primary),F(Secondary) (asynchroneous rotation)
  A-H: 2 Spots (Primary)
  I-P: 2 Spots (Secondary)
  QR: Mass, Separation
  a-l: Third Light (UBVRIJHKuvby)
  (Tolerance = 0.001 to use Simulated Annealing)
-Y[as above] Step1 Step2 Chi Square Map (2 Parameters)

About this document ...

This document was generated using the LaTeX2HTML translator Version 2008 (1.71)

Copyright © 1993, 1994, 1995, 1996, Nikos Drakos, Computer Based Learning Unit, University of Leeds.
Copyright © 1997, 1998, 1999, Ross Moore, Mathematics Department, Macquarie University, Sydney.

The command line arguments were:
latex2html -no_images -split 0 -show_section_numbers -no_navigation -link 2 UserManual.tex

The translation was initiated by rainer wichmann on 2018-06-22

rainer wichmann 2018-06-22