Skip to content
GitLab
Menu
Projects
Groups
Snippets
/
Help
Help
Support
Community forum
Keyboard shortcuts
?
Submit feedback
Contribute to GitLab
Sign in
Toggle navigation
Menu
Open sidebar
iwr
amdis
Commits
83b6de2b
Commit
83b6de2b
authored
May 25, 2009
by
Thomas Witkowski
Browse files
Adapted demos to new AMDiS (DELETE/delete and dirichlet boundary conditions).
parent
b8488200
Changes
26
Hide whitespace changes
Inline
Side-by-side
demo/src/ball.cc
View file @
83b6de2b
...
...
@@ -7,36 +7,27 @@ using namespace AMDiS;
// ===== function definitions ================================================
// ===========================================================================
/** \brief
* Dirichlet boundary function
*/
/// Dirichlet boundary function
class
G
:
public
AbstractFunction
<
double
,
WorldVector
<
double
>
>
{
public:
MEMORY_MANAGED
(
G
);
/** \brief
* Implementation of AbstractFunction::operator().
*/
double
operator
()(
const
WorldVector
<
double
>&
x
)
const
{
/// Implementation of AbstractFunction::operator().
double
operator
()(
const
WorldVector
<
double
>&
x
)
const
{
return
exp
(
-
10.0
*
(
x
*
x
));
}
};
/** \brief
* RHS function
*/
/// RHS function
class
F
:
public
AbstractFunction
<
double
,
WorldVector
<
double
>
>
{
public:
MEMORY_MANAGED
(
F
);
F
(
int
degree
)
:
AbstractFunction
<
double
,
WorldVector
<
double
>
>
(
degree
)
{};
F
(
int
degree
)
:
AbstractFunction
<
double
,
WorldVector
<
double
>
>
(
degree
)
{}
/** \brief
* Implementation of AbstractFunction::operator().
*/
double
operator
()(
const
WorldVector
<
double
>&
x
)
const
{
/// Implementation of AbstractFunction::operator().
double
operator
()(
const
WorldVector
<
double
>&
x
)
const
{
int
dim
=
x
.
getSize
();
double
r2
=
x
*
x
;
double
ux
=
exp
(
-
10.0
*
r2
);
...
...
@@ -61,7 +52,7 @@ int main(int argc, char* argv[])
// ===== create projection =====
WorldVector
<
double
>
ballCenter
;
ballCenter
.
set
(
0.0
);
NEW
BallProject
(
1
,
new
BallProject
(
1
,
BOUNDARY_PROJECTION
,
ballCenter
,
1.0
);
...
...
@@ -71,25 +62,25 @@ int main(int argc, char* argv[])
ball
.
initialize
(
INIT_ALL
);
// === create adapt info ===
AdaptInfo
*
adaptInfo
=
NEW
AdaptInfo
(
"ball->adapt"
);
AdaptInfo
*
adaptInfo
=
new
AdaptInfo
(
"ball->adapt"
);
// === create adapt ===
AdaptStationary
*
adapt
=
NEW
AdaptStationary
(
"ball->adapt"
,
AdaptStationary
*
adapt
=
new
AdaptStationary
(
"ball->adapt"
,
&
ball
,
adaptInfo
);
// ===== add boundary conditions =====
ball
.
addDirichletBC
(
1
,
NEW
G
);
ball
.
addDirichletBC
(
1
,
new
G
);
// ===== create matrix operator =====
Operator
matrixOperator
(
Operator
::
MATRIX_OPERATOR
,
ball
.
getFESpace
());
matrixOperator
.
addSecondOrderTerm
(
NEW
Laplace_SOT
);
matrixOperator
.
addSecondOrderTerm
(
new
Laplace_SOT
);
ball
.
addMatrixOperator
(
&
matrixOperator
);
// ===== create rhs operator =====
int
degree
=
ball
.
getFESpace
()
->
getBasisFcts
()
->
getDegree
();
Operator
rhsOperator
(
Operator
::
VECTOR_OPERATOR
,
ball
.
getFESpace
());
rhsOperator
.
addZeroOrderTerm
(
NEW
CoordsAtQP_ZOT
(
NEW
F
(
degree
)));
rhsOperator
.
addZeroOrderTerm
(
new
CoordsAtQP_ZOT
(
new
F
(
degree
)));
ball
.
addVectorOperator
(
&
rhsOperator
);
// ===== start adaption loop =====
...
...
demo/src/bunny.cc
View file @
83b6de2b
...
...
@@ -7,36 +7,27 @@ using namespace AMDiS;
// ===== function definitions ================================================
// ===========================================================================
/** \brief
* RHS function
*/
/// RHS function
class
F
:
public
AbstractFunction
<
double
,
WorldVector
<
double
>
>
{
public:
MEMORY_MANAGED
(
F
);
F
(
int
degree
)
:
AbstractFunction
<
double
,
WorldVector
<
double
>
>
(
degree
)
{}
F
(
int
degree
)
:
AbstractFunction
<
double
,
WorldVector
<
double
>
>
(
degree
)
{};
/** \brief
* Implementation of AbstractFunction::operator().
*/
double
operator
()(
const
WorldVector
<
double
>&
x
)
const
{
/// Implementation of AbstractFunction::operator().
double
operator
()(
const
WorldVector
<
double
>&
x
)
const
{
return
-
2
*
x
[
0
];
}
};
/** \brief
* boundary
*/
/// boundary
class
G
:
public
AbstractFunction
<
double
,
WorldVector
<
double
>
>
{
public:
MEMORY_MANAGED
(
G
);
/** \brief
* Implementation of AbstractFunction::operator().
*/
double
operator
()(
const
WorldVector
<
double
>&
x
)
const
{
/// Implementation of AbstractFunction::operator().
double
operator
()(
const
WorldVector
<
double
>&
x
)
const
{
return
10000.0
;
}
};
...
...
@@ -58,7 +49,7 @@ int main(int argc, char* argv[])
// ===== create projection =====
WorldVector
<
double
>
ballCenter
;
ballCenter
.
set
(
0.0
);
NEW
BallProject
(
1
,
new
BallProject
(
1
,
VOLUME_PROJECTION
,
ballCenter
,
1.0
);
...
...
@@ -68,20 +59,20 @@ int main(int argc, char* argv[])
bunny
.
initialize
(
INIT_ALL
);
// === create adapt info ===
AdaptInfo
*
adaptInfo
=
NEW
AdaptInfo
(
"bunny->adapt"
);
AdaptInfo
*
adaptInfo
=
new
AdaptInfo
(
"bunny->adapt"
);
// === create adapt ===
AdaptStationary
*
adapt
=
NEW
AdaptStationary
(
"bunny->adapt"
,
AdaptStationary
*
adapt
=
new
AdaptStationary
(
"bunny->adapt"
,
&
bunny
,
adaptInfo
);
// ===== add boundary conditions =====
//bunny.addDirichletBC(1111,
NEW
G);
//bunny.addDirichletBC(1111,
new
G);
// ===== create matrix operator =====
Operator
matrixOperator
(
Operator
::
MATRIX_OPERATOR
,
bunny
.
getFESpace
());
matrixOperator
.
addSecondOrderTerm
(
NEW
Laplace_SOT
);
matrixOperator
.
addSecondOrderTerm
(
new
Laplace_SOT
);
bunny
.
addMatrixOperator
(
&
matrixOperator
);
// ===== create rhs operator =====
...
...
@@ -90,7 +81,7 @@ int main(int argc, char* argv[])
int
degree
=
bunny
.
getFESpace
()
->
getBasisFcts
()
->
getDegree
();
rhsOperator
.
addZeroOrderTerm
(
NEW
CoordsAtQP_ZOT
(
NEW
F
(
degree
)));
rhsOperator
.
addZeroOrderTerm
(
new
CoordsAtQP_ZOT
(
new
F
(
degree
)));
bunny
.
addVectorOperator
(
&
rhsOperator
);
// ===== start adaption loop =====
...
...
demo/src/couple.cc
View file @
83b6de2b
...
...
@@ -6,9 +6,8 @@ using namespace AMDiS;
class
G
:
public
AbstractFunction
<
double
,
WorldVector
<
double
>
>
{
public:
MEMORY_MANAGED
(
G
);
double
operator
()(
const
WorldVector
<
double
>&
x
)
const
{
double
operator
()(
const
WorldVector
<
double
>&
x
)
const
{
return
exp
(
-
10.0
*
(
x
*
x
));
}
};
...
...
@@ -16,11 +15,10 @@ public:
class
F
:
public
AbstractFunction
<
double
,
WorldVector
<
double
>
>
{
public:
MEMORY_MANAGED
(
F
);
F
(
int
degree
)
:
AbstractFunction
<
double
,
WorldVector
<
double
>
>
(
degree
)
{};
F
(
int
degree
)
:
AbstractFunction
<
double
,
WorldVector
<
double
>
>
(
degree
)
{}
double
operator
()(
const
WorldVector
<
double
>&
x
)
const
{
double
operator
()(
const
WorldVector
<
double
>&
x
)
const
{
int
dow
=
x
.
getSize
();
double
r2
=
(
x
*
x
);
double
ux
=
exp
(
-
10.0
*
r2
);
...
...
@@ -111,9 +109,10 @@ class Identity : public AbstractFunction<double, double>
public:
MEMORY_MANAGED
(
Identity
);
Identity
(
int
degree
)
:
AbstractFunction
<
double
,
double
>
(
degree
)
{}
;
Identity
(
int
degree
)
:
AbstractFunction
<
double
,
double
>
(
degree
)
{}
double
operator
()(
const
double
&
x
)
const
{
double
operator
()(
const
double
&
x
)
const
{
return
x
;
}
};
...
...
@@ -145,35 +144,35 @@ int main(int argc, char* argv[])
adoptFlag
);
// ===== add boundary conditions for problem1 =====
problem1
.
addDirichletBC
(
1
,
NEW
G
);
problem1
.
addDirichletBC
(
1
,
new
G
);
// ===== add boundary conditions for problem1 =====
//problem2.addDirichletBC(1,
NEW
G);
//problem2.addDirichletBC(1,
new
G);
// ===== create operators for problem1 =====
Operator
matrixOperator1
(
Operator
::
MATRIX_OPERATOR
,
problem1
.
getFESpace
());
matrixOperator1
.
addSecondOrderTerm
(
NEW
Laplace_SOT
);
matrixOperator1
.
addSecondOrderTerm
(
new
Laplace_SOT
);
problem1
.
addMatrixOperator
(
&
matrixOperator1
);
int
degree
=
problem1
.
getFESpace
()
->
getBasisFcts
()
->
getDegree
();
Operator
rhsOperator1
(
Operator
::
VECTOR_OPERATOR
,
problem1
.
getFESpace
());
rhsOperator1
.
addZeroOrderTerm
(
NEW
CoordsAtQP_ZOT
(
NEW
F
(
degree
)));
rhsOperator1
.
addZeroOrderTerm
(
new
CoordsAtQP_ZOT
(
new
F
(
degree
)));
problem1
.
addVectorOperator
(
&
rhsOperator1
);
// ===== create operators for problem2 =====
Operator
matrixOperator2
(
Operator
::
MATRIX_OPERATOR
,
problem2
.
getFESpace
());
matrixOperator2
.
addZeroOrderTerm
(
NEW
Simple_ZOT
);
matrixOperator2
.
addZeroOrderTerm
(
new
Simple_ZOT
);
problem2
.
addMatrixOperator
(
&
matrixOperator2
);
Operator
rhsOperator2
(
Operator
::
VECTOR_OPERATOR
,
problem2
.
getFESpace
());
rhsOperator2
.
addZeroOrderTerm
(
NEW
VecAtQP_ZOT
(
problem1
.
getSolution
(),
NEW
Identity
(
degree
)));
rhsOperator2
.
addZeroOrderTerm
(
new
VecAtQP_ZOT
(
problem1
.
getSolution
(),
new
Identity
(
degree
)));
problem2
.
addVectorOperator
(
&
rhsOperator2
);
// ===== create adaptation loop and iteration interface =====
AdaptInfo
*
adaptInfo
=
NEW
AdaptInfo
(
"couple->adapt"
,
1
);
AdaptInfo
*
adaptInfo
=
new
AdaptInfo
(
"couple->adapt"
,
1
);
MyCoupledIteration
coupledIteration
(
&
problem1
,
&
problem2
);
AdaptStationary
*
adapt
=
NEW
AdaptStationary
(
"couple->adapt"
,
AdaptStationary
*
adapt
=
new
AdaptStationary
(
"couple->adapt"
,
&
coupledIteration
,
adaptInfo
);
...
...
demo/src/ellipt.cc
View file @
83b6de2b
...
...
@@ -3,35 +3,30 @@
using
namespace
AMDiS
;
using
namespace
std
;
// ===== function definitions //
/** \brief
* Dirichlet boundary function
*/
// ===========================================================================
// ===== function definitions ================================================
// ===========================================================================
/// Dirichlet boundary function
class
G
:
public
AbstractFunction
<
double
,
WorldVector
<
double
>
>
{
public:
MEMORY_MANAGED
(
G
);
/** \brief
* Implementation of AbstractFunction::operator().
*/
double
operator
()(
const
WorldVector
<
double
>&
x
)
const
{
/// Implementation of AbstractFunction::operator().
double
operator
()(
const
WorldVector
<
double
>&
x
)
const
{
return
exp
(
-
10.0
*
(
x
*
x
));
}
};
/** \brief
* RHS function
*/
/// RHS function
class
F
:
public
AbstractFunction
<
double
,
WorldVector
<
double
>
>
{
public:
MEMORY_MANAGED
(
F
);
/** \brief
* Implementation of AbstractFunction::operator().
*/
double
operator
()(
const
WorldVector
<
double
>&
x
)
const
{
/// Implementation of AbstractFunction::operator().
double
operator
()(
const
WorldVector
<
double
>&
x
)
const
{
int
dow
=
x
.
getSize
();
double
r2
=
(
x
*
x
);
double
ux
=
exp
(
-
10.0
*
r2
);
...
...
@@ -39,7 +34,10 @@ public:
}
};
// // ===== main program //
// ===========================================================================
// ===== main program ========================================================
// ===========================================================================
int
main
(
int
argc
,
char
*
argv
[])
{
FUNCNAME
(
"main"
);
...
...
@@ -55,24 +53,24 @@ int main(int argc, char* argv[])
ellipt
.
initialize
(
INIT_ALL
);
// === create adapt info ===
AdaptInfo
*
adaptInfo
=
NEW
AdaptInfo
(
"ellipt->adapt"
);
AdaptInfo
*
adaptInfo
=
new
AdaptInfo
(
"ellipt->adapt"
);
// === create adapt ===
AdaptStationary
*
adapt
=
NEW
AdaptStationary
(
"ellipt->adapt"
,
AdaptStationary
*
adapt
=
new
AdaptStationary
(
"ellipt->adapt"
,
&
ellipt
,
adaptInfo
);
// ===== add boundary conditions =====
ellipt
.
addDirichletBC
(
1
,
NEW
G
);
ellipt
.
addDirichletBC
(
1
,
new
G
);
// ===== create matrix operator =====
Operator
matrixOperator
(
Operator
::
MATRIX_OPERATOR
,
ellipt
.
getFESpace
());
matrixOperator
.
addSecondOrderTerm
(
NEW
Laplace_SOT
);
matrixOperator
.
addSecondOrderTerm
(
new
Laplace_SOT
);
ellipt
.
addMatrixOperator
(
&
matrixOperator
);
// ===== create rhs operator =====
Operator
rhsOperator
(
Operator
::
VECTOR_OPERATOR
,
ellipt
.
getFESpace
());
rhsOperator
.
addZeroOrderTerm
(
NEW
CoordsAtQP_ZOT
(
NEW
F
()));
rhsOperator
.
addZeroOrderTerm
(
new
CoordsAtQP_ZOT
(
new
F
()));
ellipt
.
addVectorOperator
(
&
rhsOperator
);
// ===== start adaption loop =====
...
...
demo/src/heat.cc
View file @
83b6de2b
...
...
@@ -7,38 +7,29 @@ using namespace AMDiS;
// ===== function definitions ================================================
// ===========================================================================
/** \brief
* Dirichlet boundary function
*/
/// Dirichlet boundary function
class
G
:
public
AbstractFunction
<
double
,
WorldVector
<
double
>
>
,
public
TimedObject
{
public:
MEMORY_MANAGED
(
G
);
/** \brief
* Implementation of AbstractFunction::operator().
*/
double
operator
()(
const
WorldVector
<
double
>&
x
)
const
{
/// Implementation of AbstractFunction::operator().
double
operator
()(
const
WorldVector
<
double
>&
x
)
const
{
return
sin
(
M_PI
*
(
*
timePtr
))
*
exp
(
-
10.0
*
(
x
*
x
));
}
};
/** \brief
* RHS function
*/
/// RHS function
class
F
:
public
AbstractFunction
<
double
,
WorldVector
<
double
>
>
,
public
TimedObject
{
public:
MEMORY_MANAGED
(
F
);
F
(
int
degree
)
:
AbstractFunction
<
double
,
WorldVector
<
double
>
>
(
degree
)
{}
F
(
int
degree
)
:
AbstractFunction
<
double
,
WorldVector
<
double
>
>
(
degree
)
{};
/** \brief
* Implementation of AbstractFunction::operator().
*/
double
operator
()(
const
WorldVector
<
double
>&
x
)
const
{
/// Implementation of AbstractFunction::operator().
double
operator
()(
const
WorldVector
<
double
>&
x
)
const
{
int
dim
=
x
.
getSize
();
double
r2
=
x
*
x
;
double
ux
=
sin
(
M_PI
*
(
*
timePtr
))
*
exp
(
-
10.0
*
r2
);
...
...
@@ -51,17 +42,10 @@ public:
// ===== instationary problem ================================================
// ===========================================================================
/** \brief
* Instationary problem
*/
/// Instationary problem
class
Heat
:
public
ProblemInstatScal
{
public:
MEMORY_MANAGED
(
Heat
);
/** \brief
* Constructor
*/
Heat
(
ProblemScal
*
heatSpace
)
:
ProblemInstatScal
(
"heat"
,
heatSpace
)
{
...
...
@@ -77,121 +61,98 @@ public:
theta1
=
theta
-
1
;
}
// ===== ProblemInstatBase methods ===================================
/** \brief
* set the time in all needed functions!
*/
void
setTime
(
AdaptInfo
*
adaptInfo
)
{
/// set the time in all needed functions!
void
setTime
(
AdaptInfo
*
adaptInfo
)
{
rhsTime
=
adaptInfo
->
getTime
()
-
(
1
-
theta
)
*
adaptInfo
->
getTimestep
();
boundaryTime
=
adaptInfo
->
getTime
();
tau1
=
1.0
/
adaptInfo
->
getTimestep
();
}
void
closeTimestep
(
AdaptInfo
*
adaptInfo
)
{
void
closeTimestep
(
AdaptInfo
*
adaptInfo
)
{
ProblemInstatScal
::
closeTimestep
(
adaptInfo
);
WAIT
;
}
// ===== initial problem methods =====================================
/** \brief
* Used by \ref problemInitial to solve the system of the initial problem
*/
void
solve
(
AdaptInfo
*
adaptInfo
)
{
/// Used by \ref problemInitial to solve the system of the initial problem
void
solve
(
AdaptInfo
*
adaptInfo
)
{
problemStat
->
getSolution
()
->
interpol
(
exactSolution
);
}
/** \brief
* Used by \ref problemInitial to do error estimation for the initial
* problem.
*/
void
estimate
(
AdaptInfo
*
adaptInfo
)
{
/// Used by \ref problemInitial to do error estimation for the initial problem.
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
);
}
// ===== setting methods ===============================================
/** \brief
* Sets \ref exactSolution;
*/
void
setExactSolution
(
AbstractFunction
<
double
,
WorldVector
<
double
>
>
*
fct
)
{
/// Sets \ref exactSolution;
void
setExactSolution
(
AbstractFunction
<
double
,
WorldVector
<
double
>
>
*
fct
)
{
exactSolution
=
fct
;
}
// ===== getting methods ===============================================
/** \brief
* Returns pointer to \ref theta.
*/
double
*
getThetaPtr
()
{
/// Returns pointer to \ref theta.
double
*
getThetaPtr
()
{
return
&
theta
;
}
/** \brief
* Returns pointer to \ref theta1.
*/
double
*
getTheta1Ptr
()
{
/// Returns pointer to \ref theta1.
double
*
getTheta1Ptr
()
{
return
&
theta1
;
}
/** \brief
* Returns pointer to \ref tau1
*/
double
*
getTau1Ptr
()
{
/// Returns pointer to \ref tau1
double
*
getTau1Ptr
()
{
return
&
tau1
;
}
/** \brief
* Returns pointer to \ref rhsTime.
*/
double
*
getRHSTimePtr
()
{
/// Returns pointer to \ref rhsTime.
double
*
getRHSTimePtr
()
{
return
&
rhsTime
;
}
/** \brief
* Returns pointer to \ref theta1.
*/
double
*
getBoundaryTimePtr
()
{
/// Returns pointer to \ref theta1.
double
*
getBoundaryTimePtr
()
{
return
&
boundaryTime
;
}
private:
/** \brief
* Used for theta scheme.
*/
/// Used for theta scheme.
double
theta
;
/** \brief
* theta - 1
*/
/// theta - 1
double
theta1
;
/** \brief
* 1.0 / timestep
*/
/// 1.0 / timestep
double
tau1
;
/** \brief
* time for right hand side functions.
*/
/// time for right hand side functions.
double
rhsTime
;
/** \brief
* time for boundary functions.
*/
/// time for boundary functions.
double
boundaryTime
;
/** \brief
* Pointer to boundary function. Needed for initial problem.
*/
/// Pointer to boundary function. Needed for initial problem.
AbstractFunction
<
double
,
WorldVector
<
double
>
>
*
exactSolution
;
};
...
...
@@ -209,7 +170,7 @@ int main(int argc, char** argv)
Parameters
::
readArgv
(
argc
,
argv
);
// ===== create and init stationary problem =====
ProblemScal
*
heatSpace
=
NEW
ProblemScal
(
"heat->space"
);
ProblemScal
*
heatSpace
=
new
ProblemScal
(
"heat->space"
);
heatSpace
->
initialize
(
INIT_ALL
);
// ===== create instationary problem =====
...
...
@@ -217,21 +178,21 @@ int main(int argc, char** argv)
heat
->
initialize
(
INIT_ALL
);
// create adapt info
AdaptInfo
*
adaptInfo
=
NEW
AdaptInfo
(
"heat->adapt"
);
AdaptInfo
*
adaptInfo
=
new
AdaptInfo
(
"heat->adapt"
);
// create initial adapt info
AdaptInfo
*
adaptInfoInitial
=
NEW
AdaptInfo
(
"heat->initial->adapt"
);
AdaptInfo
*
adaptInfoInitial
=
new
AdaptInfo
(
"heat->initial->adapt"
);
// create instationary adapt
AdaptInstationary
*
adaptInstat
=
NEW
AdaptInstationary
(
"heat->adapt"
,
new
AdaptInstationary
(
"heat->adapt"
,
heatSpace
,
adaptInfo
,
heat
,
adaptInfoInitial
);
// ===== create boundary functions =====
G
*
boundaryFct
=
NEW
G
;
G
*
boundaryFct
=
new
G
;
boundaryFct
->
setTimePtr
(
heat
->
getBoundaryTimePtr
());
heat
->
setExactSolution
(
boundaryFct
);
...
...
@@ -239,7 +200,7 @@ int main(int argc, char** argv)
// ===== create rhs functions =====
int
degree
=
heatSpace
->
getFESpace
()
->
getBasisFcts
()
->
getDegree
();
F
*
rhsFct
=
NEW
F
(
degree
);
F
*
rhsFct
=
new
F
(
degree
);
rhsFct
->
setTimePtr
(
heat
->
getRHSTimePtr
());
// ===== create operators =====
...
...
@@ -247,7 +208,7 @@ int main(int argc, char** argv)
double
zero
=
0.0
;
// create laplace
Operator
*
A
=
NEW
Operator
(
Operator
::
MATRIX_OPERATOR
|
Operator
*
A
=
new
Operator
(
Operator
::
MATRIX_OPERATOR
|
Operator
::
VECTOR_OPERATOR
,
heatSpace
->
getFESpace
());
...
...
@@ -255,18 +216,18 @@ int main(int argc, char** argv)
A
->
setUhOld
(
heat
->
getOldSolution
());
if
((
*
(
heat
->
getThetaPtr
()))
!=
0.0
)
if
((
*
(
heat
->
getThetaPtr
()))
!=
0.0
)
heatSpace
->
addMatrixOperator
(
A
,
heat
->
getThetaPtr
(),
&
one
);