heat.tex 17.9 KB
Newer Older
Thomas Witkowski's avatar
Thomas Witkowski committed
1
2
3
\section{Time dependent problem}
\label{s:heat}

4
5
This is an example for a time dependent scalar problem. The problem is
described by the heat equation
Thomas Witkowski's avatar
Thomas Witkowski committed
6
7
8
9
10
11
12
\begin{eqnarray}
\label{eq:heat}
\partial_t u - \Delta u &=& f ~~~~~~~ \mbox{in } \Omega \subset \mathbb{R}^{dim}\times\left({t^{begin}},t^{end}\right)\\
u &=& g ~~~~~~~ \mbox{on } \partial\Omega \times\left({t^{begin}},t^{end}\right)\\
u &=& u_0 ~~~~~ \mbox{on } \Omega \times \left({t^{begin}}\right).
\end{eqnarray}

13
14
15
16
We solve the problem in the time interval
$\left(t^{begin},t^{end}\right)$ with Dirichlet boundary conditions on
$\partial\Omega$. The problem is constructed, such that the exact
solution is $u(x,t) = \mbox{sin}(\pi t)e^{-10x^2}$. So we set
Thomas Witkowski's avatar
Thomas Witkowski committed
17
18
19
20
21
22
\begin{eqnarray}
f(x,t) &=& \pi \mbox{cos}(\pi t)e^{-10x^2}-\left(400x^2-20dow\right)\mbox{sin}(\pi t)e^{-10x^2}\\
g(x,t) &=& \mbox{sin}(\pi t)e^{-10x^2} \\
u_0(x) &=& \mbox{sin}(\pi t^{begin})e^{-10x^2}.
\end{eqnarray}

23
24
We use a variable time discretization scheme. Equation (\ref{eq:heat})
is approximated by
Thomas Witkowski's avatar
Thomas Witkowski committed
25
26
27
\begin{equation}
\frac{u^{new} - u^{old}}{\tau} - \left( \theta\Delta u^{new} + (1-\theta)\Delta u^{old} \right) = f(\cdot, t^{old}+\theta\tau).
\end{equation}
28
29
30
31
32
33
34
35
36
$\tau = t^{new} - t^{old}$ is the timestep size between the old and
the new problem time. $u^{new}$ is the (searched) solution at
$t=t^{new}$. $u^{old}$ is the solution at $t=t^{old}$, which is
already known from the last timestep. The parameter $\theta$
determines the implicit and explicit treatment of $\Delta u$. For
$\theta = 0$ we have the forward explicit Euler scheme, for $\theta =
1$ the backward implicit Euler scheme. $\theta = 0.5$ results in the
Crank-Nicholson scheme. If we bring all terms that depend on $u^{old}$
to the right hand side, the equation reads
Thomas Witkowski's avatar
Thomas Witkowski committed
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
\begin{equation}
\frac{u^{new}}{\tau} - \theta\Delta u^{new} = \frac{u^{old}}{\tau}  + (1-\theta)\Delta u^{old} + f(\cdot, t^{old}+\theta\tau).
\end{equation}

%\begin{figure}[h]
%\center
%\includegraphics[width=0.5\textwidth]{fig/theta.pdf}
%\caption{Time discretization with the $\theta$ scheme.}
%\label{f:theta}
%\end{figure}

\begin{table}
\center
\begin{tabular}{|cc|c|c|c|}
\hline 
& & {\bf 1d} & {\bf 2d} & {\bf 3d} \\
\hline 
{\bf source code} & \tt src/ & \multicolumn{3}{|c|}{\tt heat.cc} \\
\hline 
{\bf parameter file} & \tt init/ & \tt heat.dat.1d & \tt heat.dat.2d & \tt heat.dat.3d \\
\hline 
{\bf macro file} & \tt macro/ & \tt macro.stand.1d & \tt macro.stand.2d & \tt macro.stand.3d  \\
\hline 
{\bf output files} & \tt output/ & \multicolumn{3}{|c|}{\tt heat\_<t>.mesh, heat\_<t>.dat} \\
\hline 
\end{tabular}
\caption{Files of the {\tt heat} example. In the output file names, {\tt <t>} is replaced by the time.}
\end{table}

\subsection{Source code}
Now, we describe the crucial parts of the source code. First, the functions $f$ and $g$ are defined. In contrast to the ellipt example, the functions now are time dependent. This is implemented by deriving the function classes also from class \verb+TimedObject+. This class provides a pointer to the current time, as well as corresponding setting and getting methods. The usage of a pointer to a real value allows to manage the current time in one location. All objects that deal with the same time, point to the same value. In our example, $f$ is evaluated at $t=t^{old}+\theta\tau$, while $g$ (the Dirichlet boundary function for $u^{new}$) is evaluated at $t=t^{new}$. Function $g$ is implemented as follows:

\begin{lstlisting}{}
class G : public AbstractFunction<double, WorldVector<double> >,
          public TimedObject
{
public:
  double operator()(const WorldVector<double>& x) const
  {
    return sin(M_PI * (*timePtr)) * exp(-10.0 *(x * x));
  }
};
\end{lstlisting}

The variable $timePtr$ is a base class member of \verb+TimedObject+. This pointer has to be set once before $g$ is evaluated the first time. Implementation of function $f$ is done in the same way.

\begin{figure}
\center
\includegraphics[width=0.5\textwidth]{fig/heat_uml.pdf}
\caption{UML diagram for class Heat.}
\label{f:heat_uml}
\end{figure} 
Now, we begin with the implementation of class \verb+Heat+, that represents the instationary problem. In Figure \ref{f:heat_uml}, its class diagram is shown. \verb+Heat+ is derived from class \verb+ProblemInstatScal+ which leads to following properties:
\begin{description}
\item \verb+Heat+ implements the \verb+ProblemTimeInterface+, so the adaptation loop can set the current time and schedule timesteps. 
\item \verb+Heat+ implements \verb+ProblemStatBase+ in the role as initial (stationary) problem. The adaptation loop can compute the initial solution through this interface. The single iteration steps can be overloaded by sub classes of \verb+ProblemInstatScal+. Actually, the initial solution is computed through the method \verb+solveInitialProblem+ of \verb+ProblemTimeInterface+. But this method is implemented by \verb+ProblemInstatScal+ interpreting itself as initial stationary problem.  
\item \verb+Heat+ knows another implementation of \verb+ProblemStatBase+: This other implementation represents a stationary problem which is solved within each timestep.\end{description}

The first lines of class \verb+Heat+ are:

\begin{lstlisting}{}
class Heat : public ProblemInstatScal
{
public:
Thomas Witkowski's avatar
Thomas Witkowski committed
101
  Heat(ProblemScal &heatSpace)
Thomas Witkowski's avatar
Thomas Witkowski committed
102
103
104
105
106
    : ProblemInstatScal("heat", heatSpace)
  {
    theta = -1.0;
    GET_PARAMETER(0, name + "->theta", "%f", &theta);
    TEST_EXIT(theta >= 0)("theta not set!\n");
107
    theta1 = theta - 1.0;
108
  }
Thomas Witkowski's avatar
Thomas Witkowski committed
109
110
111
112
113
114
115
\end{lstlisting}

The argument \verb+heatSpace+ is a pointer to the stationary problem which is solved each timestep. It is directly handed to the base class constructor of \verb+ProblemInstatScal+. In the body of the constructor, $\theta$ is read from the parameter file and stored in a member variable. The member variable \verb+theta1+ stores the value of $\theta - 1$. A pointer to this value is used later as factor in the $\theta$-scheme. 

The next lines show the implementation of the time interface. 

\begin{lstlisting}{}
116
117
  void setTime(AdaptInfo *adaptInfo) 
  {
Thomas Witkowski's avatar
Thomas Witkowski committed
118
    rhsTime = adaptInfo->getTime() - (1 - theta) * adaptInfo->getTimestep();
Thomas Witkowski's avatar
Thomas Witkowski committed
119
    boundaryTime = adaptInfo->getTime();
120
  }
Thomas Witkowski's avatar
Thomas Witkowski committed
121

122
123
  void closeTimestep(AdaptInfo *adaptInfo) 
  {
Thomas Witkowski's avatar
Thomas Witkowski committed
124
125
    ProblemInstatScal::closeTimestep(adaptInfo);
    WAIT;
126
  }
Thomas Witkowski's avatar
Thomas Witkowski committed
127
128
\end{lstlisting}

129
130
131
The method \verb+setTime+ is called by the adaptation loop to inform
the problem about the current time. The right hand side function $f$
will be evaluated at $t^{old}+\theta\tau = t^{new} - (1-\theta)\tau$,
132
133
the Dirichlet boundary function $g$ at $t^{new}$, which is the
current time.
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150

The method \verb+closeTimestep+ is called at the end of each timestep
by the adaptation loop. In the default implementation of
\verb+ProblemInstatScal::closeTimestep+, the solution is written to
output files, if specified in the parameter file. Note that the base
class implementation of a method must be explicitly called, if the
method is overwritten in a sub class. The macro \verb+WAIT+ waits
until the \verb+return+ key is pressed by the user, if the
corresponding entry in the parameter file is set to 1. The macro
\verb+WAIT_REALLY+ would wait, independent of parameter settings. If
\verb+closeTimestep+ wouldn't be overloaded here, the default
implementation without the \verb+WAIT+ statement would be called after
each timestep.

Now, the implementation of the \verb+ProblemStatBase+ interface
begins. As mentioned above, the instationary problem plays the role of
the initial problem by implementing this interface.
Thomas Witkowski's avatar
Thomas Witkowski committed
151
152
153
154
155

\begin{lstlisting}{}
  void solve(AdaptInfo *adaptInfo) 
  {
    problemStat->getSolution()->interpol(exactSolution);
156
  }
Thomas Witkowski's avatar
Thomas Witkowski committed
157
158
159
160
161
162
163
164
165

  void estimate(AdaptInfo *adaptInfo) 
  {
    double errMax, errSum;
    errSum = Error<double>::L2Err(*exactSolution,
                                  *(problemStat->getSolution()), 
                                  0, &errMax, false);
    adaptInfo->setEstSum(errSum, 0);
    adaptInfo->setEstMax(errMax, 0);
166
  }
Thomas Witkowski's avatar
Thomas Witkowski committed
167
168
\end{lstlisting}

169
170
171
172
173
174
175
176
177
178
Here, only the solve and the estimate step are overloaded. For the
other steps, there are empty default implementations in
\verb+ProblemInstatScal+. Since the mesh is not adapted in the initial
problem, the initial adaptation loop will stop after one iteration. In
the solve step, the exact solution is interpolated on the macro mesh
and stored in the solution vector of the stationary problem. In the
estimate step, the L2 error is computed. The maximal element error and
the sum over all element errors are stored in \verb+adaptInfo+. To
make the exact solution known to the problem, we need a setting
function:
Thomas Witkowski's avatar
Thomas Witkowski committed
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203

\begin{lstlisting}{}
  void setExactSolution(AbstractFunction<double, WorldVector<double> > *fct) 
  {
    exactSolution = fct;
  } 
\end{lstlisting}

Now, we define some getting functions and the private member variables:

\begin{lstlisting}{}
  double *getThetaPtr() { return &theta; };
  double *getTheta1Ptr() { return &theta1; };
  double *getRHSTimePtr() { return &rhsTime; };
  double *getBoundaryTimePtr() { return &boundaryTime; };

private:
  double theta;
  double theta1;
  double rhsTime;
  double boundaryTime;
  AbstractFunction<double, WorldVector<double> > *exactSolution;
};
\end{lstlisting}

204
205
The definition of class \verb+Heat+ is now finished. In the following,
the main program is described.
Thomas Witkowski's avatar
Thomas Witkowski committed
206
207
208
209
210
211
212
213
214
215
216

\begin{lstlisting}{}
int main(int argc, char** argv)
{
  // ===== check for init file =====
  TEST_EXIT(argc == 2)("usage: heat initfile\n");

  // ===== init parameters =====
  Parameters::init(false, argv[1]);

  // ===== create and init stationary problem =====
Thomas Witkowski's avatar
Thomas Witkowski committed
217
218
  ProblemScal heatSpace("heat->space");
  heatSpace.initialize(INIT_ALL);
Thomas Witkowski's avatar
Thomas Witkowski committed
219
220

  // ===== create instationary problem =====
Thomas Witkowski's avatar
Thomas Witkowski committed
221
222
  Heat heat(heatSpace);
  heat.initialize(INIT_ALL);
Thomas Witkowski's avatar
Thomas Witkowski committed
223
224
\end{lstlisting}

225
226
227
228
229
230
So far, the stationary space problem \verb+heatSpace+ and the
instationary problem \verb+heat+ were created and
initialized. \verb+heatSpace+ is an instance of
\verb+ProblemScal+. \verb+heat+ is an instance of the class
\verb+Heat+ we defined above. \verb+heatSpace+ is given to \verb+heat+
as its stationary problem.
Thomas Witkowski's avatar
Thomas Witkowski committed
231

232
233
The next step is the creation of the needed \verb+AdaptInfo+ objects
and of the instationary adaptation loop:
Thomas Witkowski's avatar
Thomas Witkowski committed
234
235
236

\begin{lstlisting}{}
  // create adapt info for heat
Thomas Witkowski's avatar
Thomas Witkowski committed
237
  AdaptInfo adaptInfo("heat->adapt");
Thomas Witkowski's avatar
Thomas Witkowski committed
238
239

  // create initial adapt info
Thomas Witkowski's avatar
Thomas Witkowski committed
240
  AdaptInfo adaptInfoInitial("heat->initial->adapt");
Thomas Witkowski's avatar
Thomas Witkowski committed
241
242

  // create instationary adapt
Thomas Witkowski's avatar
Thomas Witkowski committed
243
244
245
246
247
  AdaptInstationary adaptInstat("heat->adapt",
				heatSpace,
				adaptInfo,
				heat,
				adaptInfoInitial);
Thomas Witkowski's avatar
Thomas Witkowski committed
248
249
\end{lstlisting}

250
251
252
253
254
The object \verb+heatSpace+ is handed as
\verb+ProblemIterationInterface+ (implemented by class\\
\verb+ProblemScal+) to the adaptation loop. \verb+heat+ is interpreted
as \verb+ProblemTimeInterface+ (implemented by class
\verb+ProblemInstatScal+).
Thomas Witkowski's avatar
Thomas Witkowski committed
255

256
The functions $f$ and $g$ are declared in the following way:
Thomas Witkowski's avatar
Thomas Witkowski committed
257
258
\begin{lstlisting}{}
  // ===== create boundary functions =====
259
  G *boundaryFct = new G;
Thomas Witkowski's avatar
Thomas Witkowski committed
260
261
262
  boundaryFct->setTimePtr(heat.getBoundaryTimePtr());
  heat.setExactSolution(boundaryFct);
  heatSpace.addDirichletBC(1, boundaryFct);
Thomas Witkowski's avatar
Thomas Witkowski committed
263
264

  // ===== create rhs functions =====
Thomas Witkowski's avatar
Thomas Witkowski committed
265
  int degree = heatSpace.getFESpace()->getBasisFcts()->getDegree();
266
  F *rhsFct = new F(degree);
Thomas Witkowski's avatar
Thomas Witkowski committed
267
  rhsFct->setTimePtr(heat.getRHSTimePtr());
Thomas Witkowski's avatar
Thomas Witkowski committed
268
\end{lstlisting}
269
270
271
272
The functions interpreted as \verb+TimedObject+s are linked with the
corresponding time pointers by \verb+setTimePtr+. The boundary
function is handed to \verb+heat+ as exact solution and as Dirichlet
boundary function with identifier $1$ to \verb+heatSpace+.
Thomas Witkowski's avatar
Thomas Witkowski committed
273
274
275
276
277
278
279
280
281

Now, we define the operators:

\begin{lstlisting}{}
  // ===== create operators =====
  double one = 1.0;
  double zero = 0.0;

  // create laplace
282
  Operator A(heatSpace.getFESpace());
Thomas Witkowski's avatar
Thomas Witkowski committed
283
284
285
286
287
288
  A.addSecondOrderTerm(new Laplace_SOT);
  A.setUhOld(heat.getOldSolution());
  if (*(heat.getThetaPtr()) != 0.0)
    heatSpace.addMatrixOperator(A, heat.getThetaPtr(), &one);
  if (*(heat.getTheta1Ptr()) != 0.0)
    heatSpace.addVectorOperator(A, heat.getTheta1Ptr(), &zero);
Thomas Witkowski's avatar
Thomas Witkowski committed
289
290
\end{lstlisting}

Thomas Witkowski's avatar
Thomas Witkowski committed
291
292
293
294
295
296
297
298
299
300
Operator \verb+A+ represents $-\Delta u$. It is used as matrix
operator on the left hand side with factor $\theta$ and as vector
operator on the right hand side with factor $-(1-\theta) = \theta -
1$. These assemble factors are the second arguments of
\verb+addMatrixOperator+ and \verb+addVectorOperator+. The third
argument is the factor used for estimation. In this example, the
estimator will consider the operator only on the left hand side with
factor $1$. On the right hand side the operator is applied to the
solution of the last timestep. So the old solution is handed to the
operator by \verb+setUhOld+.
Thomas Witkowski's avatar
Thomas Witkowski committed
301
302
303

\begin{lstlisting}{}
  // create zero order operator
304
  Operator C(heatSpace.getFESpace());
Thomas Witkowski's avatar
Thomas Witkowski committed
305
306
  C.addZeroOrderTerm(new Simple_ZOT);
  C.setUhOld(heat.getOldSolution());
307
308
  heatSpace.addMatrixOperator(C, heat.getInvTau(), heat.getInvTau());
  heatSpace.addVectorOperator(C, heat.getInvTau(), heat.getInvTau());
Thomas Witkowski's avatar
Thomas Witkowski committed
309
310
\end{lstlisting}

Thomas Witkowski's avatar
Thomas Witkowski committed
311
312
The \verb+Simple_ZOT+ of operator \verb+C+ represents the zero order
terms for the time discretization. On both sides of the equation $u$
313
314
315
316
are added with the factor $\frac{1}{\tau}$ for both, assembling and
error estimation. The inverse of the current timestep is returned by
the function \verb+getInvTau()+, which is a member of the class
\verb+ProblemInstat+.
Thomas Witkowski's avatar
Thomas Witkowski committed
317

Thomas Witkowski's avatar
Thomas Witkowski committed
318
319
Finally, the operator for the right hand side function $f$ is added
and the adaptation loop is started:
Thomas Witkowski's avatar
Thomas Witkowski committed
320
321
322

\begin{lstlisting}{}
  // create RHS operator
323
  Operator F(heatSpace.getFESpace());
Thomas Witkowski's avatar
Thomas Witkowski committed
324
325
  F.addZeroOrderTerm(new CoordsAtQP_ZOT(rhsFct));
  heatSpace.addVectorOperator(F);
Thomas Witkowski's avatar
Thomas Witkowski committed
326
327

  // ===== start adaption loop =====
Thomas Witkowski's avatar
Thomas Witkowski committed
328
  int errorCode = adaptInstat.adapt();
Thomas Witkowski's avatar
Thomas Witkowski committed
329
330
331
}
\end{lstlisting}

Thomas Witkowski's avatar
Thomas Witkowski committed
332
333
334
335
336
337
\verb+CoordsAtQP_ZOT+ is a zero order term that evaluates a given
function $fct$ at all needed quadrature points. At the left hand side,
it would represent the term $fct(x, t)\cdot u$, on the right hand
side, just $fct(x, t)$. Note that the old solution isn't given to the
operator here. Otherwise the term would represent $fct(x, t) \cdot
u^{old}$ on the right hand side.
Thomas Witkowski's avatar
Thomas Witkowski committed
338
339
340
341

\subsection{Parameter file}
\label{s:heat parameter}

Thomas Witkowski's avatar
Thomas Witkowski committed
342
343
In this section, we show only the relevant parts of the parameter file
\verb+heat.dat.2d+.
Thomas Witkowski's avatar
Thomas Witkowski committed
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368

First the parameter $\theta$ for the time discretization is defined:

\begin{lstlisting}{}
heat->theta:                         1.0
\end{lstlisting}

Then we define the initial timestep and the time interval:

\begin{lstlisting}{}
heat->adapt->timestep:               0.1
heat->adapt->start time:             0.0
heat->adapt->end time:               1.0
\end{lstlisting}

Now, tolerances are determined:

\begin{lstlisting}{}
heat->adapt->tolerance:              0.01
heat->adapt->rel space error:        0.5
heat->adapt->rel time error:         0.5
heat->adapt->time theta 1:           1.0
heat->adapt->time theta 2:           0.3
\end{lstlisting}

Thomas Witkowski's avatar
Thomas Witkowski committed
369
370
371
372
373
374
375
376
377
The total tolerance is divided in a space tolerance and a time
tolerance. The space tolerance is the maximal allowed space error,
given by the product of \verb+tolerance+ and \verb+rel space error+. It is reached by adaptive mesh refinements. The time tolerance
is the maximal allowed error, due to the timestep size. It is given by
the product of \verb+tolerance+ and \verb+rel time error+ and
\verb+time theta 1+. It is relevant, only if an implicit time strategy
with adaptive timestep size is used. The parameter \verb+time theta 2+
is used to enlarge the timestep, if the estimated time error falls
beneath a given threshold.
Thomas Witkowski's avatar
Thomas Witkowski committed
378
379
380
381
382
383
384

\begin{lstlisting}{}
heat->adapt->strategy:               1
heat->adapt->time delta 1:           0.7071
heat->adapt->time delta 2:           1.4142
\end{lstlisting}

Thomas Witkowski's avatar
Thomas Witkowski committed
385
386
387
388
389
390
391
392
393
394
395
If \verb+strategy+ is $0$, an explicit time strategy with fixed
timestep size is used. A value of $1$ stands for the implicit
strategy. The time tolerance is reached by successively multiplying
the timestep with \verb+time delta 1+. If the estimated timestep error
is smaller than the product of \verb+tolerance+ and
\verb+rel time error+ and \verb+time theta 2+ at the end of a
timestep, the timestep size is multiplied by \verb+time delta 2+.

The following lines determine, whether coarsening is allowed in
regions with sufficient small errors, and how many refinements or
coarsenings are performed for marked elements.
Thomas Witkowski's avatar
Thomas Witkowski committed
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
\begin{lstlisting}{}
heat->adapt->coarsen allowed:        1
heat->adapt->refine bisections:      2
heat->adapt->coarsen bisections:     2
\end{lstlisting}

Now, the output behavior is determined:

\begin{lstlisting}{}
heat->space->output->filename:       output/heat

heat->space->output->AMDiS format:   1
heat->space->output->AMDiS mesh ext: .mesh
heat->space->output->AMDiS data ext: .dat

heat->space->output->write every i-th timestep: 10

heat->space->output->append index:   1
heat->space->output->index length:   6
heat->space->output->index decimals: 3
\end{lstlisting}

Thomas Witkowski's avatar
Thomas Witkowski committed
418
419
420
421
422
423
424
In this example, all output filenames start with prefix
\verb+output/heat+ and end with the extensions \verb+.mesh+ and
\verb+.dat+. Output is written after every $10$th timestep. The time
of the single solution is added after the filename prefix with 6
letters, three of them are decimals. The solution for $t=0$ e.g. would
be written in the files \verb+output/heat00.000.mesh+ and
\verb+output/heat00.000.dat+.
Thomas Witkowski's avatar
Thomas Witkowski committed
425

Thomas Witkowski's avatar
Thomas Witkowski committed
426
427
428
Finally, we set parameter \verb+WAIT+ to $1$. So each call of the
macro \verb+WAIT+ in the application will lead to an interruption of
the program, until the \verb+return+ key is pressed.
Thomas Witkowski's avatar
Thomas Witkowski committed
429
430
431
432
433
434

\begin{lstlisting}{}
WAIT:                                1
\end{lstlisting}

\subsection{Macro file}
Thomas Witkowski's avatar
Thomas Witkowski committed
435
436
We again use the macro file \verb+macro/macro.stand.2d+, which was
described in Section \ref{s:ellipt macro}.
Thomas Witkowski's avatar
Thomas Witkowski committed
437
438

\subsection{Output}
Thomas Witkowski's avatar
Thomas Witkowski committed
439
440
441
442
443
444
As mentioned above, the output files look like
\verb+output/heat00.000.mesh+ and
\verb+output/heat00.000.dat+. Depending on the corresponding value in
the parameter file only the solution after every $i$-th timestep is
written. In Figure \ref{f:heat}, the solution at three timesteps is
visualized.
Thomas Witkowski's avatar
Thomas Witkowski committed
445
446
447
448
449
450
451
452
453
454
455
456

\begin{figure}
\center
\begin{tabular}{ccc}
\includegraphics[width=0.3\textwidth]{fig/heat1.jpg}&
\includegraphics[width=0.3\textwidth]{fig/heat2.jpg}&
\includegraphics[width=0.3\textwidth]{fig/heat3.jpg}\\
$t \approx 0.2 $ & $t \approx 0.47 $ & $t \approx 0.81$
\end{tabular}
\caption{The solution at three different timesteps.}
\label{f:heat}
\end{figure}