doc/user/escript.tex


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
% Copyright (c) 2003-2009 by University of Queensland
% Earth Systems Science Computational Center (ESSCC)
% http://www.uq.edu.au/esscc
%
% Primary Business: Queensland, Australia
% Licensed under the Open Software License version 3.0
% http://www.opensource.org/licenses/osl-3.0.php
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%


\chapter{The Module \escript}
\label{ESCRIPT CHAP}


\begin{figure}
\includegraphics[width=\textwidth]{figures/EscriptDiagram1}
\caption{\label{ESCRIPT DEP}Dependency of Function Spaces in Finley. An arrow indicates that a function in the 
function space at the starting point can be interpolated to the function space of the arrow target.
All functionspaces on the left side can be interpolated to any of the functionspaces on the right.}
\end{figure}

\escript is a Python module that allows you to represent the values of
a function at points in a \Domain in such a way that the function will
be useful for the Finite Element Method (FEM) simulation.  It also
provides what we call a function space that describes how the data is
used in the simulation.  Stored along with the data is information
about the elements and nodes which will be used by \finley.

In order to understand what we mean by the term 'function space',
consider that the solution of a partial differential equation
\index{partial differential equation} (PDE) is a function on a domain
$\Omega$.  When solving a PDE using FEM, the solution is
piecewise-differentiable but, in general, its gradient is
discontinuous.  To reflect these different degrees of smoothness,
different function spaces are used.  For instance, in FEM, the
displacement field is represented by its values at the nodes of the
mesh, and so is continuous.  The strain, which is the symmetric
part of the gradient of the displacement field, is stored on the
element centers, and so is considered to be discontinuous.

A function space is described by a \FunctionSpace object.  The
following statement generates the object \var{solution_space} which is
a \FunctionSpace object and provides access to the function space of
PDE solutions on the \Domain \var{mydomain}:

\begin{python}
  solution_space=Solution(mydomain)
\end{python}
The following generators for function spaces on a \Domain \var{mydomain} are commonly used: 
\begin{itemize}
\item \var{Solution(mydomain)}: solutions of a PDE.
\item \var{ReducedSolution(mydomain)}: solutions of a PDE with a reduced smoothness requirement.  
\item \var{ContinuousFunction(mydomain)}: continuous functions, eg. a temperature distribution.
\item \var{Function(mydomain)}: general functions which are not necessarily continuous, eg. a stress field.
\item \var{FunctionOnBoundary(mydomain)}: functions on the boundary of the domain, eg. a surface pressure.
\item \var{FunctionOnContact0(mydomain)}: functions on side $0$ of the discontinuity. 
\item \var{FunctionOnContact1(mydomain)}: functions on side $1$ of the discontinuity.
\end{itemize}
In some cases under-integration is used. For these cases the user may use a
\FunctionSpace from the floowing list:
\begin{itemize}
\item \var{ReducedFunction(mydomain)}
\item \var{ReducedFunctionOnBoundary(mydomain)}
\item \var{ReducedFunctionOnContact0(mydomain)} 
\item \var{ReducedFunctionOnContact1(mydomain)}
\end{itemize}
In comparison to the full teh corresponding full version
they use a reduced number of integration nodes (typically one only) to represent values.


The reduced smoothness for PDE solution is often used to fulfill the Ladyzhenskaya–-Babuska–-Brezzi condition \cite{LBB} when 
solving saddle point problems \index{saddle point problems}, eg. the Stokes equation.
A discontinuity \index{discontinuity} is a region within the domain across which functions may be discontinuous.  
The location of discontinuity is defined in the \Domain object.
\fig{ESCRIPT DEP} shows the dependency between the types of function spaces in Finley (other libraries may have different relationships).

The solution of a PDE is a continuous function. Any continuous function can be seen as a general function
on the domain and can be restricted to the boundary as well as to one side of
discontinuity (the result will be different depending on 
which side is chosen). Functions on any side of the  
discontinuity can be seen as a function on the corresponding other side. 

A function on the boundary or on one side of
the discontinuity cannot be seen as a general function on the domain as there are no values 
defined for the interior. For most PDE solver libraries
the space of the solution and continuous functions is identical, however in some cases, eg.
when periodic boundary conditions are used in \finley, a solution 
fulfills periodic boundary conditions while a continuous function does not have to be periodic.

The concept of function spaces describes the properties of 
functions and allows abstraction from the actual representation 
of the function in the context of a particular application. For instance, 
in the FEM context a
function of the \Function type (written as \emph{Function()} in Figure~\ref{ESCRIPT DEP})
is usually represented by its values at the element center, 
but in a finite difference scheme the edge midpoint of cells is preferred. 
By changing its function space you can use the same function in a Finite Difference
scheme instead of Finite Element scheme.
Changing the function space of a particular function
will typically lead to a change of its representation. 
So, when seen as a general function,
a continuous function which is typically represented by its values
on the node of the FEM mesh or finite difference grid 
must be interpolated to the element centers or the cell edges,
respectively. Interpolation happens automatically in \escript
whenever it is required.

In \escript the class that stores these functions is called \Data.
The function is represented through its values on \DataSamplePoints where
the \DataSamplePoints are chosen according to the function space 
of the function.  
\Data class objects are used to define the coefficients
of the PDEs to be solved by a PDE solver library 
and also to store the solutions of the PDE.

The values of the function have a rank which gives the
number of indices, and a \Shape defining the range of each index.
The rank in \escript is limited to the range $0$ through $4$ and
it is assumed that the rank and \Shape is the same for all \DataSamplePoints.
The \Shape of a \Data object is a tuple (list) \var{s} of integers. The length
of \var{s} is the rank of the \Data object and the \var{i}-th index ranges between $0$ and $\var{s[i]}-1$.
For instance, a stress field has rank $2$ and 
\Shape $(d,d)$ where $d$ is the spatial dimension.
The following statement creates the \Data object
\var{mydat} representing a 
continuous function with values 
of \Shape $(2,3)$ and rank $2$:
\begin{python}
  mydat=Data(value=1,what=ContinuousFunction(myDomain),shape=(2,3))
\end{python}
The initial value is the constant $1$ for all \DataSamplePoints and
all components.

\Data objects can also be created from any \numpy
array or any object, such as a list of floating point numbers, 
that can be converted into a \numpyNDA \cite{NUMPY}. 
The following two statements
create objects which are equivalent to \var{mydat}:
\begin{python}
  mydat1=Data(value=numpy.ones((2,3)),what=ContinuousFunction(myDomain))
  mydat2=Data(value=[[1,1],[1,1],[1,1]],what=ContinuousFunction(myDomain))
\end{python}
In the first case the initial value is \var{numpy.ones((2,3))}
which generates a $2 \times 3$ matrix as a \numpyNDA
filled with ones. The \Shape of the created \Data object
it taken from the \Shape of the array. In the second
case, the creator converts the initial value, which is a list of lists,
and converts it into a \numpyNDA before creating the actual
\Data object.      

For convenience \escript provides creators for the most common types
of \Data objects in the following forms (\var{d} defines the 
spatial dimension):
\begin{itemize}
\item \var{Scalar(0,Function(mydomain))} is the same as \var{Data(0,Function(myDomain),(,))} (each value is a scalar), 
e.g a temperature field. 
\item \var{Vector(0,Function(mydomain))} is the same as \var{Data(0,Function(myDomain),(d))} (each value is a vector), e.g
a velocity field.   
\item \var{Tensor(0,Function(mydomain))} is the same as \var{Data(0,Function(myDomain),(d,d))},
eg. a stress field.  
\item \var{Tensor4(0,Function(mydomain))} is the same as \var{Data(0,Function(myDomain),(d,d,d,d))}
eg. a Hook tensor field.   
\end{itemize}
Here the initial value is $0$ but any object that can be converted into a \numpyNDA and whose \Shape
is consistent with \Shape of the \Data object to be created can be used as the initial value.

\Data objects can be manipulated by applying unary operations (eg. cos, sin, log) point
and can be combined point-wise by applying arithmetic operations (eg. +, - ,* , /). 
It is to be emphasized that \escript itself does not handle any spatial dependencies as 
it does not know how values are interpreted by the processing PDE solver library. 
However \escript invokes interpolation if this is needed during data manipulations. 
Typically, this occurs in binary operation when both arguments belong to different
function spaces or when data are handed over to a PDE solver library 
which requires functions to be represented in a particular way. 

The following example shows the usage of {\tt Data} objects: Assume we have a
displacement field $u$ and we want to calculate the corresponding stress field
$\sigma$ using the linear--elastic isotropic material model
\begin{eqnarray}\label{eq: linear elastic stress}
\sigma\hackscore {ij}=\lambda u\hackscore {k,k} \delta\hackscore {ij} + \mu ( u\hackscore {i,j} + u\hackscore {j,i})
\end{eqnarray}
where $\delta\hackscore {ij}$ is the Kronecker symbol and 
$\lambda$ and $\mu$ are the Lame coefficients. The following function
takes the displacement {\tt u} and the Lame coefficients
\var{lam} and \var{mu} as arguments and returns the corresponding stress:
\begin{python}
  from esys.escript import *
  def getStress(u,lam,mu):
    d=u.getDomain().getDim()
    g=grad(u)
    stress=lam*trace(g)*kronecker(d)+mu*(g+transpose(g))
    return stress     
\end{python}
The variable 
\var{d} gives the spatial dimension of the 
domain on which the displacements are defined.
\var{kronecker} returns the Kronecker symbol with indexes 
$i$ and $j$ running from $0$ to \var{d}-1. The call \var{grad(u)} requires 
the displacement field \var{u} to be in the \var{Solution} or \ContinuousFunction
function space. The result \var{g} as well as the returned stress will be in the \Function function space. 
If, for example, \var{u} is the solution of a PDE then \var{getStress} might be called
in the following way:
\begin{python}
  s=getStress(u,1.,2.)
\end{python}
However \var{getStress} can also be called with \Data objects as values for
\var{lam} and \var{mu} which,
for instance in the case of a temperature dependency, are calculated by an expression. 
The following call is equivalent to the previous example:
\begin{python}
  lam=Scalar(1.,ContinuousFunction(mydomain))
  mu=Scalar(2.,Function(mydomain))
  s=getStress(u,lam,mu)
\end{python}

The function \var{lam} belongs to the \ContinuousFunction function space
but with \var{g} the function \var{trace(g)} is in the \Function function space.
In the evaluation of the product \var{lam*trace(g)} we have different function
spaces (on the nodes versus in the centers) and at first glance we have incompatible data.
\escript converts the arguments in an appropriate function space according to
Table~\ref{ESCRIPT DEP}. In this example that means
\escript sees \var{lam} as a function of the \Function function space. 
In the context of FEM this means the nodal values of 
\var{lam} are interpolated to the element centers.
The interpolation is automatic and requires no special handling.

\begin{figure}
\includegraphics[width=\textwidth]{figures/EscriptDiagram2}
\caption{\label{Figure: tag}Element Tagging. A rectangular mesh over a region with two rock types {\it white} and {\it gray}.
The number in each cell refers to the major rock type present in the cell ($1$ for {\it white} and $2$ for {\it gray}).
}
\end{figure}

Material parameters such as the Lame coefficients are typically dependent on rock types present in the 
area of interest. A common technique to handle these kinds of material parameters is "tagging"\index{tagging}, which
uses storage efficiently. Figure \fig{Figure: tag}
shows an example. In this case two rock types {\it white} and {\it gray} can be found in the domain. The domain
is subdivided into triangular shaped cells. Each 
cell has a tag indicating the rock type predominately found in this cell. Here $1$ is used to indicate
rock type {\it white} and $2$ for rock type {\it gray}. The tags are assigned at the time when the cells are generated
and stored in the \Domain class object. To allow easier usage of tags, names can be used instead of numbers. These names are typically defined 
at the time when the geometry is generated. 

The following statements show how, for the
example of \fig{Figure: tag}, the stress calculation discussed above and tagged values are used for
\var{lam}:
\begin{python}
  lam=Scalar(value=2.,what=Function(mydomain))
  insertTaggedValue(lam,white=30.,gray=5000.)
  s=getStress(u,lam,2.)
\end{python}
In this example \var{lam} is set to $30$ for those cells with tag {\it white} (=$1$) and to $5000.$ for those cells 
with tag {\it gray} (=$2$_. The initial value $2$ of \var{lam} is used as a default value for the case when a tag
is encountered which has not been linked with a value. The \var{getStress} method
does not need to be changed now that we are using tags.
\escript resolves the tags when \var{lam*trace(g)} is calculated.

This brings us to a very important point about \escript.
You can develop a simulation with constant Lame coefficients, and then later switch to tagged
Lame coefficients without otherwise changing your python script.
In short, you can use the same script to model with different domains and different types of input data.

There are three main ways in which \Data objects are represented internally: constant, tagged, and expanded.
In the constant case, the same value is used at each sample point and only a single value is stored to save memory. 
In the expanded case, each sample point has an individual value (such as for the solution of a PDE).
This is where your largest data sets will be created because the values are stored as a complete array.
The tagged case has already been discussed above.

Expanded data is created when you create a \Data object with expanded=True.
Tagged data sets are created when you use the insertTaggedValue() method as shown above.
 
Values are accessed through a sample reference number. Operations on expanded \Data
objects have to be performed for each sample point individually. When tagged values are used, the values are
held in a dictionary. Operations on tagged data require processing the set of tagged values only, rather than 
processing the value for each individual sample point. 
\escript allows any mixture of constant, tagged and expanded data in a single expression.

\Data objects can be written to disk files and read with \var{dump} and \var{load}, both of which use \netCDF\cite{NETCDF}.
Use these to save data for visualization, checkpoint/restart or simply to save and reuse data that was expensive to compute.

For instance to save the coordinates of the data points of the
\ContinuousFunction to the file {\tt x.nc} use 
\begin{python}
  x=ContinuousFunction(mydomain).getX()
  x.dump("x.nc")
  mydomain.dump(`dom.nc`)
\end{python}
To recover the object \var{x} and \var{mydomain} was a \finley mesh use
\begin{python}
  from esys.finley import LoadMesh
  mydomain=LoadMesh('dom.nc')
  x=load("x.nc", mydomain)
\end{python}
It possible to rerun the mechanism that was originally used to generates
\var{mydomain} to recreate \var{mydomain}. However in most cases using \var{dump} and
load is faster in particular if optimization has been applied. In case that
\escript is running on more than one \MPI processor the \var{dump} will create an individual file for each processor containing the local data. In order to avoid conflicts the file name is extended by the \MPI processor rank.

The function space of the \Data is stored in {\tt x.nc}, though.
If the \Data object
is expanded, the number of data points in the file and of the \Domain for the particular \FunctionSpace must match. 
Moreover, the ordering of the values is checked using the reference identifiers provided by 
\FunctionSpace on the \Domain. In some cases, data points will be re-ordered. Take care to be sure you get what you want!


\section{\escript Classes}
\declaremodule{extension}{esys.escript} 
\modulesynopsis{Data manipulation}

\subsection{\Domain class}
\begin{classdesc}{Domain}{}
A \Domain object is used to describe a geometric region together with 
a way of representing functions over this region.
The \Domain class provides an abstract interface to the domain of \FunctionSpace and \Data objects.
\Domain needs to be subclassed in order to provide a complete implementation.
\end{classdesc}
The following methods are available:
\begin{methoddesc}[Domain]{getDim}{}
returns the spatial dimension of the \Domain.
\end{methoddesc}
\begin{methoddesc}[Domain]{dump}{filename}
dumps the \Domain into the file \var{filename}.
\end{methoddesc}
\begin{methoddesc}[Domain]{getX}{}
returns the locations in the \Domain. The \FunctionSpace of the returned
\Data object is chosen by the \Domain implementation. Typically it will be
in the \Function.
\end{methoddesc}

\begin{methoddesc}[Domain]{setX}{newX}
assigns a new location to the \Domain. \var{newX} has to have \Shape $(d,)$
where $d$ is the spatial dimension of the domain. Typically \var{newX} must be
in the \ContinuousFunction but the space actually to be used depends on the \Domain implementation.
\end{methoddesc}

\begin{methoddesc}[Domain]{getNormal}{}
returns the surface normals on the boundary of the \Domain as \Data object.
\end{methoddesc}

\begin{methoddesc}[Domain]{getSize}{}
returns the local sample size, e.g. the element diameter, as \Data object.
\end{methoddesc}

\begin{methoddesc}[Domain]{setTagMap}{tag_name, tag}
defines a mapping of the tag name  \var{tag_name} to the \var{tag}. 
\end{methoddesc}
\begin{methoddesc}[Domain]{getTag}{tag_name}
returns the tag associated with the tag name \var{tag_name}.
\end{methoddesc}
\begin{methoddesc}[Domain]{isValidTagName}{tag_name}
return \True if \var{tag_name} is a valid tag name.
\end{methoddesc}

\begin{methoddesc}[Domain]{__eq__}{arg}
(python == operator) returns \True if the \Domain \var{arg} describes the same domain. Otherwise
\False is returned.
\end{methoddesc}

\begin{methoddesc}[Domain]{__ne__}{arg}
(python != operator) returns \True if the \Domain \var{arg} does not describe the same domain. 
Otherwise \False is returned.
\end{methoddesc}

\begin{methoddesc}[Domain]{__str__}{arg}
(python str() function) returns string representation of the \Domain.
\end{methoddesc}

\begin{methoddesc}[Domain]{onMasterProcessor)}{}
returns \True if the processor is the master processor within 
the \MPI processor group used by the \Domain. This is the processor with rank 0.
If \MPI support is not enabled the return value is always \True.
\end{methoddesc}

\begin{methoddesc}[Domain]{getMPISize}{}
 returns the number of \MPI processors used for this \Domain. If \MPI support is not enabled
1 is returned.
\end{methoddesc}

\begin{methoddesc}[Domain]{getMPIRank}{}
 returns the rank of the processor executing the statement 
within the  \MPI processor group used by the \Domain. 
If \MPI support is not enabled 0 is returned.
\end{methoddesc}

\begin{methoddesc}[Domain]{MPIBarrier}{}
executes barrier synchronization within 
the \MPI processor group used by the \Domain.
If \MPI support is not enabled, this command does nothing. 
\end{methoddesc}

\subsection{\FunctionSpace class}
\begin{classdesc}{FunctionSpace}{}
\FunctionSpace objects are used to define properties of \Data objects, such as continuity. \FunctionSpace objects
are instantiated by generator functions. A \Data object in a particular \FunctionSpace is 
represented by its values at \DataSamplePoints which are defined by the type and the \Domain of the
\FunctionSpace.
\end{classdesc}
The following methods are available:
\begin{methoddesc}[FunctionSpace]{getDim}{}
returns the spatial dimension of the \Domain of the \FunctionSpace.
\end{methoddesc}


\begin{methoddesc}[FunctionSpace]{getX}{}
returns the location of the \DataSamplePoints.
\end{methoddesc}

\begin{methoddesc}[FunctionSpace]{getNormal}{}
If the domain of functions in the \FunctionSpace 
is a hyper-manifold (e.g. the boundary of a domain)
the method returns the outer normal at each of the 
\DataSamplePoints. Otherwise an exception is raised.
\end{methoddesc}

\begin{methoddesc}[FunctionSpace]{getSize}{}
returns a \Data objects measuring the spacing of the \DataSamplePoints.  
The size may be zero.
\end{methoddesc}

\begin{methoddesc}[FunctionSpace]{getDomain}{}
returns the \Domain of the \FunctionSpace.
\end{methoddesc}

\begin{methoddesc}[FunctionSpace]{setTags}{new_tag, mask}
assigns a new tag \var{new_tag} to all data sample 
where \var{mask} is positive for a least one data point. 
\var{mask} must be defined on the this \FunctionSpace.
Use the \var{setTagMap} to assign a tag name to \var{new_tag}.
\end{methoddesc}

\begin{methoddesc}[FunctionSpace]{__eq__}{arg}
(python == operator) returns \True if the \Domain \var{arg} describes the same domain. Otherwise
\False is returned.
\end{methoddesc}

\begin{methoddesc}[FunctionSpace]{__ne__}{arg}
(python != operator) returns \True if the \Domain \var{arg} do not describe the same domain. 
Otherwise \False is returned.
\end{methoddesc}

\begin{methoddesc}[Domain]{__str__}{g}
(python str() function) returns string representation of the \Domain.
\end{methoddesc}

The following function provide generators for \FunctionSpace objects:
\begin{funcdesc}{Function}{domain}
returns the \Function on the \Domain \var{domain}. \Data objects in this type of \Function
are defined over the whole geometric region defined by \var{domain}. 
\end{funcdesc}

\begin{funcdesc}{ContinuousFunction}{domain}
returns the \ContinuousFunction on the \Domain domain. \Data objects in this type of \Function
are defined over the whole geometric region defined by \var{domain} and assumed to represent
a continuous function.
\end{funcdesc}

\begin{funcdesc}{FunctionOnBoundary}{domain}
returns the \ContinuousFunction on the \Domain domain. \Data objects in this type of \Function
are defined on the boundary of the geometric region defined by \var{domain}. 
\end{funcdesc}

\begin{funcdesc}{FunctionOnContactZero}{domain}
returns the \FunctionOnContactZero the \Domain domain. \Data objects in this type of \Function
are defined on side 0 of a discontinuity  within the geometric region defined by \var{domain}.
The discontinuity is defined when \var{domain} is instantiated.
\end{funcdesc}

\begin{funcdesc}{FunctionOnContactOne}{domain}
returns the \FunctionOnContactOne on the \Domain domain. 
\Data objects in this type of \Function
are defined on side 1 of a discontinuity  within the geometric region defined by \var{domain}.
The discontinuity is defined when \var{domain} is instantiated.
\end{funcdesc}

\begin{funcdesc}{Solution}{domain}
returns the \SolutionFS on the \Domain domain. \Data objects in this type of \Function
are defined on geometric region defined by \var{domain} and are solutions of
partial differential equations \index{partial differential equation}. 
\end{funcdesc}

\begin{funcdesc}{ReducedSolution}{domain}
returns the \ReducedSolutionFS on the \Domain domain. \Data objects in this type of \Function
are defined on geometric region defined by \var{domain} and are solutions of
partial differential equations \index{partial differential equation} with a reduced smoothness 
for the solution approximation.
\end{funcdesc}

\subsection{\Data Class}
\label{SEC ESCRIPT DATA}

The following table shows arithmetic operations that can be performed point-wise on
\Data objects.
\begin{tableii}{l|l}{textrm}{expression}{Description}
\lineii{+\var{arg0}} {identical to \var{arg} \index{+}}
\lineii{-\var{arg0}} {negation\index{-}}
\lineii{\var{arg0}+\var{arg1}} {adds \var{arg0} and \var{arg1} \index{+}}
\lineii{\var{arg0}*\var{arg1}} {multiplies \var{arg0} and \var{arg1} \index{*}}
\lineii{\var{arg0}-\var{arg1}} {difference \var{arg1} from\var{arg1} \index{-}}
\lineii{\var{arg0}/\var{arg1}} {divide \var{arg0} by \var{arg1} \index{/}}
\lineii{\var{arg0}**\var{arg1}} {raises \var{arg0} to the power of \var{arg1} \index{**}}
\end{tableii}
At least one of the arguments \var{arg0} or \var{arg1} must be a
\Data object.
Either of the arguments may be a \Data object, a python number or a \numpy object.

If \var{arg0} or \var{arg1} are
not defined on the same \FunctionSpace, then an attempt is made to convert \var{arg0}
to the \FunctionSpace of \var{arg1} or to convert \var{arg1} to
the \FunctionSpace of \var{arg0}. Both arguments must have the same
\Shape or one of the arguments may be of rank 0 (a constant).

The returned \Data object has the same \Shape and is defined on
the \DataSamplePoints as \var{arg0} or \var{arg1}.

The following table shows the update operations that can be applied to
\Data objects:
\begin{tableii}{l|l}{textrm}{expression}{Description}
\lineii{\var{arg0}+=\var{arg2}} {adds \var{arg0} to \var{arg2} \index{+}}
\lineii{\var{arg0}*=\var{arg2}} {multiplies \var{arg0} with \var{arg2} \index{*}}
\lineii{\var{arg0}-=\var{arg2}} {subtracts \var{arg2} from\var{arg2} \index{-}}
\lineii{\var{arg0}/=\var{arg2}} {divides \var{arg0} by \var{arg2} \index{/}}
\lineii{\var{arg0}**=\var{arg2}} {raises \var{arg0} by \var{arg2} \index{**}}
\end{tableii}
\var{arg0} must be a \Data object. \var{arg1} must be a
\Data object or an object that can be converted into a
\Data object. \var{arg1} must have the same \Shape as
\var{arg0} or have rank 0.  In the latter case it is
assumed that the values of \var{arg1} are constant for all
components. \var{arg1} must be defined in the same \FunctionSpace as
\var{arg0} or it must be possible to interpolate \var{arg1} onto the
\FunctionSpace of \var{arg0}.

The \Data class supports taking slices from a \Data object as well as assigning new values to a slice of an existing
\Data object. \index{slicing}
The following expressions for taking and setting slices are valid:
\begin{tableiii}{l|ll}{textrm}{rank of \var{arg}}{slicing expression}{\Shape of returned and assigned object}
\lineiii{0}{ no slicing }                      {-}
\lineiii{1}{\var{arg[l0:u0]}}                   {(\var{u0}-\var{l0},)}
\lineiii{2}{\var{arg[l0:u0,l1:u1]}}             {(\var{u0}-\var{l0},\var{u1}-\var{l1})}
\lineiii{3}{\var{arg[l0:u0,l1:u1,l2:u2]} }      {(\var{u0}-\var{l0},\var{u1}-\var{l1},\var{u2}-\var{l2})}
\lineiii{4}{\var{arg[l0:u0,l1:u1,l2:u2,l3:u3]}} {(\var{u0}-\var{l0},\var{u1}-\var{l1},\var{u2}-\var{l2},\var{u3}-\var{l3})}
\end{tableiii}
where \var{s} is the \Shape of \var{arg} and 
\[0 \le \var{l0} \le \var{u0} \le \var{s[0]},\]
\[0 \le \var{l1} \le \var{u1} \le \var{s[1]},\] 
\[0 \le \var{l2} \le \var{u2} \le \var{s[2]},\] 
\[0 \le \var{l3} \le \var{u3} \le \var{s[3]}.\]
Any of the lower indexes \var{l0}, \var{l1}, \var{l2} and \var{l3} may not be present in which case 
$0$ is assumed. 
Any of the upper indexes \var{u0}, \var{u1}, \var{u2} and \var{u3} may be omitted, in which case, the upper limit for that dimension is assumed. 
The lower and upper index may be identical, in which case the column and the lower or upper
index may be dropped. In the returned or in the object assigned to a slice, the corresponding component is dropped,
i.e. the rank is reduced by one in comparison to \var{arg}.
The following examples show slicing in action:
\begin{python}
  t=Data(1.,(4,4,6,6),Function(mydomain))
  t[1,1,1,0]=9.
  s=t[:2,:,2:6,5] # s has rank 3
  s[:,:,1]=1.
  t[:2,:2,5,5]=s[2:4,1,:2]
\end{python}

\subsection{Generation of \Data objects}
\begin{classdesc}{Data}{value=0,shape=(,),what=FunctionSpace(),expand=\False}
creates a \Data object with \Shape \var{shape} in the \FunctionSpace \var{what}.
The values at all \DataSamplePoints are set to the double value \var{value}. If \var{expanded} is \True
the \Data object is represented in expanded from.
\end{classdesc}

\begin{classdesc}{Data}{value,what=FunctionSpace(),expand=\False}
creates a \Data object in the \FunctionSpace \var{what}. 
The value for each \DataSamplePoints is set to \var{value}, which could be a \numpy, \Data object \var{value} or a dictionary of 
\numpy or floating point numbers. In the latter case the keys must be integers and are used
as tags.
The \Shape of the returned object is equal to the \Shape of \var{value}. If \var{expanded} is \True
the \Data object is represented in expanded form.
\end{classdesc}

\begin{classdesc}{Data}{}
creates an \EmptyData object. The \EmptyData object is used to indicate that an argument is not present
where a \Data object is required.
\end{classdesc}

\begin{funcdesc}{Scalar}{value=0.,what=FunctionSpace(),expand=\False}
returns a \Data object of rank 0 (a constant) in the \FunctionSpace \var{what}.
Values are initialized with \var{value}, a double precision quantity. If \var{expanded} is \True
the \Data object is represented in expanded from.
\end{funcdesc}

\begin{funcdesc}{Vector}{value=0.,what=FunctionSpace(),expand=\False}
returns a \Data object of \Shape \var{(d,)} in the \FunctionSpace \var{what},
where \var{d} is the spatial dimension of the \Domain of \var{what}.
Values are initialed with \var{value}, a double precision quantity. If \var{expanded} is \True
the \Data object is represented in expanded from.
\end{funcdesc}

\begin{funcdesc}{Tensor}{value=0.,what=FunctionSpace(),expand=\False}
returns a \Data object of \Shape \var{(d,d)} in the \FunctionSpace \var{what},
where \var{d} is the spatial dimension of the \Domain of \var{what}.
Values are initialed with \var{value}, a double precision quantity. If \var{expanded} is \True
the \Data object is represented in expanded from.
\end{funcdesc}

\begin{funcdesc}{Tensor3}{value=0.,what=FunctionSpace(),expand=\False}
returns a \Data object of \Shape \var{(d,d,d)} in the \FunctionSpace \var{what},
where \var{d} is the spatial dimension of the \Domain of \var{what}.
Values are initialed with \var{value}, a double precision quantity. If \var{expanded} is \True
the \Data object is re\var{arg}presented in expanded from.
\end{funcdesc}

\begin{funcdesc}{Tensor4}{value=0.,what=FunctionSpace(),expand=\False}
returns a \Data object of \Shape \var{(d,d,d,d)} in the \FunctionSpace \var{what},
where \var{d} is the spatial dimension of the \Domain of \var{what}.
Values are initialized with \var{value}, a double precision quantity. If \var{expanded} is \True
the \Data object is represented in expanded from.
\end{funcdesc}

\begin{funcdesc}{load}{filename,domain}
recovers a \Data object on \Domain \var{domain} from the file \var{filename}, which was created by \function{dump}.
\end{funcdesc}

\subsection{\Data methods}
These are the most frequently-used methods of the 
\Data class. A complete list of methods can be found on \ReferenceGuide.
\begin{methoddesc}[Data]{getFunctionSpace}{}
returns the \FunctionSpace of the object.
\end{methoddesc}

\begin{methoddesc}[Data]{getDomain}{}
returns the \Domain  of the object.
\end{methoddesc}

\begin{methoddesc}[Data]{getShape}{}
returns the \Shape  of the object as a \class{tuple} of
integers.
\end{methoddesc}

\begin{methoddesc}[Data]{getRank}{}
returns the rank of the data on each data point. \index{rank}
\end{methoddesc}

\begin{methoddesc}[Data]{isEmpty}{}
returns \True id the \Data object is the \EmptyData object.
Otherwise \False is returned.
Note that this is not the same as asking if the object contains no \DataSamplePoints.
\end{methoddesc}

\begin{methoddesc}[Data]{setTaggedValue}{tag_name,value}
assigns the \var{value} to all \DataSamplePoints which have the tag
assigned to \var{tag_name}. \var{value} must be an object of class
\class{numpy.ndarray} or must be convertible into a
\class{numpy.ndarray} object. \var{value} (or the corresponding
\class{numpy.ndarray} object) must be of rank $0$ or must have the
same rank like the object.
If a value has already be defined for tag \var{tag_name} within the object
it is overwritten by the new \var{value}.  If the object is expanded,
the value assigned to \DataSamplePoints with tag \var{tag_name} is replaced by
\var{value}. If no tag is assigned tag name \var{tag_name}, no value is set.
\end{methoddesc}

\begin{methoddesc}[Data]{dump}{filename}
dumps the \Data object to the file \var{filename}. The file stores the
function space but not the \Domain. It is in the responsibility of the user to
save the \Domain. 
\end{methoddesc}

\begin{methoddesc}[Data]{__str__}{}
returns a string representation of the object.
\end{methoddesc}

\subsection{Functions of \Data objects}
This section lists the most important functions for \Data class objects \var{a}.
A complete list and a more detailed description of the functionality can be found on \ReferenceGuide.
\begin{funcdesc}{saveVTK}{filename,**kwdata}
writes \Data defined by keywords in the file with \var{filename} using the 
vtk file format \VTK file format. The key word is used as an identifier. The statement
\begin{python}
  saveVTK("out.xml",temperature=T,velocity=v)
\end{python} 
will write the scalar \var{T} as \var{temperature} and the vector \var{v} as  \var{velocity} into the 
file \file{out.xml}. Restrictions on the allowed combinations of \FunctionSpace apply.
\end{funcdesc}
\begin{funcdesc}{saveDX}{filename,**kwdata}
writes \Data defined by keywords in the file with \var{filename} using the 
vtk file format \OpenDX file format. The key word is used as an identifier. The statement
\begin{python}
  saveDX("out.dx",temperature=T,velocity=v)
\end{python} 
will write the scalar \var{T} as \var{temperature} and the vector \var{v} as  \var{velocity} into the 
file \file{out.dx}. Restrictions on the allowed combinations of \FunctionSpace apply.
\end{funcdesc}
\begin{funcdesc}{kronecker}{d}
returns a \RankTwo \Data object in \FunctionSpace \var{d} such that
\begin{equation}
\code{kronecker(d)}\left[ i,j\right] = \left\{ 
\begin{array}{cc}
1 & \mbox{ if } i=j \\
0 & \mbox{ otherwise }
\end{array}
\right.
\end{equation}
If \var{d} is an integer a $(d,d)$ \numpy array is returned.
\end{funcdesc}
\begin{funcdesc}{identityTensor}{d}
is a synonym for \code{kronecker} (see above).
% returns a \RankTwo \Data object in \FunctionSpace \var{d} such that
% \begin{equation}
% \code{identityTensor(d)}\left[ i,j\right] = \left\{ 
% \begin{array}{cc}
% 1 & \mbox{ if } i=j \\
% 0 & \mbox{ otherwise }
% \end{array}
% \right.
% \end{equation}
% If \var{d} is an integer a $(d,d)$ \numpy array is returned.
\end{funcdesc}
\begin{funcdesc}{identityTensor4}{d}
returns a \RankFour \Data object in \FunctionSpace \var{d} such that
\begin{equation}
\code{identityTensor(d)}\left[ i,j,k,l\right] = \left\{ 
\begin{array}{cc}
1 & \mbox{ if } i=k \mbox{ and } j=l\\
0 & \mbox{ otherwise }
\end{array}
\right.
\end{equation}
If \var{d} is an integer a $(d,d,d,d)$ \numpy array is returned.
\end{funcdesc}
\begin{funcdesc}{unitVector}{i,d}
returns a \RankOne \Data object in \FunctionSpace \var{d} such that
\begin{equation}
\code{identityTensor(d)}\left[ j \right] = \left\{ 
\begin{array}{cc}
1 & \mbox{ if } j=i\\
0 & \mbox{ otherwise }
\end{array}
\right.
\end{equation}
If \var{d} is an integer a $(d,)$ \numpy array is returned.

\end{funcdesc}

\begin{funcdesc}{Lsup}{a}
returns the $L^{sup}$ norm of \var{arg}. This is the maximum of the absolute values
 over all components and all \DataSamplePoints of \var{a}. 
\end{funcdesc}

\begin{funcdesc}{sup}{a}
returns the maximum value over all components and all \DataSamplePoints of \var{a}.
\end{funcdesc}

\begin{funcdesc}{inf}{a}
returns the minimum value over all components and all \DataSamplePoints of \var{a}
\end{funcdesc}


\begin{funcdesc}{minval}{a}
returns at each \DataSamplePoints the minimum value over all components.
\end{funcdesc}

\begin{funcdesc}{maxval}{a}
returns at each \DataSamplePoints the maximum value over all components.
\end{funcdesc}

\begin{funcdesc}{length}{a}
returns at Euclidean norm at each \DataSamplePoints. For a \RankFour \var{a} this is
\begin{equation}
\code{length(a)}=\sqrt{\sum\hackscore{ijkl} \var{a} \left[i,j,k,l\right]^2}
\end{equation} 
\end{funcdesc}
\begin{funcdesc}{trace}{a\optional{,axis_offset=0}}
returns the trace of \var{a}. This is the sum over components \var{axis_offset} and \var{axis_offset+1} with the same index. For instance in the
case of a \RankTwo function and this is 
\begin{equation}
\code{trace(a)}=\sum\hackscore{i} \var{a} \left[i,i\right]
\end{equation} 
and for a \RankFour function and  \code{axis_offset=1} this is
\begin{equation}
\code{trace(a,1)}\left[i,j\right]=\sum\hackscore{k} \var{a} \left[i,k,k,j\right]
\end{equation} 
\end{funcdesc}

\begin{funcdesc}{transpose}{a\optional{, axis_offset=None}}
returns the transpose of \var{a}. This swaps the first \var{axis_offset} components of \var{a} with the rest. If \var{axis_offset} is not
present \code{int(r/2)} is used where \var{r} is the rank of \var{a}. 
 the sum over components \var{axis_offset} and \var{axis_offset+1} with the same index. For instance in the
case of a \RankTwo function and this is 
\begin{equation}
\code{transpose(a)}\left[i,j\right]=\var{a} \left[j,i\right]
\end{equation} 
and for a \RankFour function and  \code{axis_offset=1} this is
\begin{equation}
\code{transpose(a,1)}\left[i,j,k,l\right]=\var{a} \left[j,k,l,i\right]
\end{equation} 
\end{funcdesc}

\begin{funcdesc}{swap_axes}{a\optional{, axis0=0 \optional{, axis1=1 }}}
returns \var{a} but with swapped components \var{axis0} and  \var{axis1}. The argument \var{a} must be
at least of \RankTwo. For instance in the 
for a \RankFour argument, \code{axis0=1} and \code{axis1=2} this is
\begin{equation}
\code{swap_axes(a,1,2)}\left[i,j,k,l\right]=\var{a} \left[i,k,j,l\right]
\end{equation} 
\end{funcdesc}

\begin{funcdesc}{symmetric}{a}
returns the symmetric part of \var{a}. This is \code{(a+transpose(a))/2}.
\end{funcdesc}
\begin{funcdesc}{nonsymmetric}{a}
returns the non--symmetric part of \var{a}. This is \code{(a-transpose(a))/2}.
\end{funcdesc}
\begin{funcdesc}{inverse}{a}
return the inverse of \var{a}. This is 
\begin{equation}
\code{matrix_mult(inverse(a),a)=kronecker(d)}
\end{equation} 
if \var{a} has shape \code{(d,d)}. The current implementation is restricted to arguments of shape 
\code{(2,2)} and \code{(3,3)}.
\end{funcdesc}
\begin{funcdesc}{eigenvalues}{a}
return the eigenvalues of \var{a}. This is 
\begin{equation}
\code{matrix_mult(a,V)=e[i]*V}
\end{equation} 
where \code{e=eigenvalues(a)} and \var{V} is suitable non--zero vector \var{V}. 
The eigenvalues are ordered in increasing size.
The argument \var{a} has to be the symmetric, ie. \code{a=symmetric(a)}.  
The current implementation is restricted to arguments of shape 
\code{(2,2)} and \code{(3,3)}.
\end{funcdesc}
\begin{funcdesc}{eigenvalues_and_eigenvectors}{a}
return the eigenvalues and eigenvectors of \var{a}. This is 
\begin{equation}
\code{matrix_mult(a,V[:,i])=e[i]*V[:,i]}
\end{equation} 
where \code{e,V=eigenvalues_and_eigenvectors(a)}. The eigenvectors \var{V} are orthogonal and normalized, ie.
\begin{equation}
\code{matrix_mult(transpose(V),V)=kronecker(d)}
\end{equation} 
if \var{a} has shape \code{(d,d)}. The eigenvalues are ordered in increasing size.
The argument \var{a} has to be the symmetric, ie. \code{a=symmetric(a)}.  
The current implementation is restricted to arguments of shape 
\code{(2,2)} and \code{(3,3)}.
\end{funcdesc}
\begin{funcdesc}{maximum}{*a}
returns the maximum value over all arguments at all \DataSamplePoints and for each component.
For instance 
\begin{equation}
\code{maximum(a0,a1)}\left[i,j\right]=max(\var{a0} \left[i,j\right],\var{a1} \left[i,j\right])
\end{equation}
at all \DataSamplePoints.
\end{funcdesc}
\begin{funcdesc}{minimum}{*a}
returns the minimum value over all arguments at all \DataSamplePoints and for each component.
For instance 
\begin{equation}
\code{minimum(a0,a1)}\left[i,j\right]=min(\var{a0} \left[i,j\right],\var{a1} \left[i,j\right])
\end{equation}
at all \DataSamplePoints.
\end{funcdesc}

\begin{funcdesc}{clip}{a\optional{, minval=0.}\optional{, maxval=1.}}
cuts back \var{a} into the range between \var{minval} and \var{maxval}. A value in the returned object equals 
\var{minval} if the corresponding value of \var{a} is less than \var{minval}, equals \var{maxval} if the 
 corresponding value of \var{a} is greater than \var{maxval}
or corresponding value of \var{a} otherwise.
\end{funcdesc}
\begin{funcdesc}{inner}{a0,a1}
returns the inner product of \var{a0} and \var{a1}. For instance in the
case of \RankTwo arguments and this is 
\begin{equation}
\code{inner(a)}=\sum\hackscore{ij}\var{a0} \left[j,i\right]  \cdot \var{a1} \left[j,i\right]
\end{equation} 
and for a \RankFour arguments this is
\begin{equation}
\code{inner(a)}=\sum\hackscore{ijkl}\var{a0} \left[i,j,k,l\right]  \cdot \var{a1} \left[j,i,k,l\right]
\end{equation} 
\end{funcdesc}

\begin{funcdesc}{matrix_mult}{a0,a1}
returns the matrix product of \var{a0} and \var{a1}. If \var{a1} is \RankOne this is
\begin{equation}
\code{matrix_mult(a)}\left[i\right]=\sum\hackscore{k}\var{a0}  \cdot \left[i,k\right]\var{a1} \left[k\right]
\end{equation} 
and if \var{a1} is \RankTwo this is
\begin{equation}
\code{matrix_mult(a)}\left[i,j\right]=\sum\hackscore{k}\var{a0}  \cdot \left[i,k\right]\var{a1} \left[k,j\right]
\end{equation} 
\end{funcdesc}

\begin{funcdesc}{transposed_matrix_mult}{a0,a1}
returns the matrix product of the transposed of \var{a0} and \var{a1}. The function is equivalent to 
\code{matrix_mult(transpose(a0),a1)}.
If \var{a1} is \RankOne this is
\begin{equation}
\code{transposed_matrix_mult(a)}\left[i\right]=\sum\hackscore{k}\var{a0}  \cdot \left[k,i\right]\var{a1} \left[k\right]
\end{equation} 
and if \var{a1} is \RankTwo this is
\begin{equation}
\code{transposed_matrix_mult(a)}\left[i,j\right]=\sum\hackscore{k}\var{a0}  \cdot \left[k,i\right]\var{a1} \left[k,j\right]
\end{equation} 
\end{funcdesc}

\begin{funcdesc}{matrix_transposed_mult}{a0,a1}
returns the matrix product of \var{a0} and the transposed of \var{a1}.
The function is equivalent to
\code{matrix_mult(a0,transpose(a1))}.  
If \var{a1} is \RankTwo this is
\begin{equation}
\code{matrix_transposed_mult(a)}\left[i,j\right]=\sum\hackscore{k}\var{a0}  \cdot \left[i,k\right]\var{a1} \left[j,k\right]
\end{equation} 
\end{funcdesc}

\begin{funcdesc}{outer}{a0,a1}
returns the outer product of \var{a0} and \var{a1}. For instance if \var{a0} and \var{a1} both are \RankOne then
\begin{equation}
\code{outer(a)}\left[i,j\right]=\var{a0} \left[i\right]  \cdot  \var{a1}\left[j\right]
\end{equation} 
and if \var{a0} is \RankOne and \var{a1} is \RankThree
\begin{equation}
\code{outer(a)}\left[i,j,k\right]=\var{a0} \left[i\right] \cdot \var{a1}\left[j,k\right]
\end{equation} 
\end{funcdesc}

\begin{funcdesc}{tensor_mult}{a0,a1}
returns the tensor product of \var{a0} and \var{a1}. If \var{a1} is \RankTwo this is
\begin{equation}
\code{tensor_mult(a)}\left[i,j\right]=\sum\hackscore{kl}\var{a0}\left[i,j,k,l\right] \cdot \var{a1} \left[k,l\right]
\end{equation} 
and if \var{a1} is \RankFour this is
\begin{equation}
\code{tensor_mult(a)}\left[i,j,k,l\right]=\sum\hackscore{mn}\var{a0} \left[i,j,m,n\right] \cdot \var{a1} \left[m,n,k,l\right]
\end{equation} 
\end{funcdesc}

\begin{funcdesc}{transposed_tensor_mult}{a0,a1}
returns the tensor product of the transposed of \var{a0} and \var{a1}. The function is equivalent to
\code{tensor_mult(transpose(a0),a1)}.
If \var{a1} is \RankTwo this is
\begin{equation}
\code{transposed_tensor_mult(a)}\left[i,j\right]=\sum\hackscore{kl}\var{a0}\left[k,l,i,j\right] \cdot \var{a1} \left[k,l\right]
\end{equation} 
and if \var{a1} is \RankFour this is
\begin{equation}
\code{transposed_tensor_mult(a)}\left[i,j,k,l\right]=\sum\hackscore{mn}\var{a0} \left[m,n,i,j\right] \cdot \var{a1} \left[m,n,k,l\right]
\end{equation} 
\end{funcdesc}

\begin{funcdesc}{tensor_transposed_mult}{a0,a1}
returns the tensor product of \var{a0} and the transposed of \var{a1}. 
The function is equivalent to
\code{tensor_mult(a0,transpose(a1))}.
If \var{a1} is \RankTwo this is
\begin{equation}
\code{tensor_transposed_mult(a)}\left[i,j\right]=\sum\hackscore{kl}\var{a0}\left[i,j,k,l\right] \cdot \var{a1} \left[l,k\right]
\end{equation} 
and if \var{a1} is \RankFour this is
\begin{equation}
\code{tensor_transposed_mult(a)}\left[i,j,k,l\right]=\sum\hackscore{mn}\var{a0} \left[i,j,m,n\right] \cdot \var{a1} \left[k,l,m,n\right]
\end{equation} 
\end{funcdesc}

\begin{funcdesc}{grad}{a\optional{, where=None}}
returns the gradient of \var{a}. If \var{where} is present the gradient will be calculated in \FunctionSpace \var{where} otherwise a 
default \FunctionSpace is used. In case that \var{a} has \RankTwo one has
\begin{equation}
\code{grad(a)}\left[i,j,k\right]=\frac{\partial \var{a} \left[i,j\right]}{\partial x\hackscore{k}}
\end{equation} 
\end{funcdesc}
\begin{funcdesc}{integrate}{a\optional{ ,where=None}}
returns the integral of \var{a} where the domain of integration is defined by the \FunctionSpace of \var{a}. If \var{where} is 
present the argument is interpolated into \FunctionSpace \var{where} before integration. For instance in the case of 
a \RankTwo argument in \ContinuousFunction it is
\begin{equation}
\code{integrate(a)}\left[i,j\right]=\int\hackscore{\Omega}\var{a} \left[i,j\right] \; d\Omega
\end{equation} 
where $\Omega$ is the spatial domain and $d\Omega$ volume integration. To integrate over the boundary of the domain one uses  
\begin{equation}
\code{integrate(a,where=FunctionOnBoundary(a.getDomain))}\left[i,j\right]=\int\hackscore{\partial \Omega} a\left[i,j\right] \; ds
\end{equation} 
where $\partial \Omega$ is the surface of the spatial domain and $ds$ area or line integration.
\end{funcdesc}
\begin{funcdesc}{interpolate}{a,where}
interpolates argument \var{a} into the \FunctionSpace \var{where}. 
\end{funcdesc}
\begin{funcdesc}{div}{a\optional{ ,where=None}}
returns the divergence of \var{a}. This 
\begin{equation}
\code{div(a)}=trace(grad(a),where)
\end{equation}
\end{funcdesc}
\begin{funcdesc}{jump}{a\optional{ ,domain=None}}
returns the jump of \var{a} over the discontinuity in its domain or if \Domain \var{domain} is present
in \var{domain}.
\begin{equation}
\begin{array}{rcl}
\code{jump(a)}& = &\code{interpolate(a,FunctionOnContactOne(domain))} \\
              &   & \hfill - \code{interpolate(a,FunctionOnContactZero(domain))}
\end{array}
\end{equation}
\end{funcdesc}
\begin{funcdesc}{L2}{a}
returns the $L^2$-norm of \var{a} in its function space. This is 
\begin{equation}
\code{L2(a)=integrate(length(a)}^2\code{)} \; .
\end{equation} 
\end{funcdesc}

The following functions operate ``point-wise''. That is, the operation is applied to each component of each point
individually.

\begin{funcdesc}{sin}{a}
applies sine function to \var{a}.
\end{funcdesc}

\begin{funcdesc}{cos}{a}
applies cosine function to \var{a}.
\end{funcdesc}

\begin{funcdesc}{tan}{a}
applies tangent function to \var{a}.
\end{funcdesc}

\begin{funcdesc}{asin}{a}
applies arc (inverse) sine function to \var{a}.
\end{funcdesc}

\begin{funcdesc}{acos}{a}
applies arc (inverse) cosine function to \var{a}.
\end{funcdesc}

\begin{funcdesc}{atan}{a}
applies arc (inverse) tangent function to \var{a}.
\end{funcdesc}

\begin{funcdesc}{sinh}{a}
applies hyperbolic sine function to \var{a}.
\end{funcdesc}

\begin{funcdesc}{cosh}{a}
applies hyperbolic cosine function to \var{a}.
\end{funcdesc}

\begin{funcdesc}{tanh}{a}
applies hyperbolic tangent function to \var{a}.
\end{funcdesc}

\begin{funcdesc}{asinh}{a}
applies arc (inverse) hyperbolic sine function to \var{a}.
\end{funcdesc}

\begin{funcdesc}{acosh}{a}
applies arc (inverse) hyperbolic cosine function to \var{a}.
\end{funcdesc}

\begin{funcdesc}{atanh}{a}
applies arc (inverse) hyperbolic tangent function to \var{a}.
\end{funcdesc}

\begin{funcdesc}{exp}{a}
applies exponential function to \var{a}.
\end{funcdesc}

\begin{funcdesc}{sqrt}{a}
applies square root function to \var{a}.
\end{funcdesc}

\begin{funcdesc}{log}{a}
applies the natural logarithm to \var{a}.
\end{funcdesc}

\begin{funcdesc}{log10}{a}
applies the base-$10$ logarithm to \var{a}.
\end{funcdesc}

\begin{funcdesc}{sign}{a}
applies the sign function to \var{a}, that is $1$ where \var{a} is positive,
$-1$ where \var{a} is negative and $0$ otherwise.
\end{funcdesc}

\begin{funcdesc}{wherePositive}{a}
returns a function which is $1$ where \var{a} is positive and $0$ otherwise.
\end{funcdesc}

\begin{funcdesc}{whereNegative}{a}
returns a function which is $1$ where \var{a} is negative and $0$ otherwise.
\end{funcdesc}

\begin{funcdesc}{whereNonNegative}{a}
returns a function which is $1$ where \var{a} is non--negative and $0$ otherwise.
\end{funcdesc}

\begin{funcdesc}{whereNonPositive}{a}
returns a function which is $1$ where \var{a} is non--positive and $0$ otherwise.
\end{funcdesc}

\begin{funcdesc}{whereZero}{a\optional{, tol=None, \optional{, rtol=1.e-8}}}
returns a function which is $1$ where \var{a} equals zero with tolerance \var{tol} and $0$ otherwise. If \var{tol} is not present, the absolute maximum value of C{a} times C{rtol} is used.
\end{funcdesc}

\begin{funcdesc}{whereNonZero}{a, \optional{, tol=None, \optional{, rtol=1.e-8}}}
returns a function which is $1$ where \var{a} different from zero with tolerance \var{tol} and $0$ otherwise. If \var{tol} is not present, the absolute maximum value of C{a} times C{rtol} is used.
\end{funcdesc}

\subsection{Interpolating Data}
\index{interpolateTable}
In some cases, it may be useful to produce Data objects which fit some user defined function.
Manually modifying each value in the Data object is not a good idea since it depends on 
knowing the location and order of each datapoint in the domain.
Instead \escript can use an interpolation table to produce a Data object.

The following example is available as \file{int_save.py} in the examples directory.
We will produce a \Data object which aproximates a sine curve.

\begin{python}
from esys.escript import saveDataCSV, sup
import numpy
from esys.finley import Rectangle

n=4
r=Rectangle(n,n)
x=r.getX()
x0=x[0]
x1=x[1]    #we'll use this later
toobig=100  
\end{python}

First we produce an interpolation table.
\begin{python}
sine_table=[0, 0.70710678118654746, 1, 0.70710678118654746, 0, 
           -0.70710678118654746, -1, -0.70710678118654746, 0]
\end{python}

We wish to identify $0$ and $1$ with the ends of the curve.
That is, with the first and eighth values in the table. 

\begin{python}
numslices=len(sine_table)-1

minval=0
maxval=1

step=sup(maxval-minval)/numslices
\end{python}

So the values $v$ from the input lie in the interval minval$\leq v < $maxval.
\var{step} represents the gap (in the input range) between entries in the table.
By default values of $v$ outside the table argument range (minval, maxval) will
be pushed back into the range, ie. if $v <$ minval the value minval will be used to
evaluate the table. Simularly, for values $v>$ maxval the value maxval is used.

Now we produce our new \Data object.

\begin{python}
result=x0.interpolateTable(sine_table, minval, step, toobig)
\end{python}
Any values which interpolate to larger than \var{toobig} will raise an exception. You can 
switch on boundary checking by adding ''check_boundaries=True`` the argument list.


Now for a 2D example.
We will interpolate a surface such that the bottom edge is the sine curve described above.
The amplitude of the curve decreases as we move towards the top edge.

Our interpolation table will have three rows.
\begin{python}
st=numpy.array(sine_table)

table=[st, 0.5*st, 0*st ]
\end{python}

The use of numpy and multiplication here is just to save typing.

\begin{python}
result2=x1.interpolateTable(table, 0, 0.55, x0, minval, step, toobig)
\end{python}

In the 2D case, the params for the x1 direction (min=0, step=0.55) come first followed by the x0 data object and 
its params.
By default, if a point is specified which is outside the boundary, then \var{interpolateTable} will operate
as if the point was on the boundary.
Passing \var{check_boundaries}=\var{True} will \var{interpolateTable} to reject any points outside the boundaries.

\subsection{Saving Data as CSV}
\index{saveDataCSV}
\index{CSV}
For simple post-processing, \Data objects can be saved in comma separated value format.

If \var{mydata1} and \var{mydata2} are scalar data, the following command:
\begin{python}
saveDataCSV('output.csv',U=mydata1, V=mydata2)
\end{python}
will record the values of mydata in \texttt{output.csv} in the following format:
\begin{verbatim}
U, V
1.0000000e+0, 2.0000000e-1
5.0000000e-0, 1.0000000e+1
...
\end{verbatim}

The names of the keyword parameters form the names of columns in the ouput.
If the data objects are over different function spaces, then saveDataCSV will attempt to
interpolate to a common function space.
If this is not possible, then an exception will be raised.

Output can be restricted using a scalar mask.
\begin{python}
saveDataCSV('outfile.csv', U=mydata1, V=mydata2, mask=myscalar)
\end{python}
Will only output those rows which correspond to to positive values of \var{myscalar}.
Some aspects of the output can be tuned using additional params.
\begin{python}
saveDataCSV('data.csv', append=True, sep=' ', csep='/', mask=mymask, e=mat1)
\end{python}

\begin{itemize}
 \item \var{append} - specifies that the output should be written to the end of an existing file.
\item \var{sep} - defines the separator between fields.
\item \var{csep} - defines the separator between components in the header line. For example between the components of a matrix. 
\end{itemize}

The above command would produce output like this:
\begin{verbatim}
e/0/0 e/1/0 e/0/1 e/1/1
1.0000000000e+00 2.0000000000e+00 3.0000000000e+00 4.0000000000e+00
... 
\end{verbatim}


\subsection{\Operator Class}
The \Operator class provides an abstract access to operators build
within the \LinearPDE class. \Operator objects are created 
when a PDE is handed over to a PDE solver library and handled
by the \LinearPDE object defining the PDE. The user can gain access
to the \Operator of a \LinearPDE object through the \var{getOperator}
method.

\begin{classdesc}{Operator}{}
creates an empty \Operator object.
\end{classdesc}

\begin{methoddesc}[Operator]{isEmpty}{fileName}
returns \True is the object is empty. Otherwise \True is returned.
\end{methoddesc}

\begin{methoddesc}[Operator]{setValue}{value}
resets all entries in the object representation to \var{value}
\end{methoddesc}

\begin{methoddesc}[Operator]{solves}{rhs}
solves the operator equation with right hand side \var{rhs}
\end{methoddesc}

\begin{methoddesc}[Operator]{of}{u}
applies the operator to the \Data object \var{u}
\end{methoddesc}

\begin{methoddesc}[Operator]{saveMM}{fileName}
saves the object to a matrix market format file of name
\var{fileName}, see
\url{http://maths.nist.gov/MatrixMarket}
% \ulink{maths.nist.gov/MatrixMarket}{\url{http://maths.nist.gov/MatrixMarket}}.
\index{Matrix Market}
\end{methoddesc}

\section{Physical Units}
\escript provides support for physical units in the SI system \index{SI units} including unit conversion. So the
user can define variables in the form
\begin{python}
from esys.escript.unitsSI import *
l=20*m
w=30*kg
w2=40*lb
T=100*Celsius
\end{python}
In the two latter cases an conversion from pounds\index{pounds} and degree Celsius\index{Celsius} is performed into the appropriate SI units kg and Kelvin is performed. In addition 
composed units can be used, for instance
\begin{python}
from esys.escript.unitsSI import *
rho=40*lb/cm**3
\end{python}
to define the density in the units of pounds per cubic centimeter. The value $40$ will be converted 
into SI units, in this case kg per cubic meter.
Moreover unit prefixes are supported:
\begin{python}
from esys.escript.unitsSI import *
p=40*Mega*Pa
\end{python}
to the the pressure to 40 Mega Pascal. Units can also be converted back from the SI system into
a desired unit, e.g
\begin{python}
from esys.escript.unitsSI import *
print p/atm
\end{python}
can be used print the pressure in units of atmosphere\index{atmosphere}.  

This is an incomplete list of supported physical units:

\begin{datadesc}{km}
unit of kilo meter
\end{datadesc}

\begin{datadesc}{m}
unit of meter
\end{datadesc}

\begin{datadesc}{cm}
unit of centi meter
\end{datadesc}

\begin{datadesc}{mm}
  unit of milli meter
\end{datadesc}

\begin{datadesc}{sec}
 unit of second
\end{datadesc}

\begin{datadesc}{minute}
  unit of minute
\end{datadesc}

\begin{datadesc}{h}
unit of hour 
\end{datadesc}
\begin{datadesc}{day}
unit of day 
\end{datadesc}
\begin{datadesc}{yr}
 unit of year
\end{datadesc}

\begin{datadesc}{gram}
unit of gram
\end{datadesc}
\begin{datadesc}{kg}
unit of kilo gram
 \end{datadesc}
\begin{datadesc}{lb}
unit of pound 
\end{datadesc}
\begin{datadesc}{ton}
 metric ton
\end{datadesc}

\begin{datadesc}{A}
 unit of Ampere
\end{datadesc}

\begin{datadesc}{Hz}
 unit of Hertz
\end{datadesc}

\begin{datadesc}{N}
 unit of Newton
\end{datadesc}
\begin{datadesc}{Pa}
unit of Pascal 
\end{datadesc}
\begin{datadesc}{atm}
unit of atmosphere 
\end{datadesc}
\begin{datadesc}{J}
unit of Joule 
\end{datadesc}

\begin{datadesc}{W}
unit of Watt 
\end{datadesc}

\begin{datadesc}{C}
unit of Coulomb 
\end{datadesc}
\begin{datadesc}{V}
unit of Volt 
\end{datadesc}
\begin{datadesc}{F}
 unit of Farad 
\end{datadesc}

\begin{datadesc}{Ohm}
 unit of Ohm
\end{datadesc}
\begin{datadesc}{K}
unit of Kelvin 
\end{datadesc}
\begin{datadesc}{Celsius}
 unit of Celsius
\end{datadesc}

\begin{datadesc}{Fahrenheit}
unit of Fahrenheit 
\end{datadesc}

Moreover unit prefixes are supported:

\begin{datadesc}{Yotta}
prefix yotta = $10^{24}$.
 
\end{datadesc}

\begin{datadesc}{Zetta}
prefix zetta= $10^{21}$.
\end{datadesc}

\begin{datadesc}{Exa}
prefix exa= $10^{18}$.
 \end{datadesc}

\begin{datadesc}{Peta}
prefix peta= $10^{15}$.
 \end{datadesc}

\begin{datadesc}{Tera}
prefix tera= $10^{12}$.
 \end{datadesc}

\begin{datadesc}{Giga}
prefix giga= $10^9$.
 \end{datadesc}

\begin{datadesc}{Mega}
prefix mega= $10^6$.
 \end{datadesc}

\begin{datadesc}{Kilo}
prefix kilo= $10^3$.
 \end{datadesc}

\begin{datadesc}{Hecto}
prefix hecto= $10^2$.
 \end{datadesc}

\begin{datadesc}{Deca}
prefix deca= $10^1$.
 \end{datadesc}

\begin{datadesc}{Deci}
prefix deci= $10^{-1}$.
 \end{datadesc}

\begin{datadesc}{Centi}
prefix centi= $10^{-2}$. 
\end{datadesc}

\begin{datadesc}{Milli}
prefix milli= $10^{-3}$.
\end{datadesc}

\begin{datadesc}{Micro}
prefix micro= $10^{-6}$.
 \end{datadesc}

\begin{datadesc}{Nano}
prefix nano= $10^{-9}$.
 \end{datadesc}

\begin{datadesc}{Pico}
prefix pico= $10^{-12}$.
 \end{datadesc}

\begin{datadesc}{Femto}
prefix femto= $10^{-15}$.
 \end{datadesc}

\begin{datadesc}{Atto}
prefix atto= $10^{-18}$.
 \end{datadesc}

\begin{datadesc}{Zepto}
prefix zepto= $10^{-21}$.
 \end{datadesc}

\begin{datadesc}{Yocto}
prefix yocto= $10^{-24}$.
 \end{datadesc}


\section{Utilities}

The \class{FileWriter} provides a mechanism to write data to a file. 
In essence, this class wraps the standard \class{file} class to write data 
that are global in MPI to a file. In fact, data are written on the processor 
with \MPI rank 0 only. It is recommended to use \class{FileWriter} 
rather than \class{open} in order to write code that is running 
with and without \MPI. It is save to use \class{open} under MPI to read data which are global under \MPI.

\begin{classdesc}{FileWriter}{fn\optional{,append=\False, \optional{createLocalFiles=\False}})}
Opens a file of name \var{fn} for writing. If \var{append} is set to \True
written data are append at the end of the file.
If running under \MPI only the first processor with rank==0
will open the file and write to it. 
If \var{createLocalFiles} is set each individual processor will create a file
where for any processor with rank>0 the file name is extended by its rank. This option is normally used for debug purposes only.
\end{classdesc}

The following methods are available:
\begin{methoddesc}[FileWriter]{close}{}
closes the file.
\end{methoddesc}
\begin{methoddesc}[FileWriter]{flush}{}
flushes the internal buffer to disk.
\end{methoddesc}
\begin{methoddesc}[FileWriter]{write}{txt}
Write string \var{txt} to file.
Note that newline is not added.
\end{methoddesc}
\begin{methoddesc}[FileWriter]{writelines}{txts}
Write the list \var{txts} of strings to the file..
Note that newlines are not added. 
This method is equivalent to call write() for each string.
\end{methoddesc}
\begin{memberdesc}[FileWriter]{closed}
\True if file is closed.
\end{memberdesc}
\begin{memberdesc}[FileWriter]{mode}
access mode. 
\end{memberdesc}
\begin{memberdesc}[FileWriter]{name}
file name.
\end{memberdesc}
\begin{memberdesc}[FileWriter]{newlines}
line separator
\end{memberdesc}


\begin{funcdesc}{setEscriptParamInt}{name,value}
assigns the integer value \var{value} to the parameter \var{name}.
If \var{name}="TOO_MANY_LINES" conversion of any \Data object to a string switches to a 
condensed format if more than \var{value} lines would be created.
\end{funcdesc}

\begin{funcdesc}{getEscriptParamInt}{name}
returns the current value of integer parameter \var{name}. 
\end{funcdesc}

\begin{funcdesc}{listEscriptParams}{a}
returns a list of valid parameters and their description.
\end{funcdesc}

\begin{funcdesc}{getMPISizeWorld}{}
returns the number of \MPI processors in use in the \env{MPI_COMM_WORLD} processor group.
If \MPI is not used 1 is returned.
\end{funcdesc}
\begin{funcdesc}{getMPIRankWorld}{}
returns the rank of the process within the \env{MPI_COMM_WORLD} processor group.
If \MPI is not used 0 is returned.
\end{funcdesc}
\begin{funcdesc}{MPIBarrierWorld}{}
performs a barrier synchronization across all processors within \env{MPI_COMM_WORLD}
processor group.
\end{funcdesc}
\begin{funcdesc}{getMPIWorldMax}{a}
returns the maximum value of the integer \var{a} across all 
processors within \env{MPI_COMM_WORLD}.
\end{funcdesc}

Name	Value
svn:eol-style	native
svn:keywords	Author Date Id Revision