# Contents of /trunk/doc/user/escript.tex

Revision 6923 - (show annotations)
Mon Dec 2 05:57:07 2019 UTC (5 months, 3 weeks ago) by aellery
File MIME type: application/x-tex
File size: 89629 byte(s)
- Amended the userguide to include ComplexScalar, ComplexVector etc.
- Fixed a bug in run_comm1 & run_comm4 (failure if escript was compiled without scipy)
- Temporarily removed the ability to interpolate from ReducedFunction to Function in dudley and finley


 1 2 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 3 % Copyright (c) 2003-2018 by The University of Queensland 4 5 % 6 % Primary Business: Queensland, Australia 7 % Licensed under the Apache License, version 2.0 8 9 % 10 % Development until 2012 by Earth Systems Science Computational Center (ESSCC) 11 % Development 2012-2013 by School of Earth Sciences 12 % Development from 2014 by Centre for Geoscience Computing (GeoComp) 13 % 14 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 15 16 \chapter{The \escript Module}\label{ESCRIPT CHAP} 17 18 \section{Concepts} 19 \escript is a \PYTHON module that allows you to represent the values of 20 a function at points in a \Domain in such a way that the function will 21 be useful for the Finite Element Method (FEM) simulation. It also 22 provides what we call a function space that describes how the data is 23 used in the simulation. Stored along with the data is information 24 about the elements and nodes which will be used by the domain (e.g. \finley). 25 26 \subsection{Function spaces} 27 In order to understand what we mean by the term 'function space', 28 consider that the solution of a partial differential 29 equation\index{partial differential equation} (PDE) is a function on a domain 30 $\Omega$. When solving a PDE using FEM, the solution is 31 piecewise-differentiable but, in general, its gradient is discontinuous. 32 To reflect these different degrees of smoothness, different function spaces 33 are used. 34 For instance, in FEM, the displacement field is represented by its values at 35 the nodes of the mesh, and so is continuous. 36 The strain, which is the symmetric part of the gradient of the displacement 37 field, is stored on the element centers, and so is considered to be 38 discontinuous. 39 40 A function space is described by a \FunctionSpace object. 41 The following statement generates the object \var{solution_space} which is 42 a \FunctionSpace object and provides access to the function space of 43 PDE solutions on the \Domain \var{mydomain}: 44 45 \begin{python} 46 solution_space=Solution(mydomain) 47 \end{python} 48 The following generators for function spaces on a \Domain \var{mydomain} are commonly used: 49 \begin{itemize} 50 \item \var{Solution(mydomain)}: solutions of a PDE 51 \item \var{ReducedSolution(mydomain)}: solutions of a PDE with a reduced 52 smoothness requirement, e.g. using a lower order approximation on the same 53 element or using macro elements\index{macro elements} 54 \item \var{ContinuousFunction(mydomain)}: continuous functions, e.g. a temperature distribution 55 \item \var{Function(mydomain)}: general functions which are not necessarily continuous, e.g. a stress field 56 \item \var{FunctionOnBoundary(mydomain)}: functions on the boundary of the domain, e.g. a surface pressure 57 \item \var{DiracDeltaFunctions(mydomain)}: functions defined on a set of points 58 \item \var{FunctionOnContact0(mydomain)}: functions on side $0$ of a discontinuity 59 \item \var{FunctionOnContact1(mydomain)}: functions on side $1$ of a discontinuity 60 \end{itemize} 61 In some cases under-integration is used. For these cases the user may use a 62 \FunctionSpace from the following list: 63 \begin{itemize} 64 \item \var{ReducedFunction(mydomain)} 65 \item \var{ReducedFunctionOnBoundary(mydomain)} 66 \item \var{ReducedFunctionOnContact0(mydomain)} 67 \item \var{ReducedFunctionOnContact1(mydomain)} 68 \end{itemize} 69 In comparison to the corresponding full version they use a reduced number of 70 integration nodes (typically one only) to represent values. 71 72 \begin{figure} 73 \centering 74 \scalebox{0.97}{\includegraphics{EscriptDiagram1}} 75 \caption{\label{ESCRIPT DEP}Dependency of function spaces in \finley. 76 An arrow indicates that a function in the \FunctionSpace at the starting point 77 can be interpolated to the \FunctionSpace of the arrow target. 78 All function spaces above the dotted line can be interpolated to any of 79 the function spaces below the line. See also \Sec{SEC Projection}.} 80 \end{figure} 81 82 The reduced smoothness for a PDE solution is often used to fulfill the 83 Ladyzhenskaya-Babuska-Brezzi condition~\cite{LBB} when solving saddle point 84 problems\index{saddle point problems}, e.g. the Stokes equation. 85 A discontinuity\index{discontinuity} is a region within the domain across 86 which functions may be discontinuous. 87 The location of a discontinuity is defined in the \Domain object. 88 \fig{ESCRIPT DEP} shows the dependency between the types of function spaces 89 in \finley (other libraries may have different relationships). 90 91 The solution of a PDE is a continuous function. Any continuous function can 92 be seen as a general function on the domain and can be restricted to the 93 boundary as well as to one side of a discontinuity (the result will be 94 different depending on which side is chosen). Functions on any side of the 95 discontinuity can be seen as a function on the corresponding other side. 96 97 A function on the boundary or on one side of the discontinuity cannot be seen 98 as a general function on the domain as there are no values defined for the 99 interior. For most PDE solver libraries the space of the solution and 100 continuous functions is identical, however in some cases, for example when 101 periodic boundary conditions are used in \finley, a solution fulfills periodic 102 boundary conditions while a continuous function does not have to be periodic. 103 104 The concept of function spaces describes the properties of functions and 105 allows abstraction from the actual representation of the function in the 106 context of a particular application. For instance, in the FEM context a 107 function of the \Function type (written as \emph{Function()} in \fig{ESCRIPT DEP}) 108 is usually represented by its values at the element center, 109 but in a finite difference scheme the edge midpoint of cells is preferred. 110 By changing its function space you can use the same function in a Finite 111 Difference scheme instead of Finite Element scheme. 112 Changing the function space of a particular function will typically lead to 113 a change of its representation. 114 So, when seen as a general function, a continuous function which is typically 115 represented by its values on the nodes of the FEM mesh or finite difference 116 grid must be interpolated to the element centers or the cell edges, 117 respectively. Interpolation happens automatically in \escript whenever it is 118 required\index{interpolation}. The user needs to be aware that an 119 interpolation is not always possible, see \fig{ESCRIPT DEP} for \finley. 120 An alternative approach to change the representation (=\FunctionSpace) is 121 projection\index{projection}, see \Sec{SEC Projection}. 122 123 \subsection{\Data Objects} 124 In \escript the class that stores these functions is called \Data. 125 The function is represented through its values on \DataSamplePoints where 126 the \DataSamplePoints are chosen according to the function space of the 127 function. 128 \Data class objects are used to define the coefficients of the PDEs to be 129 solved by a PDE solver library and also to store the solutions of the PDE. 130 131 The values of the function have a rank which gives the number of indices, 132 and a \Shape defining the range of each index. 133 The rank in \escript is limited to the range 0 through 4 and it is assumed 134 that the rank and \Shape is the same for all \DataSamplePoints. 135 The \Shape of a \Data object is a tuple (list) \var{s} of integers. 136 The length of \var{s} is the rank of the \Data object and the \var{i}-th 137 index ranges between 0 and $\var{s[i]}-1$. 138 For instance, a stress field has rank 2 and \Shape $(d,d)$ where $d$ is the 139 number of spatial dimensions. 140 The following statement creates the \Data object \var{mydat} representing a 141 continuous function with values of \Shape $(2,3)$ and rank $2$: 142 \begin{python} 143 mydat=Data(value=1, what=ContinuousFunction(myDomain), shape=(2,3)) 144 \end{python} 145 The initial value is the constant 1 for all \DataSamplePoints and all 146 components. 147 148 \Data objects can also be created from any \numpy array or any object, such 149 as a list of floating point numbers, that can be converted into 150 a \numpyNDA\cite{NUMPY}. 151 The following two statements create objects which are equivalent 152 to \var{mydat}: 153 \begin{python} 154 mydat1=Data(value=numpy.ones((2,3)), what=ContinuousFunction(myDomain)) 155 mydat2=Data(value=[[1,1], [1,1], [1,1]], what=ContinuousFunction(myDomain)) 156 \end{python} 157 In the first case the initial value is \var{numpy.ones((2,3))} which generates 158 a $2 \times 3$ matrix as an instance of \numpyNDA filled with ones. 159 The \Shape of the created \Data object is taken from the \Shape of the array. 160 In the second case, the creator converts the initial value, which is a list of 161 lists, into a \numpyNDA before creating the actual \Data object. 162 163 For convenience \escript provides creators for the most common types 164 of \Data objects in the following forms (\var{d} defines the spatial 165 dimensionality): 166 \begin{itemize} 167 \item \code{Scalar(0, Function(mydomain))} is the same as \code{Data(0, Function(myDomain),(,))}\\ 168 (each value is a scalar), e.g. a temperature field 169 \item \code{Vector(0, Function(mydomain))} is the same as \code{Data(0, Function(myDomain),(d,))}\\ 170 (each value is a vector), e.g. a velocity field 171 \item \code{Tensor(0, Function(mydomain))} equals \code{Data(0, Function(myDomain), (d,d))}, 172 e.g. a stress field 173 \item \code{Tensor4(0,Function(mydomain))} equals \code{Data(0,Function(myDomain), (d,d,d,d))}, 174 e.g. a Hook tensor field 175 \item \code{ComplexScalar(0+0j, Function(mydomain))} is the same as \\ 176 \code{ComplexData(0+0j, Function(myDomain),(,))} 177 (each value is a complex scalar), e.g. a temperature field 178 \item \code{ComplexVector(0+0j, Function(mydomain))} is the same as \\ 179 \code{ComplexData(0+0j, Function(myDomain),(d,))} 180 (each value is a complex vector), e.g. a velocity field 181 \item \code{ComplexTensor(0+0j, Function(mydomain))} is the same as \\ 182 \code{ComplexData(0+0j, Function(myDomain), (d,d))}, e.g. a stress field 183 \item \code{ComplexTensor4(0+0j,Function(mydomain))} is the same as \\ 184 \code{ComplexData(0+0j,Function(myDomain), (d,d,d,d))}, e.g. a Hook tensor field 185 \end{itemize} 186 Here the initial value is 0 but any object that can be converted into 187 a \numpyNDA and whose \Shape is consistent with \Shape of the \Data object to 188 be created can be used as the initial value. 189 190 \Data objects can be manipulated by applying unary operations (e.g. cos, sin, 191 log), and they can be combined point-wise by applying arithmetic operations 192 (e.g. +, - ,* , /). 193 We emphasize that \escript itself does not handle any spatial dependencies as 194 it does not know how values are interpreted by the processing PDE solver library. 195 However \escript invokes interpolation if this is needed during data manipulations. 196 Typically, this occurs in binary operations when the arguments belong to 197 different function spaces or when data are handed over to a PDE solver library 198 which requires functions to be represented in a particular way. 199 200 The following example shows the usage of \Data objects. Assume we have a 201 displacement field $u$ and we want to calculate the corresponding stress field 202 $\sigma$ using the linear-elastic isotropic material model 203 \begin{eqnarray}\label{eq: linear elastic stress} 204 \sigma_{ij}=\lambda u_{k,k} \delta_{ij} + \mu ( u_{i,j} + u_{j,i}) 205 \end{eqnarray} 206 where $\delta_{ij}$ is the Kronecker symbol and 207 $\lambda$ and $\mu$ are the Lam\'e coefficients. The following function 208 takes the displacement \var{u} and the Lam\'e coefficients \var{lam} and \var{mu} 209 as arguments and returns the corresponding stress: 210 \begin{python} 211 from esys.escript import * 212 def getStress(u, lam, mu): 213 d=u.getDomain().getDim() 214 g=grad(u) 215 stress=lam*trace(g)*kronecker(d)+mu*(g+transpose(g)) 216 return stress 217 \end{python} 218 The variable \var{d} gives the spatial dimensionality of the domain on which 219 the displacements are defined. 220 The \code{kronecker(d)} call, returns the Kronecker symbol with indices $i$ and $j$ running 221 from 0 to \var{d}-1. 222 The \var{grad(u)} call, requires the displacement field \var{u} to be in 223 the \var{Solution} or \ContinuousFunction. 224 The result \var{g} as well as the returned stress will be in the \Function. 225 If, for example, \var{u} is the solution of a PDE then \code{getStress} might 226 be called in the following way: 227 \begin{python} 228 s=getStress(u, 1., 2.) 229 \end{python} 230 However \code{getStress} can also be called with \Data objects as values for 231 \var{lam} and \var{mu} which, for instance in the case of a temperature 232 dependency, are calculated by an expression. 233 The following call is equivalent to the previous example: 234 \begin{python} 235 lam=Scalar(1., ContinuousFunction(mydomain)) 236 mu=Scalar(2., Function(mydomain)) 237 s=getStress(u, lam, mu) 238 \end{python} 239 % 240 The function \var{lam} belongs to the \ContinuousFunction but with \var{g} the 241 function \var{trace(g)} is in the \Function. 242 In the evaluation of the product \var{lam*trace(g)} we have different function 243 spaces (on the nodes versus in the centers) and at first glance we have incompatible data. 244 \escript converts the arguments into an appropriate function space according 245 to \fig{ESCRIPT DEP}. 246 In this example that means \escript sees \var{lam} as a function of the \Function. 247 In the context of FEM this means the nodal values of \var{lam} are 248 interpolated to the element centers. 249 The interpolation is automatic and requires no special handling. 250 251 \begin{figure} 252 \centering 253 \includegraphics{EscriptDiagram2} 254 \caption{\label{Figure: tag}Element Tagging. A rectangular mesh over a region 255 with two rock types {\it white} and {\it gray} is shown. 256 The number in each cell refers to the major rock type present in the cell 257 ($1$ for {\it white} and $2$ for {\it gray}).} 258 \end{figure} 259 260 \subsection{Tagged, Expanded and Constant Data} 261 Material parameters such as the Lam\'e coefficients are typically dependent on 262 rock types present in the area of interest. 263 A common technique to handle these kinds of material parameters is 264 \emph{tagging}\index{tagging}, which uses storage efficiently. 265 \fig{Figure: tag} shows an example. In this case two rock types {\it white} 266 and {\it gray} can be found in the domain. 267 The domain is subdivided into triangular shaped cells. 268 Each cell has a tag indicating the rock type predominantly found in this cell. 269 Here $1$ is used to indicate rock type {\it white} and $2$ for rock type {\it gray}. 270 The tags are assigned at the time when the cells are generated and stored in 271 the \Domain class object. To allow easier usage of tags, names can be used 272 instead of numbers. These names are typically defined at the time when the 273 geometry is generated. 274 275 The following statements show how to use tagged values for \var{lam} as shown 276 in \fig{Figure: tag} for the stress calculation discussed above: 277 \begin{python} 278 lam=Scalar(value=2., what=Function(mydomain)) 279 insertTaggedValue(lam, white=30., gray=5000.) 280 s=getStress(u, lam, 2.) 281 \end{python} 282 In this example \var{lam} is set to $30$ for those cells with tag {\it white} 283 (=$1$) and to $5000$ for cells with tag {\it gray} (=$2$). 284 The initial value $2$ of \var{lam} is used as a default value for the case 285 when a tag is encountered which has not been linked with a value. 286 The \code{getStress} method does not need to be changed now that we are using tags. 287 \escript resolves the tags when \var{lam*trace(g)} is calculated. 288 289 This brings us to a very important point about \escript. 290 You can develop a simulation with constant Lam\'e coefficients, and then later 291 switch to tagged Lam\'e coefficients without otherwise changing your \PYTHON script. 292 In short, you can use the same script for models with different domains and 293 different types of input data. 294 295 There are three main ways in which \Data objects are represented internally -- 296 constant, tagged, and expanded. 297 In the constant case, the same value is used at each sample point while only a 298 single value is stored to save memory. 299 In the expanded case, each sample point has an individual value (such as for the solution of a PDE). 300 This is where your largest data sets will be created because the values are 301 stored as a complete array. 302 The tagged case has already been discussed above. 303 Expanded data is created when specifying \code{expanded=True} in the \Data 304 object constructor, while tagged data requires calling the \member{insertTaggedValue} 305 method as shown above. 306 307 Values are accessed through a sample reference number. 308 Operations on expanded \Data objects have to be performed for each sample 309 point individually. 310 When tagged values are used, the values are held in a dictionary. 311 Operations on tagged data require processing the set of tagged values only, 312 rather than processing the value for each individual sample point. 313 \escript allows any mixture of constant, tagged and expanded data in a single expression. 314 315 \subsection{Saving and Restoring Simulation Data} 316 \Data objects can be written to disk files with the \member{dump} method and 317 read back using the \member{load} method, both of which use the 318 \netCDF\cite{NETCDF} file format. 319 Use these to save data for checkpoint/restart or simply to save and reuse data 320 that was expensive to compute. 321 For instance, to save the coordinates of the data points of a 322 \ContinuousFunction to the file \file{x.nc} use 323 \begin{python} 324 x=ContinuousFunction(mydomain).getX() 325 x.dump("x.nc") 326 mydomain.dump("dom.nc") 327 \end{python} 328 To recover the object \var{x}, and you know that \var{mydomain} was an \finley 329 mesh, use 330 \begin{python} 331 from esys.finley import LoadMesh 332 mydomain=LoadMesh("dom.nc") 333 x=load("x.nc", mydomain) 334 \end{python} 335 Obviously, it is possible to execute the same steps that were originally used 336 to generate \var{mydomain} to recreate it. However, in most cases using 337 \member{dump} and \member{load} is faster, particularly if optimization has 338 been applied. 339 If \escript is running on more than one \MPI process \member{dump} will create 340 an individual file for each process containing the local data. 341 In order to avoid conflicts the \MPI processor 342 rank is appended to the file names. 343 That is instead of one file \file{dom.nc} you would get 344 \file{dom.nc.0000}, \file{dom.nc.0001}, etc. 345 You still call \code{LoadMesh("dom.nc")} to load the domain but you have to 346 make sure that the appropriate file is accessible from the corresponding rank, 347 and loading will only succeed if you run with as many processes as were used 348 when calling \member{dump}. 349 350 The function space of the \Data is stored in \file{x.nc}. 351 If the \Data object is expanded, the number of data points in the file and of 352 the \Domain for the particular \FunctionSpace must match. 353 Moreover, the ordering of the values is checked using the reference 354 identifiers provided by the \FunctionSpace on the \Domain. 355 In some cases, data points will be reordered so be aware and confirm that you 356 get what you wanted. 357 358 A more flexible way of saving and restoring \escript simulation data 359 is through an instance of the \class{DataManager} class. 360 It has the advantage of allowing to save and load not only a \Domain and 361 \Data objects but also other values\footnote{The \PYTHON \emph{pickle} module 362 is used for other types.} you compute in your simulation script. 363 Further, \class{DataManager} objects can simultaneously create files for 364 visualization so no extra calls to \code{saveVTK} etc. are needed. 365 366 The following example shows how the \class{DataManager} class can be used. 367 For an explanation of all member functions and options see the class reference 368 Section \ref{sec:datamanager}. 369 \begin{python} 370 from esys.escript import DataManager, Scalar, Function 371 from esys.finley import Rectangle 372 373 dm = DataManager(formats=[DataManager.RESTART, DataManager.VTK]) 374 if dm.hasData(): 375 mydomain=dm.getDomain() 376 val=dm.getValue("val") 377 t=dm.getValue("t") 378 t_max=dm.getValue("t_max") 379 else: 380 mydomain=Rectangle() 381 val=Function(mydomain).getX() 382 t=0. 383 t_max=2.5 384 385 while t$\var{maxval} the value \var{maxval} is used. 1449 1450 Now we produce our new \Data object: 1451 1452 \begin{python} 1453 result=interpolateTable(sine_table, x[0], minval, step, toobig) 1454 \end{python} 1455 Any values which interpolate to larger than \var{toobig} will raise an 1456 exception. You can switch on boundary checking by adding 1457 \code{check_boundaries=True} to the argument list. 1458 1459 Now consider a 2D example. We will interpolate from a plane where$\forall x,y\in[0,9]:(x,y)=x+y\cdot10$. 1460 1461 \begin{python} 1462 from esys.escript import whereZero 1463 table2=[] 1464 for y in range(0,10): 1465 r=[] 1466 for x in range(0,10): 1467 r.append(x+y*10) 1468 table2.append(r) 1469 xstep=(maxval-minval)/(10-1) 1470 ystep=(maxval-minval)/(10-1) 1471 1472 xmin=minval 1473 ymin=minval 1474 1475 result2=interpolateTable(table2, x2, (xmin, ymin), (xstep, ystep), toobig) 1476 \end{python} 1477 1478 We can check the values using \function{whereZero}. 1479 For example, for$x=0$: 1480 \begin{python} 1481 print(result2*whereZero(x[0])) 1482 \end{python} 1483 1484 Finally let us look at a 3D example. Note that the parameter tuples should be 1485$(x,y,z)$but that in the interpolation table,$x$is the innermost dimension. 1486 \begin{python} 1487 b=Brick(n,n,n) 1488 x3=b.getX() 1489 toobig=1000000 1490 1491 table3=[] 1492 for z in range(0,10): 1493 face=[] 1494 for y in range(0,10): 1495 r=[] 1496 for x in range(0,10): 1497 r.append(x+y*10+z*100) 1498 face.append(r) 1499 table3.append(face); 1500 1501 zstep=(maxval-minval)/(10-1) 1502 1503 zmin=minval 1504 1505 result3=interpolateTable(table3, x3, (xmin, ymin, zmin), 1506 (xstep, ystep, zstep), toobig) 1507 \end{python} 1508 1509 1510 \subsubsection{Non-uniform Interpolation} 1511 Non-uniform interpolation is also supported for the one dimensional case. 1512 \begin{python} 1513 Data.nonuniformInterpolate(in, out, check_boundaries) 1514 Data.nonuniformSlope(in, out, check_boundaries) 1515 \end{python} 1516 1517 Will produce a new \Data object by mapping the given \Data object through the user-defined function 1518 specified by \texttt{in} and \texttt{out}. 1519 The \ldots Interpolate version gives the value of the function at the specified point and the 1520 \ldots Slope version gives the slope at those points. 1521 The check_boundaries boolean argument specifies what the function should do if the \Data object contains 1522 values outside the range specified by the \texttt{in} parameter. 1523 If the argument is \texttt{False}, then those datapoints will be interpolated to the value of the edge 1524 they are closest to (or assigned a slope of zero). 1525 If the argument is \texttt{True}, then an exception will be thrown if out of bounds values are detected. 1526 Note that the values given by the \texttt{in} parameter must be monotonically increasing. 1527 1528 \noindent For example:\\ 1529 If \texttt{d} contains the values \texttt{\{1,2,3,4,5\}}, then 1530 \begin{python} 1531 d.nonuniformInterpolate([1.5, 2, 2.8, 4.6], [4, 5, -1, 1], False) 1532 \end{python} 1533 would produce a \Data object containing \texttt{\{4, 5, -0.7777, 0.3333, 1\}}.\\ 1534 A similar call to \texttt{nonuniformSlope} would produce a \Data object containing \texttt{\{0, 2, 1.1111, 1.1111, 0\}}. 1535 % 1536 % 1537 % We will interpolate a surface such that the bottom 1538 % edge is the sine curve described above. 1539 % The amplitude of the curve decreases as we move towards the top edge. 1540 % Our interpolation table will have three rows: 1541 % 1542 % \begin{python} 1543 % st=numpy.array(sine_table) 1544 % table=[st, 0.5*st, 0*st] 1545 % \end{python} 1546 % % 1547 % The use of \numpy and multiplication here is just to save typing. 1548 % 1549 % % result2=x1.interpolateTable(table, 0, 0.55, x0, minval, step, toobig) 1550 % \begin{python} 1551 % result=interpolateTable(table, x (minval,0), (0.55, step), toobig) 1552 % \end{python} 1553 % 1554 % In the 2D case the start and step parameters are tuples$(x,y)$. 1555 % By default, if a point is specified which is outside the boundary, then 1556 % \var{interpolateTable} will operate as if the point was on the boundary. 1557 % Passing \code{check_boundaries=True} will lead to the rejection of any points 1558 % outside the boundaries by \var{interpolateTable}. 1559 % 1560 % This method can also be called with three dimensional tables and \Data objects. 1561 % Tuples should be ordered$(x,y,z)$. 1562 1563 \subsection{The \var{DataManager} Class} 1564 \label{sec:datamanager} 1565 1566 The \var{DataManager} class can be used to conveniently add checkpoint/restart 1567 functionality to \escript simulations. 1568 Once an instance is created \Data objects and other values can be added and 1569 dumped to disk by a single method call. 1570 If required the object can be set up to also save the data in a format suitable 1571 for visualization. 1572 Internally the \var{DataManager} interfaces with \weipa for this. 1573 1574 \begin{classdesc}{DataManager}{formats=[RESTART], work_dir=".", restart_prefix="restart", do_restart=\True} 1575 initializes a new \var{DataManager} object which can be used to save, 1576 restore and export simulation data in a number of formats. 1577 All files and directories saved or restored by this object are located 1578 under the directory specified by \var{work_dir}. 1579 If \var{RESTART} is specified in \var{formats}, the \var{DataManager} will 1580 look for directories whose name starts with \var{restart_prefix}. 1581 In case \var{do_restart} is \True, the last of these directories is used 1582 to restore simulation data while all others are deleted. 1583 If \var{do_restart} is \False, then all of those directories are deleted. 1584 The \var{restart_prefix} and \var{do_restart} parameters are ignored if 1585 \var{RESTART} is not specified in \var{formats}. 1586 \end{classdesc} 1587 1588 \noindent Valid values for the \var{formats} parameter are: 1589 \begin{memberdesc}[DataManager]{RESTART} 1590 enables writing of checkpoint files to be able to continue simulations 1591 as explained in the class description. 1592 \end{memberdesc} 1593 \begin{memberdesc}[DataManager]{SILO} 1594 exports simulation data in the \SILO file format. \escript must have 1595 been compiled with \SILO support for this to work. 1596 \end{memberdesc} 1597 \begin{memberdesc}[DataManager]{VISIT} 1598 enables the \VisIt simulation interface which allows connecting to and 1599 interacting with the running simulation from a compatible \VisIt client. 1600 \escript must have been compiled with \VisIt (version 2) support and the 1601 version of the client has to match the version used at compile time. 1602 In order to connect to the simulation the client needs to have access and 1603 load the file \file{escriptsim.sim2} located under the work directory. 1604 \end{memberdesc} 1605 \begin{memberdesc}[DataManager]{VTK} 1606 exports simulation data in the \VTK file format. 1607 \end{memberdesc} 1608 1609 \noindent The \var{DataManager} class has the following methods: 1610 \begin{methoddesc}[DataManager]{addData}{**data} 1611 adds \Data objects and other data to the manager. Calling this method does 1612 not save or export the data yet so it is allowed to incrementally add data 1613 at various points in the simulation script if required. 1614 Note, that only a single domain is supported so all \Data objects have to 1615 be defined on the same one or an exception is raised. 1616 \end{methoddesc} 1617 1618 \begin{methoddesc}[DataManager]{setDomain}{domain} 1619 explicitly sets the domain for this manager. 1620 It is generally not required to call this method directly. 1621 Instead, the \var{addData} method will set the domain used by the \Data 1622 objects. 1623 An exception is raised if the domain was set to a different domain before 1624 (explicitly or implicitly). 1625 \end{methoddesc} 1626 1627 \begin{methoddesc}[DataManager]{hasData}{} 1628 returns \True if the manager has loaded simulation data for a restart. 1629 \end{methoddesc} 1630 1631 \begin{methoddesc}[DataManager]{getDomain}{} 1632 returns the domain as recovered from a restart. 1633 \end{methoddesc} 1634 1635 \begin{methoddesc}[DataManager]{getValue}{value_name} 1636 returns a \Data object or other value with the name \var{value_name} that 1637 has been recovered after a restart. 1638 \end{methoddesc} 1639 1640 \begin{methoddesc}[DataManager]{getCycle}{} 1641 returns the export cycle, i.e. the number of times \var{export()} has been 1642 called. 1643 \end{methoddesc} 1644 1645 \begin{methoddesc}[DataManager]{setCheckpointFrequency}{freq} 1646 sets the frequency with which checkpoint files are created. This is only 1647 useful if the \var{DataManager} object was created with at least one other 1648 format next to \var{RESTART}. The frequency is 1 by default which means 1649 that checkpoint files are created every time \var{export()} is called. 1650 Unlike visualization output, a simulation checkpoint is usually not 1651 required at every time step. Thus, the frequency can be decreased by 1652 calling this method with$\var{freq}>1$which would then create restart 1653 files every \var{freq} times \var{export()} is called. 1654 \end{methoddesc} 1655 1656 \begin{methoddesc}[DataManager]{setTime}{time} 1657 sets the simulation time stamp. This floating point number is stored in 1658 the metadata of exported data but not used by \var{RESTART}. 1659 \end{methoddesc} 1660 1661 \begin{methoddesc}[DataManager]{setMeshLabels}{x, y, z=""} 1662 sets labels for the mesh axes. These are currently only used by the \SILO 1663 exporter. 1664 \end{methoddesc} 1665 1666 \begin{methoddesc}[DataManager]{setMeshUnits}{x, y, z=""} 1667 sets units for the mesh axes. These are currently only used by the \SILO 1668 exporter. 1669 \end{methoddesc} 1670 1671 \begin{methoddesc}[DataManager]{setMetadataSchemaString}{schema, metadata=""} 1672 sets metadata namespaces and the corresponding metadata. These are 1673 currently only used by the \VTK exporter. 1674 \var{schema} is a dictionary that maps prefixes to namespace names, e.g.\\ 1675 \code{\{"gml": "http://www.opengis.net/gml"\}} and \var{metadata} is a 1676 string with the actual content which will be enclosed in \var{} 1677 tags. 1678 \end{methoddesc} 1679 1680 \begin{methoddesc}[DataManager]{export}{} 1681 executes the actual data export. Depending on the \var{formats} parameter 1682 used in the constructor all data added by \var{addData()} is written to 1683 disk (\var{RESTART,SILO,VTK}) or made available through the \VisIt 1684 simulation interface (\var{VISIT}). 1685 At least the domain must be set for something to be exported. 1686 \end{methoddesc} 1687 1688 \subsection{Saving Data as CSV} 1689 \label{sec:savedatacsv} 1690 \index{saveDataCSV}\index{CSV} 1691 For simple post-processing, \Data objects can be saved in comma separated 1692 value (\emph{CSV}) format. 1693 If \var{mydata1} and \var{mydata2} are scalar data, the command 1694 \begin{python} 1695 saveDataCSV('output.csv', U=mydata1, V=mydata2) 1696 \end{python} 1697 will record the values in \file{output.csv} in the following format: 1698 \begin{verbatim} 1699 U, V 1700 1.0000000e+0, 2.0000000e-1 1701 5.0000000e-0, 1.0000000e+1 1702 ... 1703 \end{verbatim} 1704 1705 The names of the keyword parameters form the names of columns in the output. 1706 If the data objects are over different function spaces, then \var{saveDataCSV} 1707 will attempt to interpolate to a common function space. 1708 If this is not possible, then an exception is raised. 1709 1710 Output can be restricted using a scalar mask as follows: 1711 \begin{python} 1712 saveDataCSV('outfile.csv', U=mydata1, V=mydata2, mask=myscalar) 1713 \end{python} 1714 This command will only output those rows which correspond to to positive 1715 values of \var{myscalar}. 1716 Some aspects of the output can be tuned using additional parameters: 1717 \begin{python} 1718 saveDataCSV('data.csv', refid=True, append=True, sep=' ', csep='/', mask=mymask, e=mat1) 1719 \end{python} 1720 1721 \begin{itemize} 1722 \item \var{refid} -- specifies that the output should include the reference IDs of the elements or nodes 1723 \item \var{append} -- specifies that the output should be written to the end of an existing file 1724 \item \var{sep} -- defines the separator between fields 1725 \item \var{csep} -- defines the separator between components in the header 1726 line. For example between the components of a matrix. 1727 \end{itemize} 1728 % 1729 The above command would produce output like this: 1730 \begin{verbatim} 1731 refid e/0/0 e/1/0 e/0/1 e/1/1 1732 0 1.0000000000e+00 2.0000000000e+00 3.0000000000e+00 4.0000000000e+00 1733 ... 1734 \end{verbatim} 1735 1736 Note that while the order in which rows are output can vary, all the elements 1737 in a given row always correspond to the same input. 1738 1739 \subsection{Converting \Data to a Numpy Array} 1740 \label{sec:getnumpy} 1741 \index{getNumpy}\index{GN} 1742 \Data objects can be converted into a numpy structured array using the commands \var{getNumpy} and \var{convertNumpy}. 1743 \subsubsection{getNumpy} 1744 If \var{mydata1} and \var{mydata2} are scalar \Data, then the command 1745 \begin{python} 1746 a,b = getNumpy(U=mydata1, V=mydata2) 1747 \end{python} 1748 will return two structured ndarrays with the names '\emph{U}' and '\emph{V}'. 1749 \begin{verbatim} 1750 a['U'] = [1.0000000e+0, 2.0000000e-1, ... 1751 b['V'] = [2.0000000e+0, 3.0000000e-1, ... 1752 \end{verbatim} 1753 1754 Up to five \Data objects can be passed to \var{getNumpy} at the time. These objects can be scalar, vector or tensor \Data objects. The names of the keyword parameters form the names of the returned arrays. 1755 If the data objects are over different function spaces, then \var{getNumpy} 1756 will attempt to interpolate to a common function space. 1757 If this is not possible, then an exception is raised. 1758 1759 Output can be restricted using a scalar mask as follows: 1760 \begin{python} 1761 a,b,c = getNumpy(U=mydata1, V=mydata2, W=mydata3, mask=myscalar) 1762 \end{python} 1763 This command will only output those rows which correspond to to positive 1764 values of \var{myscalar}. 1765 1766 Note that while the order in which output rows are output can vary, all the elements 1767 in a given row always correspond to the same input. 1768 1769 \subsubsection{convertNumpy} 1770 \Data objects can also be converted into a numpy structured array using the command \var{convertNumpy}. 1771 If \var{mydata1} is a \Data object, then the command 1772 \begin{python} 1773 a = convertNumpy(mydata1) 1774 \end{python} 1775 will return a structured ndarray containing all of the data in \var{mydata1}. Unlike \var{getNumpy}, this function 1776 does not support the use of masks and does not use MPI. 1777 1778 \subsection{The \Operator Class} 1779 The \Operator class provides an abstract access to operators built 1780 within the \LinearPDE class. \Operator objects are created 1781 when a PDE is handed over to a PDE solver library and handled 1782 by the \LinearPDE object defining the PDE. The user can gain access 1783 to the \Operator of a \LinearPDE object through the \var{getOperator} 1784 method. 1785 1786 \begin{classdesc}{Operator}{} 1787 creates an empty \Operator object. 1788 \end{classdesc} 1789 1790 \begin{methoddesc}[Operator]{isEmpty}{fileName} 1791 returns \True is the object is empty, \False otherwise. 1792 \end{methoddesc} 1793 1794 \begin{methoddesc}[Operator]{resetValues}{} 1795 resets all entries in the operator. 1796 \end{methoddesc} 1797 1798 \begin{methoddesc}[Operator]{solve}{rhs} 1799 returns the solution \var{u} of: operator * \var{u} = \var{rhs}. 1800 \end{methoddesc} 1801 1802 \begin{methoddesc}[Operator]{of}{u} 1803 applies the operator to the \Data object \var{u}, i.e. performs a matrix-vector 1804 multiplication. 1805 \end{methoddesc} 1806 1807 \begin{methoddesc}[Operator]{saveMM}{fileName}\index{Matrix Market} 1808 saves the object to a Matrix Market format file with name \var{fileName}, see 1809 \url{http://math.nist.gov/MatrixMarket} 1810 \end{methoddesc} 1811 1812 \section{Physical Units} 1813 \escript provides support for physical units in the SI system\index{SI units} 1814 including unit conversion. So the user can define variables in the form 1815 \begin{python} 1816 from esys.escript.unitsSI import * 1817 l=20*m 1818 w=30*kg 1819 w2=40*lb 1820 T=100*Celsius 1821 \end{python} 1822 In the two latter cases a conversion from pounds\index{pounds} and degrees 1823 Celsius\index{Celsius} is performed into the appropriate SI units \emph{kg} 1824 and \emph{Kelvin}. 1825 In addition, composed units can be used, for instance 1826 \begin{python} 1827 from esys.escript.unitsSI import * 1828 rho=40*lb/cm**3 1829 \end{python} 1830 defines the density in the units of pounds per cubic centimeter. 1831 The value$40$will be converted into SI units, in this case kg per cubic 1832 meter. Moreover unit prefixes are supported: 1833 \begin{python} 1834 from esys.escript.unitsSI import * 1835 p=40*Mega*Pa 1836 \end{python} 1837 The pressure \var{p} is set to 40 Mega Pascal. Units can also be converted 1838 back from the SI system into a desired unit, e.g. 1839 \begin{python} 1840 from esys.escript.unitsSI import * 1841 print(p/atm) 1842 \end{python} 1843 can be used print the pressure in units of atmosphere\index{atmosphere}. 1844 1845 The following is an incomplete list of supported physical units: 1846 1847 \begin{datadesc}{km} 1848 unit of kilometer 1849 \end{datadesc} 1850 1851 \begin{datadesc}{m} 1852 unit of meter 1853 \end{datadesc} 1854 1855 \begin{datadesc}{cm} 1856 unit of centimeter 1857 \end{datadesc} 1858 1859 \begin{datadesc}{mm} 1860 unit of millimeter 1861 \end{datadesc} 1862 1863 \begin{datadesc}{sec} 1864 unit of second 1865 \end{datadesc} 1866 1867 \begin{datadesc}{minute} 1868 unit of minute 1869 \end{datadesc} 1870 1871 \begin{datadesc}{h} 1872 unit of hour 1873 \end{datadesc} 1874 1875 \begin{datadesc}{day} 1876 unit of day 1877 \end{datadesc} 1878 1879 \begin{datadesc}{yr} 1880 unit of year 1881 \end{datadesc} 1882 1883 \begin{datadesc}{gram} 1884 unit of gram 1885 \end{datadesc} 1886 1887 \begin{datadesc}{kg} 1888 unit of kilogram 1889 \end{datadesc} 1890 1891 \begin{datadesc}{lb} 1892 unit of pound 1893 \end{datadesc} 1894 1895 \begin{datadesc}{ton} 1896 metric ton 1897 \end{datadesc} 1898 1899 \begin{datadesc}{A} 1900 unit of Ampere 1901 \end{datadesc} 1902 1903 \begin{datadesc}{Hz} 1904 unit of Hertz 1905 \end{datadesc} 1906 1907 \begin{datadesc}{N} 1908 unit of Newton 1909 \end{datadesc} 1910 1911 \begin{datadesc}{Pa} 1912 unit of Pascal 1913 \end{datadesc} 1914 1915 \begin{datadesc}{atm} 1916 unit of atmosphere 1917 \end{datadesc} 1918 1919 \begin{datadesc}{J} 1920 unit of Joule 1921 \end{datadesc} 1922 1923 \begin{datadesc}{W} 1924 unit of Watt 1925 \end{datadesc} 1926 1927 \begin{datadesc}{C} 1928 unit of Coulomb 1929 \end{datadesc} 1930 1931 \begin{datadesc}{V} 1932 unit of Volt 1933 \end{datadesc} 1934 1935 \begin{datadesc}{F} 1936 unit of Farad 1937 \end{datadesc} 1938 1939 \begin{datadesc}{Ohm} 1940 unit of Ohm 1941 \end{datadesc} 1942 1943 \begin{datadesc}{K} 1944 unit of degrees Kelvin 1945 \end{datadesc} 1946 1947 \begin{datadesc}{Celsius} 1948 unit of degrees Celsius 1949 \end{datadesc} 1950 1951 \begin{datadesc}{Fahrenheit} 1952 unit of degrees Fahrenheit 1953 \end{datadesc} 1954 1955 \noindent Supported unit prefixes: 1956 1957 \begin{datadesc}{Yotta} 1958 prefix yotta =$10^{24}$1959 \end{datadesc} 1960 1961 \begin{datadesc}{Zetta} 1962 prefix zetta =$10^{21}$1963 \end{datadesc} 1964 1965 \begin{datadesc}{Exa} 1966 prefix exa =$10^{18}$1967 \end{datadesc} 1968 1969 \begin{datadesc}{Peta} 1970 prefix peta =$10^{15}$1971 \end{datadesc} 1972 1973 \begin{datadesc}{Tera} 1974 prefix tera =$10^{12}$1975 \end{datadesc} 1976 1977 \begin{datadesc}{Giga} 1978 prefix giga =$10^9$1979 \end{datadesc} 1980 1981 \begin{datadesc}{Mega} 1982 prefix mega =$10^6$1983 \end{datadesc} 1984 1985 \begin{datadesc}{Kilo} 1986 prefix kilo =$10^3$1987 \end{datadesc} 1988 1989 \begin{datadesc}{Hecto} 1990 prefix hecto =$10^2$1991 \end{datadesc} 1992 1993 \begin{datadesc}{Deca} 1994 prefix deca =$10^1$1995 \end{datadesc} 1996 1997 \begin{datadesc}{Deci} 1998 prefix deci =$10^{-1}$1999 \end{datadesc} 2000 2001 \begin{datadesc}{Centi} 2002 prefix centi =$10^{-2}$2003 \end{datadesc} 2004 2005 \begin{datadesc}{Milli} 2006 prefix milli =$10^{-3}$2007 \end{datadesc} 2008 2009 \begin{datadesc}{Micro} 2010 prefix micro =$10^{-6}$2011 \end{datadesc} 2012 2013 \begin{datadesc}{Nano} 2014 prefix nano =$10^{-9}$2015 \end{datadesc} 2016 2017 \begin{datadesc}{Pico} 2018 prefix pico =$10^{-12}$2019 \end{datadesc} 2020 2021 \begin{datadesc}{Femto} 2022 prefix femto =$10^{-15}$2023 \end{datadesc} 2024 2025 \begin{datadesc}{Atto} 2026 prefix atto =$10^{-18}$2027 \end{datadesc} 2028 2029 \begin{datadesc}{Zepto} 2030 prefix zepto =$10^{-21}$2031 \end{datadesc} 2032 2033 \begin{datadesc}{Yocto} 2034 prefix yocto =$10^{-24}$2035 \end{datadesc} 2036 2037 \section{Utilities} 2038 The \class{FileWriter} class provides a mechanism to write data to a file. 2039 In essence, this class wraps the standard \PYTHON \class{file} class to write 2040 data that are global in \MPI to a file. In fact, data are written on the 2041 processor with \MPI rank 0 only. It is recommended to use \class{FileWriter} 2042 rather than \class{open} in order to write code that will run with and without 2043 \MPI. It is safe to use \class{open} under \MPI to \emph{read} data which are 2044 global under \MPI. 2045 2046 \begin{classdesc}{FileWriter}{fn\optional{,append=\False, \optional{createLocalFiles=\False}})} 2047 Opens a file with name \var{fn} for writing. If \var{append} is set to \True 2048 data are appended at the end of the file. 2049 If running under \MPI, only the first processor (rank==0) will open the file 2050 and write to it. 2051 If \var{createLocalFiles} is set each individual processor will create a file 2052 where for any processor with rank$> 0\$ the file name is extended by its rank. 2053 This option is normally used for debugging purposes only. 2054 \end{classdesc} 2055 2056 \vspace{1em}\noindent The following methods are available: 2057 \begin{methoddesc}[FileWriter]{close}{} 2058 closes the file. 2059 \end{methoddesc} 2060 \begin{methoddesc}[FileWriter]{flush}{} 2061 flushes the internal buffer to disk. 2062 \end{methoddesc} 2063 \begin{methoddesc}[FileWriter]{write}{txt} 2064 writes string \var{txt} to the file. Note that a newline is not added. 2065 \end{methoddesc} 2066 \begin{methoddesc}[FileWriter]{writelines}{txts} 2067 writes the list \var{txts} of strings to the file. 2068 Note that newlines are not added. 2069 This method is equivalent to calling \var{write()} for each string. 2070 \end{methoddesc} 2071 \begin{memberdesc}[FileWriter]{closed} 2072 this member is \True if the file is closed. 2073 \end{memberdesc} 2074 \begin{memberdesc}[FileWriter]{mode} 2075 holds the access mode. 2076 \end{memberdesc} 2077 \begin{memberdesc}[FileWriter]{name} 2078 holds the file name. 2079 \end{memberdesc} 2080 \begin{memberdesc}[FileWriter]{newlines} 2081 holds the line separator. 2082 \end{memberdesc} 2083 2084 \noindent The following additional functions are available in the \escript 2085 module: 2086 \begin{funcdesc}{setEscriptParamInt}{name,value} 2087 assigns the integer value \var{value} to the internal Escript parameter 2088 \var{name}. This should be considered an advanced feature and it is generally 2089 not required to call this function. One parameter worth mentioning is 2090 \var{name}="TOO_MANY_LINES" which affects the conversion of \Data objects to a 2091 string. If more than \var{value} lines would be created, a condensed format is 2092 used instead which reports the minimum and maximum values and general 2093 information about the \Data object rather than all values. 2094 \end{funcdesc} 2095 2096 \begin{funcdesc}{getEscriptParamInt}{name} 2097 returns the current value of internal Escript parameter \var{name}. 2098 \end{funcdesc} 2099 2100 \begin{funcdesc}{listEscriptParams}{a} 2101 returns a list of valid Escript parameters and their description. 2102 \end{funcdesc} 2103 2104 \begin{funcdesc}{getMPISizeWorld}{} 2105 returns the number of \MPI processes in use in the \env{MPI_COMM_WORLD} 2106 process group. If \MPI is not used 1 is returned. 2107 \end{funcdesc} 2108 2109 \begin{funcdesc}{getMPIRankWorld}{} 2110 returns the rank of the current process within the \env{MPI_COMM_WORLD} 2111 process group. If \MPI is not used 0 is returned. 2112 \end{funcdesc} 2113 2114 \begin{funcdesc}{MPIBarrierWorld}{} 2115 performs a barrier synchronization across all processes within the 2116 \env{MPI_COMM_WORLD} process group. 2117 \end{funcdesc} 2118 2119 \begin{funcdesc}{getMPIWorldMax}{a} 2120 returns the maximum value of the integer \var{a} across all processes within 2121 \env{MPI_COMM_WORLD}. 2122 \end{funcdesc} 2123 2124 \section{Lazy Evaluation of Data} 2125 \label{sec:lazy} 2126 Constant and Tagged representations of Data are relatively small but Expanded\footnote{Separate values stored for each point of the FunctionSpace.} are larger and 2127 will not entirely fit in CPU cache. 2128 2129 Escript's lazy evaluation features record operations performed on Data objects but do not actually carry them out until the Data is resolved''. 2130 2131 Consider the following code: 2132 \begin{python} 2133 from esys.escript import * 2134 from esys.dudley import Rectangle 2135 x=Rectangle(3,3) 2136 x=Rectangle(3,3).getX() 2137 c=Data((1.5, 1), x.getFunctionSpace()) 2138 t=Data(((1,1),(0,1)), x.getFunctionSpace()) 2139 t.tag() 2140 \end{python} 2141 2142 The variables \var{c}, \var{t}, \var{x} are stored as \texttt{constant}, \texttt{tagged} and \texttt{expanded} Data respectively. 2143 Printing those variables will show the values stored (or if we were to use a larger Rectangle, a summary). 2144 2145 \begin{python} 2146 v = matrix_mult(t,x) + c 2147 print(v.isExpanded()) 2148 print(v) 2149 \end{python} 2150 2151 Will output \texttt{True} followed by all of the values for \var{v}. 2152 Now we'll introduce lazy evaluation: 2153 2154 \begin{python} 2155 xx = x.delay() 2156 print(xx.isExpanded(), xx.isLazy()) 2157 print(x.isExpanded(), x.isLazy()) 2158 print(xx) 2159 \end{python} 2160 2161 The first print will show that \var{xx} is not considered to be expanded'', while the second print shows that \var{x} is unaffected. 2162 The last print will produce something like: 2163 \begin{python} 2164 Lazy Data: [depth=0] E@0x55ed512ad760 2165 \end{python} 2166 The \texttt{E} before the \verb|@| shows that this lazy Data is wrapping expanded'' Data. 2167 Calling \texttt{.delay()} on constant or tagged Data results in \verb|C@...| and \verb|T@...| respectively. 2168 2169 If an input to an operation is lazy, then the result will be lazy as well\footnote{Matrix inverse is an exception to this.}: 2170 \begin{python} 2171 res = matrix_mult(t,-xx) + c 2172 print(res) 2173 \end{python} 2174 Will produce: 2175 \begin{python} 2176 Lazy Data: [depth=3] (prod(T@0x..., neg(E@...)) + C@0x...) 2177 \end{python} 2178 Depth indicates the largest number of operators from the top of the expression to the bottom. 2179 2180 To actually find the value of this lazy Data object, we need to resolve it: 2181 \begin{python} 2182 res.resolve() 2183 \end{python} 2184 Note that \texttt{resolve()} doesn't return a new object, but transforms the object it is called on. 2185 Printing, \var{res} now will show the values at each point. 2186 2187 \subsection{Lazyness and non-expanded Data} 2188 While it is possible to call delay on constant or tagged Data, escript will not build expressions consisting solely of such Data. 2189 \begin{python} 2190 cx=c.delay() 2191 res=cx+cx 2192 print(res) 2193 \end{python} 2194 would output: 2195 \begin{python} 2196 Lazy Data: [depth=0] C@0x55ed512cc7c0 2197 # Not 2198 Lazy Data: [depth=1] (C@0x... + C@0x...) 2199 \end{python} 2200 2201 2202 \subsection{When to resolve} 2203 2204 You are never \emph{required} to manually resolve lazy Data in \texttt{escript}. 2205 Any operations which need the actual values of an expression will either 2206 \begin{itemize} 2207 \item compute the values without resolving the whole Data object at once (solvers assembling FEM matrices) 2208 \item resolve the data automatically (everthing else) 2209 \end{itemize} 2210 2211 \noindent Escript will automatically resolve lazy Data: 2212 \begin{enumerate} 2213 \item If a matrix inversion operation is applied to the Data. 2214 \item If the expression tree becomes too deep\footnote{At time of writing, this threshold is somewhat arbitrarily set at \texttt{depth>9}, but this is configurable.}. 2215 \end{enumerate} 2216 Note, the second point is important when writing loops like this: 2217 \begin{python} 2218 # x is initial guess 2219 while err > tol: 2220 construct PDE coefficients involving x 2221 solve PDE 2222 calculate err 2223 update x 2224 \end{python} 2225 2226 After a few iterations of the loop, \var{x} may be something like \texttt{x=F(F(F(F(originalX))))}. 2227 So it will probably be better to \texttt{resolve} \var{x} at the end of each loop iteration. 2228 Alternatively, if \var{x} is included in many expressions in the loop, it may be better to resolve it earlier. 2229 2230 \subsection{Options for using lazy evaluation} 2231 2232 There are two ways to enable lazy evaluation: 2233 \begin{enumerate} 2234 \item Any escript script can make use of lazy evaluation by \texttt{delay()}-ing one of its expanded Data variables. 2235 Any expressions including that delayed variable (directly or indirectly) will be lazy until resolved. 2236 \item Setting the \texttt{AUTOLAZY} parameter for \texttt{escript} to \texttt{1}. 2237 In this case, most escript operation which would normally produce extended Data, will produce lazy Data instead. 2238 In general, this option is not recommended for two reasons: 2239 \begin{itemize} 2240 \item AUTOLAZY uses the \texttt{setEscriptParamInt()} which is not guaranteed to have continued support. 2241 \item Making everything lazy instead of just more complex objects is not likely to give significant efficiency improvements. 2242 \end{itemize} 2243 \end{enumerate} 2244 2245 \subsection{When to use lazy evaluation?} 2246 Exactly when using lazy evaluation will be more efficient is still an open question. 2247 When the objects being manipulated are large (eg 4-Tensors in Drucker-Prager), significant memory and runtime improvements can be achieved. 2248 See~\cite{lazyauspdc}. 2249 2250 Our best advice is to experiment with it.

## Properties

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