Main function
This page has not been updated. However almost all
explanations are fully relevant.
track.cpp
The main source is in the file track.cpp.
The simulation is performed by calling an instance of the class Csimulator.
This is the last function called in the main. This simulation
(which uses templates) takes different arguments that are
instances of Multivac classes (e.g. SpeedType, MeshType...). These
instances (by the way, all the classes use templates) are to be
initialized, and that is precisely what is done in the first part
of the main. We have for example to choose the type of speed
function we want to use, the size of the mesh, which initial curve
we want, whether we want to save every curve or just the last one
and whether we want to use the regular Level set method, the
narrow band Level set method (this will be explained further) or
the Fast marching method...
The main is split in two parts:
- initialization of all variables and Multivac classes
- given this data, we launch the simulation
So each new set of values, written in track.cpp, will trigger a
new and different simulation (depending on these new
arguments).
From now on, the variables that can be modified to create a new
simulation will be followed, in the source file, by:
@***** simulation variable ******@.
The easiest way to proceed is to take the main and explain it in a
linear way, which is the natural "program order".
Note: A general explanation about the different types of
changes that have to be made by the user to initialize the
different class instances used by Multivac is given in DOMAIN MESH(which is used as an
example). This procedure, explained with DOMAIN MESH is quite the
same for SPEED FUNCTION, INITIAL CURVE, LEVEL SET FUNCTION,
INITIALIZER, UPDATER and SAVER.
Includes and options
Time
Domain mesh
Speed function
Initial curve
Level set function
Initializer
Updater
Saver
Simulation
INCLUDES AND OPTIONS
00001
00002
00003
00004
00005
00006
00007
00008 #define SELDON_DEBUG_LEVEL_1
00009 #define MULTIVAC_DEBUG_LEVEL_1
00010
00011
00012
00013 #define MULTIVAC_REPORT_NO
00014
00015 #define ERR(x) cout << "Hermes - " #x << endl
00016 #define DISPLAY(x) cout << #x ": " << x << endl
00017 #define DISP(x) cout << #x ": " << x << endl
00018
00019
00020 #include "seldon.hxx"
00021
00022
00023 #include "includes.hxx"
00024
00025 using namespace Multivac;
This is the header of the file track.cpp, before the definition
of the main function.
One notices the use of the namespace Multivac and the two
inclusions seldon.hxx and includes.hxx. These two header files
include functions and classes that will be used in the
main.
As for the namespace, Multivac uses two namespaces Multivac and
Seldon.
Actually, Seldon is another "little library" that defines
different useful tools such as a class vector, a class matrix,
a class list... which are used by Multivac.
There are also three macros defined: ERR(x), DISPLAY(x) and
DISP(x).
Finally, there are three other items "defined":
#define MULTIVAC_DEBUG_LEVEL_1
#define SELDON_DEBUG_LEVEL_1
#define MULTIVAC_REPORT_NO
#define MULTIVAC_DEBUG_LEVEL_1 & #define
SELDON_DEBUG_LEVEL_1
(MULTIVAC_DEBUG_LEVEL_1 implies SELDON_DEBUG_LEVEL_1 so here, it
is redundant)
Multivac and Seldon Libraries are provided with different levels
of error management (four for each):
MULTIVAC_DEBUG_LEVEL_1, MULTIVAC_DEBUG_LEVEL_2,
MULTIVAC_DEBUG_LEVEL_3, MULTIVAC_DEBUG_LEVEL_4
SELDON_DEBUG_LEVEL_1, SELDON_DEBUG_LEVEL_2,
SELDON_DEBUG_LEVEL_3, SELDON_DEBUG_LEVEL_4.
The user is then able to use 4 levels of errors handling (more
or less strict).
The level 1 is the least "strict" and the level 4 in the most
"complete" one. Each level has the same features than the
previous levels plus new ones. (level 2 is level 1 plus new
functions). In this example, for this track.cpp, the debug level
chosen is the number 1.
Let k be an integer between 1 and 4.
if MULTIVAC_DEBUG_LEVEL_k is defined, it automatically defines
SELDON_DEBUG_LEVEL_k.
MULTIVAC_DEBUG_LEVEL_k and SELDON_DEBUG_LEVEL_k handle the same
types of errors (for example, they may both take care of memory
management) but respectively for the MULTIVAC part and the
SELDON part of the program.
What are the features specific to each level?
Level 1
The exception handling mode is "off". No "throws" can be
used. All exception specifications are disabled.
Level 2
Memory management & Input/Output operations are
checked. (example: creating a new variable, reading/writing
a file...). If an error of one of these types occurs, then the
program catches the exception and processes a return (an error
message is also sent).
Level 3
Memory management & Input/Output operations are still handled
(cf Level 2)
Indices are also checked. Multivac & Seldon use matrices
and vectors and this Level enable the program to check and
report any problem of indices (negative or greater than the
size).
Level 4
This level checks everythingBut for the moment, it
does exactly the same thing as Level 3.
#define MULTIVAC_REPORT_NO, as mentioned in the comments
is to assess the duration of the updates and initializations and
display it on the screen.
#define MULTIVAC_REPORT could have been another possible
"define" here, it is used to monitor the whole process:
(initMesh, initphiandF, saveatcurrentiteration...). It records
any operation of Multivac and prints the corresponding times on
the screen.
Note about a programming technique: The definition of
the content of each level or of MULTIVAC_REPORT can be found in
the files includes.hxx and seldon/seldon.hxx. Then in the
other files, with a #ifdef ... #endif or a #ifndef ... #endif
structure, we can decide to tell the
compiler to include
that or that action or to disregard it depending on wether or not
an item has been "defined". That is why defining
MULTIVAC_REPORT_NO or MULTIVAC_REPORT for example, will not give
the same results. This use of the tools #define, #ifdef, #ifndef
... which are commands treated by the compiler (because of the
"#") is a C++ programming tip. When the compiler sees a #ifdef
NAME ... #endif for example, it will take into account the define
situation,i.e wether or not NAME has been previously defined by
the command #define NAME, if so, it will include the code between
#idef NAME and #endif in the compilation, ortherwise it won't. The
program can then have different features just by changing the
tools defined.
TIME
00028
00029
00030
00031
00032 int main()
00033 {
00034
00035 clock_t timer = clock();
00036
00037
00038 TRY
00039
00040
00041 typedef double real;
00042
00043
00045
00047
00048 real FinalTime = 0.1; @***** simulation variable ******@
00049 real Delta_t = 0.0001; @***** simulation variable ******@
00050
00051
00052 int NbIterations = int (FinalTime / Delta_t);
We can first notice the timer which will enable us to clock the
program, the typedef (double will be referred to as real) and
the TRY to catch the exceptions.
Note: In C++ the try
syntax is the following: try{...} catch(), but a lighter
notation: TRY ... END has been defined in Multivac in the file
includes.hxx.
Then FinalTime, which is the time of the simulation and Delta_t
the time step have to be set to a new value by the user
depending on the simulation he wishes to do. Delta_t is the step
used by the numerical schemes. The program will compute the
values of our functions (here, it is the function phi, our
"Level set function" which has the value 0 on the front, and
which is evaluated at each time step). The simulation stops when
the "FinalTime" is reached.
So the number of iterations is known when we have the time step
and the final time. The user doesn't have to change it.
DOMAIN MESH
00056
00058
00059
00060
00061
00062 typedef COrthogonalMesh<real> MeshType;
00063
00064
00065 real Xmin = 0.0; @***** simulation variable ******@
00066 real Xmax = 3.0; @***** simulation variable ******@
00067 real Ymin = 0.0; @***** simulation variable ******@
00068 real Ymax = 3.0; @***** simulation variable ******@
00069
00070
00071
00072 int Nx = 601; @***** simulation variable ******@
00073 int Ny = Nx; @***** simulation variable ******@
00074
00075
00076 MeshType Mesh(Xmin, Xmax, Ymin, Ymax, Nx, Ny);
This first part of the initialization aims at creating the
mesh. For the moment,
only one type of mesh can be used:
an orthogonal one. Xmin, Xmax, Ymin, Ymax which defines the
domain and Nx, Ny (number of grid points) have to be defined by
the user.
The structure "typedef COrthogonalMesh<real> MeshType;" is
used here so that MeshType will represent, in the rest of the
file, COrthogonalMesh<real>.
If another mesh type
called, let say, CTriangleMesh<real> had been available
for example, we would have written "typedef
CTriangleMesh<real> MeshType;" to actually use the
triangle mesh instead of "typedef COrthogonalMesh<real>
MeshType;". The use of this typedef is really convenient because
in the rest of the document, the mesh type class is referred to
as
MeshType (and stands for COrthogonalMesh or
CTriangleMesh in our example) and we only have to change one
line of the code to select which one we need and change the mesh
type.Then we create a variable Mesh of this class Meshtype
(MeshType Mesh(Xmin, Xmax, Ymin, Ymax, Nx, Ny)) with a set of
arguments. Indeed, the mesh is stored and saved in an instance,
an object of the class MeshType.
General Note: This technique will be used for each class
we will study further in this document. Each time a typedef is
used to associate a specific class to a generic name, name which
is used in the rest of the file.
Then we update a set of variables which are given as an argument
to a contructor of the class, to create an instance with the
properties needed.
In Multivac, most properties and datas are set in specific
"Multivac class instances". For example, the mesh is stored in
an instance of MeshType (we have seen previously that MeshType
refers to one of the mesh classes because there could be
different meshes possible, these different mesh classes are in
fact sub-classes of a mother mesh class: C++ heritage). The
files containing the definitions of these classes and their
sub-classes can be found in the directory
includes of
Multivac and have
explicit names to enable the user to
find his way in the program easily. Then for each mother
class, let say for example the speed function class (cf next
section), the user just has to select which version (i.e which
sub-class) he wants to use, in this precise case CConstantSpeed,
CSpeedFunction1 or CFireModel and write it in the typedef
(typedef sub-class Speedtype). Then he has to initialize the
different variables given as an argument to the instance of the
class created. These arguments can be different(most of the
time, their number is different) from one sub-class to
another.
Finally, depending on which sub-class is used, the user has to
give the constructor of the class the corresponding set of
arguments. In our example of the speed function, we have to
choose between:
1) F(SpeedRate) if the speed type is CConstantSpeed
2) F(SpeedRate, SpeedRate0, Limit) if the speed type is
CSpeedFunction1
3) F(U, a, b, epsilon0, epsilon1, c1, mu, CFL, coeff) if the
speed type is CFireModel
If the third one is chosen, "SpeedType F(U, a, b, epsilon0,
epsilon1, c1, mu, CFL, coeff)" creates an instance named "F" of
SpeedType (which stands for CFireModel here) with all the
paramaters it needs (9 here).
If the first one, a constant speed had been chosen we have
had:
typedef CConstantSpeed SpeedType;
initialization of SpeedRate
SpeedType F(SpeedRate);
Of course, in this case we just need one argument SpeedRate so
we don't have to change or write the others (unused). Each time,
the user just has to initialize the parameters needed (however
in this version of track.cpp, all parameters for every cases are
written to help the user, we then just need to modify the
parameters we want). Moreover, the different classes and
sub-classes available are written in the source code as
comments.
Note: We will discuss the nature of these "speed function"
parameters in the next section.
SPEED FUNCTION
00080
00082
00083
00084
00085
00086
00087 typedef CFireModel<real> SpeedType;
00088
00089
00090 real SpeedRate = 0.5; @***** simulation variable ******@
00091
00092
00093 real SpeedRate0 = 0.2; @***** simulation variable ******@
00094
00095
00096 real Limit = 0.9; @***** simulation variable ******@
00097
00098
00099
00100
00101
00102 real U = 100.0; @***** simulation variable ******@
00103 real a = 0.1; @***** simulation variable ******@
00104 real b = 1.0; @***** simulation variable ******@
00105 real epsilon0 = 0.02; @***** simulation variable ******@
00106 real epsilon1 = 0.003; @***** simulation variable ******@
00107 real c1 = 0.5; @***** simulation variable ******@
00108 real mu = 2.0; @***** simulation variable ******@
00109
00110 real CFL = FinalTime / NbIterations
00111 / min( (Xmax - Xmin) / (Nx-1), (Ymax - Ymin) / (Ny-1));
00112 real coeff = 2.0; @***** simulation variable ******@
00113 real MaxGrad = 0.5 * sqrt( pow(Xmax - Xmin, 2.0)
00114 + pow(Ymax - Ymin, 2.0) );
00115
00116
00117
00118
00119
00120 CFL = 0.9; @***** simulation variable ******@
00121 SpeedType F(U, a, b, epsilon0, epsilon1, c1, mu, CFL, coeff);
This part deals with the speed function, (the speed is known for
every point of the curve).
Programming note: One more time, the speed values (on the
grid) are saved in an instance of the speed class. But before
being saved, these values have to be assessed depending on which
speed model we want to use. This job is done thanks to heritage,
virtual functions and sub-classes of the speed class.
There are, as mentionned in the previous section, 3 possible
sub-classes defined for the moment:
- a constant speed :
This speed function, as its name indicates it, is a constant
speed function. The only parameter to initialize is
SpeedRate, the constant speed value.
- a speed function function1 :
We need to initialize 3 parameters: SpeedRate, SpeedRate0,
Limit.
function1 actually defines a piecewise constant speed
function. The Domain is split in two parts,
horizontally. Limit is the ordinate of the "cut". The speed
value is SpeedRate0 in the lower part (y < Limit) and
SpeedRate in the upper part (y > Limit).
- a fire model speed :
This speed function is more complicated. It is, as
mentionned in the subtitle a fire model speed. The speed
rate is a simple fire model, with wind.
Note: There are in fact 4 fire models defined, plus the
possibility to rotate the wind speed along with the time (the
wind speed is a vector, with a norm equal to the parameter U,
which has otherwise a constant direction and is then rotated of
an angle of 10 x current time). As seen in the previous section,
the use of #ifdef enable us here to define several models and
select just one during the compilation using the key word
#define NAME where NAME belongs to:
MULTIVAC_FIREMODEL4, MULTIVAC_FIREMODEL3, MULTIVAC_FIREMODEL2 or
MULTIVAC_FIREMODEL1
The same thing can be applied to MULTIVAC_ROTATING_WIND0
(constant wind direction), MULTIVAC_ROTATING_WIND (rotating
wind).
But this time, in this version of the code, #define
MULTIVAC_FIREMODEL4 and #define MULTIVAC_ROTATING_WIND0 are
called in the head of the file
/userpath/includes/speedfunctions/firemodel.cxx. So if we want
to try to change it, we have to modify that in this file. These
different models are defined in the end of the file
/userpath/includes/speedfunctions/firemodel.cxx.
More information can be found on this subject on the website of
Vivien Mallet:
http://spacetown.free.fr/fronts/
For any point M, let theta be the angle between the direction of
the wind and the normal of the curve at this point M. Let cos2 =
cos(theta)*cos(theta) and sin2 =
sin(theta)*sin(theta).
Speed function : F
MULTIVAC_FIREMODEL1
F(theta) = epsilon0 * sin2 + a * U * sin2 Exp(-b * U * sin2) +
G(theta)
where G(theta) is defined as:
If |theta|<pi/2 then G(theta) = epsilon0 * cos2 + c1 *
srqt(U) * costheta^(5/2)
else G(theta) = epsilon0 * cos2 * Exp(-epsilon1 * U *
cos2)
MULTIVAC_FIREMODEL2
F(theta) = epsilon0 * sin2 + a * U * sin2 Exp(-b * U * sin2) +
G(theta)
where G(theta) is defined as:
If |theta|<pi/2 then G(theta) = epsilon0 * cos2 + c1 *
srqt(U) * costheta^(3/2)
else G(theta) = epsilon0 * cos2 * Exp(-epsilon1 * U *
cos2)
MULTIVAC_FIREMODEL3
F(theta) = epsilon0 * sin2 + a * U * sin2 Exp(-b * U * sin2) +
G(theta)
where G(theta) is defined as:
If |theta|<pi/2 then G(theta) = epsilon0 * cos2 + c1 *
srqt(U) * costheta^(5/2)
else G(theta) = epsilon0 * cos2 * Exp(-epsilon1 * U *
cos2)
MULTIVAC_FIREMODEL4
Model used here.
if (costheta > 0)
F(theta) = epsilon0 + c1 * sqrt(U) * pow(costheta, 1.5) + a * U
* sin2 * exp(-b * U * sin2) );
else
F(theta) = epsilon0 * sin2 + a * U * sin2 * exp(-b * U * sin2)+
epsilon0 * cos2 * exp(-epsilon1 * U * cos2) );
These models calculate the speed knowing the wind speed, a few
parameters and the normal of the curve.
- U magnitude of wind velocity.
- a parameter (beta).
- b parameter (beta).
- epsilon0 parameter (v_f, epsilon, beta).
- epsilon1 parameter (epsilon).
- c1 parameter (v_f).
- mu parameter (beta).
- CFL CFL number. The CFL is the Courant-Lax-Friedrichs
condition.
- coeff correction to the CFL number.
With:
- v_f(x) = epsilon0 + c1 * sqrt(x). (head)
- beta(x) = epsilon0 + a * x * exp(- b * x). (flank)
- epsilon(x) = epsilon0 * exp(- epsilon1 * x). (rear)
Note:
In the previous versions of Multivac, the speed rate was defined
as followed:
If |theta| < pi / 2
F(U, theta) = vm_f(U cos(theta)) + beta(U sin^mu(theta))
else
F(U, theta) = beta(U sin^mu(theta)) + epsilon(U
cos^2(theta)).
This definition is written in the function "Value_Normal" in the
file "userpath/includes/speedfunction/firemodel.cxx". But this
function is not called anymore in the new versions, instead, we
use the 4 models previously defined (which are written in the
function: "Model" in the same file).
That's why in the new version, the parameter mu is not used
anymore (it is included indirectly in the different
models).
INITIAL CURVE
00125
00127
00128
00129
00130
00131
00132
00133
00134
00135 typedef CCircle<real> InitialCurveType;
00136
00137
00138 real CircleCenterX = 1.0; @***** simulation variable ******@
00139 real CircleCenterY = 1.5; @***** simulation variable ******@
00140 real CircleRadius = 0.5; @***** simulation variable ******@
00141
00142
00143 real CircleCenterX0 = 1.65; @***** simulation variable ******@
00144 real CircleCenterY0 = 1.6; @***** simulation variable ******@
00145 real CircleRadius0 = 0.3; @***** simulation variable ******@
00146
00147
00148 real CircleCenterX1 = 0.50; @***** simulation variable ******@
00149 real CircleCenterY1 = 1.0; @***** simulation variable ******@
00150 real CircleRadius1 = 0.25; @***** simulation variable ******@
00151
00152
00153
00154
00155 string InitialFrontFile = "F:/Fronts/multivac/includes/initialcurves/M.pts";
@***** simulation variable ******@
00156
00157
00158
00159
00160
00161
00162
00163
00164
00165
00166
00167 InitialCurveType InitialCurve(CircleCenterX, CircleCenterY, CircleRadius);
The notations are explicit. This class creates the first curve
(from a set of datas).
- CCircle creates on circle.
- CTwoCircles creates two circles.
- CThreeCircles creates three circles.
- CIsland creates two circles: the second circle is inside
the first one. The second circle is called "island".
- CIsland0: the initial front consists of three circles. The
third one is inside one of the two previous one (it is also
called "island").
- CSetOfPoints enable the user to write a set of points in a
file and create a curve with these points. The only thing to
do is to give CsetOfPoint contructor the complete name of the
file (with the full path). We then have to replace
F:/Fronts/multivac/includes/initialcurves/M.pts by our new
file. In the file, we have to give the coordinates of each
point P (one line per point) in this format: XP YP . The
points have to be given "sorted" in either trigonometrical or
clockwise way.
Example: The file M.pts creates the shape of the letter
M.
M.pts
0.51 0.2
0.5 0.21
0.5 2.79
0.51 2.8
0.9 2.8
1.49 2.0
1.51 2.0
2.1 2.8
2.49 2.8
2.5 2.79
2.5 0.21
2.49 0.2
2.11 0.2
2.1 0.21
2.1 1.6
2.09 1.61
1.51 1.2
1.49 1.2
0.91 1.61
0.9 1.6
0.9 0.21
0.89 0.2
Initial M curve:
LEVEL SET FUNCTION
00171
00173
00174
00175
00176
00177 typedef COrthogonalLevelSet<real> LevelSetType;
00178
00179
00180
00181 LevelSetType Phi;
We only have one choice, we have to create phi on an orthogonal
mesh because we have seen previously that it is the only type of
mesh for the moment. By the way, creating a 3D simulation or
using a mesh using triangles could be future improvements.
So, there is nothing to be done with this class.
INITIALIZER
00185
00187
00188
00189
00190
00191
00192
00193
00194
00195
00196 typedef CNarrowBandNeverInit<real> InitializerType;
00197
00198 InitializerType Initializer;
From now on, a new feature is now introduced: the method
used to realize the simulation.
In the previous initializations, it had no inluence.
The Level Set Method and the Fast Marching Method are both
implemented under Multivac.
By the way, the Level Set Method is implemented "twice". In
fact, a more efficient version of the Level Set Method has been
implemented as well: this is the Narrow Band Level Set
Method. Phi is computed within a narrow band instead of being
assessed on the whole grid. That is why this method, developped
for programming, is better than the regular algorithm. In
Multivac, the term NarrowBand is used in the name to refer to
this method. We also refer to the Fast Marching Method by using
FastMarching in the names as well.
More information about these different techniques can be found
in the book of J.A Sethian
Level Set Methods and Fast
Marching Methods Cambridge University press.
In the initializer, if we want to use the regular Level Set
Method or the Fast Marching Method we have to use respectively
the first (CNeverInit) and the fourth class
(CFastMarchingNeverInit).
These initializers only provide first initializations. It
doesn't provide any reinitialization or update.
The mesh is initialized. The level set is initialized according
to the initial curve position. Then, the speed function is
initialized according to the initial level set function. The
reinitialization only calls the speed function
updater.
But if we use the narrow band version, then we have two choices:
"never init" and "init domain".
This initializer can be used to perform initializations and
reinitializations of the level set function and the tube
involved in the narrow band level set method.
First, the mesh is initialized. The level set is initialized
according to the initial curve position. Then, the speed
function is initialized according to the initial level set
function. Reinitializations build the tube and init the
corresponding level set function. As for the speed function,
reinitializations only call the speed function updater.
In the narrow band method, the tube has to be recalculated any
time a point of the front enters a limit zone defined in the
updater (cf next section). In fact, with this method, the front
evolves in a certain tube and sometimes tends to go out of
it. In such situations, a new tube is recomputed with its
associated phi (this is handled by the initializer).
The only difference between "never init" and "init domain" is
that "init domain", reinitializes a new domain (newXmin,
newXmax, newYmin, newYmax) from the initial domain (Xmin, Xmax,
Ymin, Ymax) anytime to fit the current situation. Indeed, we
only work with a narrow band, so shrinking the domain to a
smaller one (of course still containing the curve and the band)
and using the same number of grid points (Nx x Ny) is more
accurate since the new space steps (newDelta_x & newDelta_y) are
smaller.
For "never init", the initial domain and the initial grid
pattern are kept.
UPDATER
00202
00204
00205
00206
00207
00208
00209
00210
00211
00212
00213
00214
00215
00216 typedef CNarrowBandFirstOrderEngquistOsher<real> UpdaterType;
00217
00218
00219
00220
00221
00222
00223
00224
00225
00226
00227
00228
00229
00230 UpdaterType Updater(6, 3, 1); @***** simulation variable ******@
The updater is the class that handles the numerical scheme and
"updates" phi (i.e: evaluates phi at the time t+dt knowing phi
at the time t). It also handles the tube in the Narrow Band
version (checks wether the tube's "limit barrier" is reached or
not) and the binary tree in the Fast Marching Method.
We have to choose our numerical scheme between EngquistOsher,
LaxFriedrichs and Godunov. These are all first order
schemes. Ones again, if we use the regular Level Set Method we
have to choose among the three first: 1) 2) or 3), if we use the
Narrow Band version we can use 4) 5) or 6) and we only have one
possibilty for the Fast Marching Method 7).
As for the arguments, the regular Level Set Method takes none
and the Fast Marching Method only needs "a time greater than all
arrival times that will be computed".
But for the Narrow Band Algorithm, three arguments are required:
TubeSemiWidth, BarrierWidth, OutSpaceWidth.As mentionned as
a comment in the code:
- TubeSemiWidth: number of grid points on each side
of the front. (zone 1)
- BarrierWidth: number of grid points (on each side
of the front) that imply tube reconstruction when reached.
(zone 2)
- OutSpaceWidth: number of grid points (on each side
of the front) that must not be reached. (zone 3)
Here is a scheme of the situation: (the bold line is the
front)
The user has to choose these values carefully. If the tube is
too large, the computation time increases. But if it is too
small, then it has to be reinitialized a lot, because the curve
reaches the second zone too often, and the cost of such
operation is quite high.
In this reason the choice of these parameters is an important
issue, as far as computational time is concerned. A
reconstruction of the tube every 6 to 10 iterations seems to be
a fine compromise and gives good results. That's why the
parameters should be chosen (knowing the grid and having an idea
of the speed) so that this range (6-10 iterations) could be
possible. It would give the best results.
SAVER
00234
00236
00237
00238
00239
00240
00241
00242
00243
00244
00245
00246
00247
00248
00249
00250
00251 typedef CTrack<real> SaverType;
00252
00253
00254
00255 int NbCurves = 10; @***** simulation variable ******@
00256
00257
00258 char flag = 'f';
00259
00260 int width = 10; @***** simulation variable ******@
00261
00262 int precision = 10; @***** simulation variable ******@
00263
00264
00265 string Directory = "/home/vivien/Computations/Fronts/"; @***** simulation variable ******@
00266
00267
00268 string TimeFile = Directory + "Time";
00269
00270
00271 string CurvesFile = Directory + "Curves";
00272
00273 string CurveLengthsFile = Directory + "CurveLengths";
00274
00275 string PhiFile = Directory + "Phi";
00276
00277 string FFile = Directory + "F";
00278
00279
00280 string XFile = Directory + "X";
00281
00282 string YFile = Directory + "Y";
00283
00284 string PointsFile = Directory + "Points";
00285
00286 string EdgesFile = Directory + "Edges";
00287
00288 string TrianglesFile = Directory + "Triangles";
00289
00290
00291 int Period = NbIterations / (NbCurves==0?1:NbCurves);
00292
00293 SaverType Saver(flag, width, precision, TimeFile, CurvesFile,
00294 CurveLengthsFile, PhiFile, FFile, XFile, YFile,
00295 PointsFile, EdgesFile, TrianglesFile, Period);
As we work with the level set function phi, the front is not
always reconstructed and if we want to save the front we
sometimes have to build it first. That is what the Saver class
does. We can first notice that there is a variable called
NbCurves that enable the user to choose how many curves he wants
to save. A value 0 means that all the curves will be
saved. Otherwise, if we choose to save let say 4 curves, it
won't be the first 4 curves, but a curve will be chosen every
period (the period is defined as the total number of iterations
divided by the number of curves to save). string Directory =
"/home/vivien/Computations/Fronts/" has also to be changed and
set by the user to the output directory he wants. Example: In
the "organisation tip section", the output advised is:
string Directory = "userpath/runs/runk/output/".
- CNeverSave: This saver doesn't do anything, it is
used to clock the calculations.
- CCurvesSaver: This saver saves the curves every
period. It also saves the different iteration times (of the
whole simulation), the initial mesh, the last level set
function and the last speed function.
- CSaveLastCurve: Saves the last mesh, the last speed
function, the last level set function and the last curve.
- CSaveAtTheEnd: Saves the last mesh, the last speed
function and the last level set function.
- CTrack: It is used with the optimization and
doesn't save anything.
The user just has to set the path in which he wants the output
files to be stored. These output files have explicit names and
contain the different values saved (curves, phi,
mesh...).
- flag: type of values to be saved (example: 'f'
prints a decimal floating point output -- cf 'sprintf'
argument 'format'). Here, "f" has to be used.
- width: minimum number of characters to be printed.
- precision: number of digits to be printed after de
decimal point (for
floating types) or minimum number of decimal digits to be
printed (for integers).
- TimeFileName: file in which iteration times will be
stored.
- CurvesFileName: file in which fronts will be
stored.
- CurveLengthsFileName: file in which numbers of
points (on fronts) will be stored.
- PhiFileName: file in which the level set function
will be stored.
- FFileName: file in which the speed function will be
stored.
- XFileName: file in which grid point abscissae will
be stored.
- YFileName: file in which grid point ordinates will
be stored.
- PointsFileName: file in which mesh points will be
stored.
- EdgesFileName: file in which mesh edges will be
stored.
- TrianglesFileName: file in which mesh triangles
will be stored.
- Period: Data will be saved if the current iteration
is a multiple of the period.
SIMULATION
To spawn a simulation, an instance called Simulator of the class
CSimulator, taking as parameters the different instances of
every Multivac class previously initialized (MeshType,
SpeedType, Phi...), is created.
Then, Simulator is initialized using the member function Init:
"Simulator.Init()" and the simulation is lauched with the member
function run: "Simulator.Run()".
Note about the optimization:
Multivac has also an optimization class that we can use instead
of the simulation.
The optimization part is still being improved. More informations
about the algorithm chosen to compute the optimization can be
found on the following website:
http://spacetown.free.fr/fronts/.
Both optimization and simulation are called with the same set of
parameters. Everything we have seen in this document is
unchanged when we use the optimization function.
Coptimizer2nd is the name of the class chosen if we want an
optimization.
So, the user just has to select wether he wants to do a
simulation or an optimization by creating the instance Simulator
from the right class.
We have only two choices.
Note: the structure is the following: the name of the
class (e.g. Csimulator), a list of templates used by the class
(between <...>), the name of the instance created
(Simulator) and the list of arguments.
- CSimulator<real, MeshType, SpeedType,
InitialCurveType,LevelSetType, InitializerType, UpdaterType,
SaverType> Simulator(Mesh, F, InitialCurve, Phi,
Initializer, Updater, Saver, NbIterations,
FinalTime);
- COptimizer2nd<real, MeshType, SpeedType,
InitialCurveType,LevelSetType, InitializerType, UpdaterType,
SaverType> Simulator(Mesh, F, InitialCurve, Phi,
Initializer, Updater, Saver, NbIterations, FinalTime);