Page Menu
Home
c4science
Search
Configure Global Search
Log In
Files
F65813676
influence.hh
No One
Temporary
Actions
Download File
Edit File
Delete File
View Transforms
Subscribe
Mute Notifications
Award Token
Subscribers
None
File Metadata
Details
File Info
Storage
Attached
Created
Thu, Jun 6, 09:59
Size
16 KB
Mime Type
text/x-c++
Expires
Sat, Jun 8, 09:59 (2 d)
Engine
blob
Format
Raw Data
Handle
18131296
Attached To
rTAMAAS tamaas
influence.hh
View Options
/**
* @file
*
* @author Lucas Frérot <lucas.frerot@epfl.ch>
*
* @section LICENSE
*
* Copyright (©) 2017 EPFL (Ecole Polytechnique Fédérale de
* Lausanne) Laboratory (LSMS - Laboratoire de Simulation en Mécanique des
* Solides)
*
* Tamaas is free software: you can redistribute it and/or modify it under the
* terms of the GNU Lesser General Public License as published by the Free
* Software Foundation, either version 3 of the License, or (at your option) any
* later version.
*
* Tamaas is distributed in the hope that it will be useful, but WITHOUT ANY
* WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
* A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more
* details.
*
* You should have received a copy of the GNU Lesser General Public License
* along with Tamaas. If not, see <http://www.gnu.org/licenses/>.
*
*/
/* -------------------------------------------------------------------------- */
#ifndef INFLUENCE_HH
#define INFLUENCE_HH
/* -------------------------------------------------------------------------- */
#include "loop.hh"
#include "static_types.hh"
#include <type_traits>
/* -------------------------------------------------------------------------- */
__BEGIN_TAMAAS__
/* -------------------------------------------------------------------------- */
namespace
influence
{
namespace
{
const
Complex
I
(
0
,
1
);
/// < imaginary unit
constexpr
inline
Real
sign
(
bool
upper
)
{
return
(
upper
)
?
1.
:
-
1.
;
}
}
template
<
bool
conjugate
,
UInt
dim_q
>
// dim + 1 to help dim template deduction
inline
Vector
<
Complex
,
dim_q
+
1
>
computeD
(
VectorProxy
<
const
Real
,
dim_q
>&
q
)
{
constexpr
Real
s
=
-
sign
(
conjugate
);
const
auto
q_norm
=
q
.
l2norm
();
return
{{{
s
*
I
*
q
(
0
)
/
q_norm
,
s
*
I
*
q
(
1
)
/
q_norm
,
1
}}};
}
template
<
UInt
dim_q
>
// dim + 1 to help dim template deduction
inline
Vector
<
Complex
,
dim_q
+
1
>
computeD2
(
VectorProxy
<
const
Real
,
dim_q
>&
q
)
{
const
auto
q_norm
=
q
.
l2norm
();
return
{{{
I
*
q
(
0
)
/
q_norm
,
I
*
q
(
1
)
/
q_norm
,
0
}}};
}
/// Class for the Kelvin tensor
template
<
UInt
dim
,
UInt
derivative_order
>
class
Kelvin
;
/// See details in notebook and Marc Bonnet's notes
template
<>
class
Kelvin
<
3
,
0
>
{
protected
:
constexpr
static
UInt
dim
=
3
;
public
:
Kelvin
(
Real
mu
,
Real
nu
)
:
mu
(
mu
),
b
(
4
*
(
1
-
nu
))
{}
template
<
bool
upper
,
typename
ST
>
inline
Vector
<
Complex
,
dim
>
applyU0
(
VectorProxy
<
const
Real
,
dim
-
1
>&
q
,
StaticVector
<
Complex
,
ST
,
dim
>&
f
)
const
;
template
<
bool
upper
,
typename
ST
>
inline
Vector
<
Complex
,
dim
>
applyU1
(
VectorProxy
<
const
Real
,
dim
-
1
>&
q
,
StaticVector
<
Complex
,
ST
,
dim
>&
f
)
const
;
protected
:
const
Real
mu
,
b
;
};
/* -------------------------------------------------------------------------- */
/* Base Kelvin tensor */
/* -------------------------------------------------------------------------- */
template
<
bool
upper
,
typename
ST
>
//< this is ignored
inline
Vector
<
Complex
,
Kelvin
<
3
,
0
>::
dim
>
Kelvin
<
3
,
0
>::
applyU0
(
VectorProxy
<
const
Real
,
dim
-
1
>&
q
,
StaticVector
<
Complex
,
ST
,
dim
>&
f
)
const
{
auto
d2
=
computeD2
(
q
);
auto
tmp
=
d2
.
dot
(
f
)
*
d2
;
tmp
+=
b
*
f
;
tmp
(
2
)
-=
f
(
2
);
tmp
*=
1.
/
(
2.
*
mu
*
q
.
l2norm
()
*
b
);
return
tmp
;
}
template
<
bool
upper
,
typename
ST
>
inline
Vector
<
Complex
,
Kelvin
<
3
,
0
>::
dim
>
Kelvin
<
3
,
0
>::
applyU1
(
VectorProxy
<
const
Real
,
dim
-
1
>&
q
,
StaticVector
<
Complex
,
ST
,
dim
>&
f
)
const
{
auto
d
=
computeD
<
upper
>
(
q
);
auto
tmp
=
d
.
dot
(
f
)
*
d
;
tmp
*=
sign
(
upper
)
/
(
2.
*
mu
*
q
.
l2norm
()
*
b
);
return
tmp
;
}
/* -------------------------------------------------------------------------- */
/* Kelvin first derivative */
/* -------------------------------------------------------------------------- */
template
<>
class
Kelvin
<
3
,
1
>
:
protected
Kelvin
<
3
,
0
>
{
using
parent
=
Kelvin
<
3
,
0
>
;
protected
:
static
constexpr
UInt
dim
=
parent
::
dim
;
public
:
using
parent
::
parent
;
template
<
bool
upper
,
typename
ST
>
inline
Vector
<
Complex
,
dim
>
applyU0
(
VectorProxy
<
const
Real
,
dim
-
1
>&
q
,
StaticMatrix
<
Complex
,
ST
,
dim
,
dim
>&
f
)
const
;
template
<
bool
upper
,
typename
ST
>
inline
Vector
<
Complex
,
dim
>
applyU1
(
VectorProxy
<
const
Real
,
dim
-
1
>&
q
,
StaticMatrix
<
Complex
,
ST
,
dim
,
dim
>&
f
)
const
;
};
template
<
bool
upper
,
typename
ST
>
inline
Vector
<
Complex
,
Kelvin
<
3
,
1
>::
dim
>
Kelvin
<
3
,
1
>::
applyU0
(
VectorProxy
<
const
Real
,
dim
-
1
>&
q
,
StaticMatrix
<
Complex
,
ST
,
dim
,
dim
>&
f
)
const
{
auto
tmp
=
f
*
computeD
<
upper
>
(
q
);
auto
res
=
Kelvin
<
3
,
0
>::
applyU0
<
false
>
(
q
,
tmp
);
res
*=
-
sign
(
upper
);
Vector
<
Real
,
dim
>
e3
{{{
0
,
0
,
1
}}};
tmp
.
template
mul
<
false
>
(
f
,
e3
);
res
+=
Kelvin
<
3
,
0
>::
applyU1
<
upper
>
(
q
,
tmp
);
res
*=
q
.
l2norm
();
return
res
;
}
template
<
bool
upper
,
typename
ST
>
inline
Vector
<
Complex
,
Kelvin
<
3
,
1
>::
dim
>
Kelvin
<
3
,
1
>::
applyU1
(
VectorProxy
<
const
Real
,
dim
-
1
>&
q
,
StaticMatrix
<
Complex
,
ST
,
dim
,
dim
>&
f
)
const
{
auto
tmp
=
f
*
computeD
<
upper
>
(
q
);
auto
res
=
Kelvin
<
3
,
0
>::
applyU1
<
upper
>
(
q
,
tmp
);
res
*=
-
sign
(
upper
)
*
q
.
l2norm
();
return
res
;
}
/* -------------------------------------------------------------------------- */
/* Kelvin second derivative */
/* -------------------------------------------------------------------------- */
template
<>
class
Kelvin
<
3
,
2
>
:
protected
Kelvin
<
3
,
1
>
{
using
parent
=
Kelvin
<
3
,
1
>
;
static
constexpr
UInt
dim
=
parent
::
dim
;
public
:
using
parent
::
parent
;
template
<
bool
upper
,
typename
ST
>
inline
Matrix
<
Complex
,
dim
,
dim
>
applyU0
(
VectorProxy
<
const
Real
,
dim
-
1
>&
q
,
StaticMatrix
<
Complex
,
ST
,
dim
,
dim
>&
f
)
const
;
template
<
bool
upper
,
typename
ST
>
inline
Matrix
<
Complex
,
dim
,
dim
>
applyU1
(
VectorProxy
<
const
Real
,
dim
-
1
>&
q
,
StaticMatrix
<
Complex
,
ST
,
dim
,
dim
>&
f
)
const
;
template
<
typename
ST
>
inline
Matrix
<
Complex
,
dim
,
dim
>
applyDiscontinuityTerm
(
VectorProxy
<
const
Real
,
dim
-
1
>&
q
,
StaticMatrix
<
Complex
,
ST
,
dim
,
dim
>&
f
)
const
;
};
template
<
bool
upper
,
typename
ST
>
inline
Matrix
<
Complex
,
Kelvin
<
3
,
2
>::
dim
,
Kelvin
<
3
,
2
>::
dim
>
Kelvin
<
3
,
2
>::
applyU0
(
VectorProxy
<
const
Real
,
dim
-
1
>&
q
,
StaticMatrix
<
Complex
,
ST
,
dim
,
dim
>&
f
)
const
{
auto
tmp
=
Kelvin
<
3
,
1
>::
applyU0
<
upper
>
(
q
,
f
);
tmp
*=
-
sign
(
upper
);
auto
res
=
outer
(
computeD
<
upper
>
(
q
),
tmp
);
Vector
<
Real
,
dim
>
e3
{{{
0
,
0
,
1
}}};
res
+=
outer
(
e3
,
Kelvin
<
3
,
1
>::
applyU1
<
upper
>
(
q
,
f
));
res
*=
-
1.
*
q
.
l2norm
();
// << -1 for grad_x
return
res
;
}
template
<
bool
upper
,
typename
ST
>
inline
Matrix
<
Complex
,
Kelvin
<
3
,
2
>::
dim
,
Kelvin
<
3
,
2
>::
dim
>
Kelvin
<
3
,
2
>::
applyU1
(
VectorProxy
<
const
Real
,
dim
-
1
>&
q
,
StaticMatrix
<
Complex
,
ST
,
dim
,
dim
>&
f
)
const
{
auto
tmp
=
Kelvin
<
3
,
1
>::
applyU1
<
upper
>
(
q
,
f
);
auto
res
=
outer
(
computeD
<
upper
>
(
q
),
tmp
);
res
*=
-
1.
*
-
sign
(
upper
)
*
q
.
l2norm
();
// << -1 for grad_x
return
res
;
}
template
<
typename
ST
>
inline
Matrix
<
Complex
,
Kelvin
<
3
,
2
>::
dim
,
Kelvin
<
3
,
2
>::
dim
>
Kelvin
<
3
,
2
>::
applyDiscontinuityTerm
(
VectorProxy
<
const
Real
,
dim
-
1
>&
q
,
StaticMatrix
<
Complex
,
ST
,
dim
,
dim
>&
f
)
const
{
Vector
<
Real
,
dim
>
e3
{{{
0
,
0
,
1
}}};
auto
res
=
(
f
*
e3
);
res
*=
-
this
->
b
;
res
(
2
)
+=
2.
*
f
(
2
,
2
);
e3
(
2
)
*=
-
1
/
(
this
->
mu
*
this
->
b
);
// << -1 for grad_x
return
outer
(
e3
,
res
);
}
/* -------------------------------------------------------------------------- */
/* Implementation */
/* -------------------------------------------------------------------------- */
/* -------------------------------------------------------------------------- */
template
<
UInt
order
>
struct
KelvinIntegrator
{
static
inline
Real
int_exp
(
Real
a
,
Real
z
)
{
return
std
::
exp
(
a
*
z
)
/
a
;
}
static
inline
Real
int_exp_z
(
Real
a
,
Real
z
)
{
const
Real
a_inv
=
1.
/
a
;
return
std
::
exp
(
a
*
z
)
*
(
z
-
a_inv
)
*
a_inv
;
}
static
inline
Real
int_exp_zz
(
Real
a
,
Real
z
)
{
const
Real
a_inv
=
1.
/
a
;
return
std
::
exp
(
a
*
z
)
*
(
z
*
z
-
2
*
a_inv
*
z
+
2
*
a_inv
*
a_inv
)
*
a_inv
;
}
template
<
bool
upper
>
static
inline
Real
g0
(
Real
z
)
{
return
std
::
exp
(
-
sign
(
upper
)
*
z
);
}
template
<
bool
upper
>
static
inline
Real
int_g0
(
const
thrust
::
pair
<
Real
,
Real
>&
bounds
,
Real
alpha
)
{
alpha
*=
-
sign
(
upper
);
return
int_exp
(
alpha
,
bounds
.
second
)
-
int_exp
(
alpha
,
bounds
.
first
);
}
template
<
bool
upper
>
static
inline
Real
int_g0_z
(
const
thrust
::
pair
<
Real
,
Real
>&
bounds
,
Real
alpha
)
{
alpha
*=
-
sign
(
upper
);
return
int_exp_z
(
alpha
,
bounds
.
second
)
-
int_exp_z
(
alpha
,
bounds
.
first
);
}
template
<
bool
upper
>
static
inline
Real
g1
(
Real
z
)
{
return
g0
<
upper
>
(
z
)
*
z
;
}
template
<
bool
upper
>
static
inline
Real
int_g1
(
const
thrust
::
pair
<
Real
,
Real
>&
bounds
,
Real
alpha
)
{
return
alpha
*
int_g0_z
<
upper
>
(
bounds
,
alpha
);
}
template
<
bool
upper
>
static
inline
Real
int_g1_z
(
const
thrust
::
pair
<
Real
,
Real
>&
bounds
,
Real
alpha
)
{
return
alpha
*
(
int_exp_zz
(
-
alpha
*
sign
(
upper
),
bounds
.
second
)
-
int_exp_zz
(
-
alpha
*
sign
(
upper
),
bounds
.
first
));
}
template
<
Int
yj_xi
,
UInt
shape_num
,
UInt
dim_q
,
typename
OutType
,
typename
InType
,
typename
KelvinType
>
static
inline
void
integrate
(
OutType
&&
out
,
InType
&&
in
,
KelvinType
&&
kelvin
,
VectorProxy
<
const
Real
,
dim_q
>&
q_
,
Real
dij
,
Real
dl
);
};
// Integrate against a constant shape function
// template <>
// template <Int yj_xi, UInt shape_num, UInt dim_q, typename OutType,
// typename InType, typename KelvinType>
// inline void KelvinIntegrator<0>::integrate(OutType&& out, InType&& in,
// KelvinType&& kelvin,
// VectorProxy<const Real, dim_q>&
// q_,
// Real dij, Real dl) {
// const auto q = q_.l2norm();
// const auto alpha = q * dl, beta = q * dij;
// // Integrate on [-1, 0]
// const auto bounds_low = thrust::make_pair(-0.5, 0);
// constexpr bool domain_low = yj_xi > 0;
// auto g0_mul{kelvin.template applyU0<domain_low>(q_, in)};
// auto g1_mul{kelvin.template applyU1<domain_low>(q_, in)};
// decltype(g0_mul) res;
// res = g0_mul * (g0<domain_low>(beta) * int_g0<domain_low>(bounds_low,
// alpha));
// res +=
// g1_mul * (g1<domain_low>(beta) * int_g0<domain_low>(bounds_low, alpha)
// +
// g0<domain_low>(beta) * int_g1<domain_low>(bounds_low,
// alpha));
// // Integrate on [0, 1]
// const auto bounds_up = thrust::make_pair(0, 0.5);
// constexpr bool domain_up = yj_xi >= 0;
// // Avoid computing
// if (yj_xi == 0) {
// g0_mul = kelvin.template applyU0<true>(q_, in);
// g1_mul = kelvin.template applyU1<true>(q_, in);
// }
// res += g0_mul * (g0<domain_up>(beta) * int_g0<domain_up>(bounds_up,
// alpha));
// res += g1_mul * (g1<domain_up>(beta) * int_g0<domain_up>(bounds_up, alpha)
// +
// g0<domain_up>(beta) * int_g1<domain_up>(bounds_up,
// alpha));
// res *= dl;
// out += res;
// }
namespace
detail
{
template
<
UInt
shape_num
>
struct
LinearShape
{
using
integrator
=
KelvinIntegrator
<
1
>
;
template
<
bool
domain
,
typename
T
>
static
inline
T
integrateG0
(
T
&
g0_mul
,
Real
alpha
,
Real
beta
);
template
<
bool
domain
,
typename
T
>
static
inline
T
integrateG1
(
T
&
g1_mul
,
Real
alpha
,
Real
beta
);
};
template
<
UInt
shape_num
>
template
<
bool
domain
,
typename
T
>
inline
T
LinearShape
<
shape_num
>::
integrateG0
(
T
&
g0_mul
,
Real
alpha
,
Real
beta
)
{
constexpr
Int
mul
=
(
shape_num
==
0
)
?
-
1
:
1
;
const
auto
bounds
=
thrust
::
make_pair
(
std
::
min
(
-
mul
,
0
),
std
::
max
(
-
mul
,
0
));
return
g0_mul
*
(
integrator
::
g0
<
domain
>
(
beta
)
*
(
integrator
::
int_g0
<
domain
>
(
bounds
,
alpha
)
+
mul
*
integrator
::
int_g0_z
<
domain
>
(
bounds
,
alpha
)));
}
template
<
UInt
shape_num
>
template
<
bool
domain
,
typename
T
>
inline
T
LinearShape
<
shape_num
>::
integrateG1
(
T
&
g1_mul
,
Real
alpha
,
Real
beta
)
{
constexpr
Int
mul
=
(
shape_num
==
0
)
?
-
1
:
1
;
const
auto
bounds
=
thrust
::
make_pair
(
std
::
min
(
-
mul
,
0
),
std
::
max
(
-
mul
,
0
));
return
g1_mul
*
(
integrator
::
g1
<
domain
>
(
beta
)
*
(
integrator
::
int_g0
<
domain
>
(
bounds
,
alpha
)
+
mul
*
integrator
::
int_g0_z
<
domain
>
(
bounds
,
alpha
))
+
integrator
::
g0
<
domain
>
(
beta
)
*
(
integrator
::
int_g1
<
domain
>
(
bounds
,
alpha
)
+
mul
*
integrator
::
int_g1_z
<
domain
>
(
bounds
,
alpha
)));
}
}
// Integrate against a linear shape function
template
<>
template
<
Int
yj_xi
,
UInt
shape_num
,
UInt
dim_q
,
typename
OutType
,
typename
InType
,
typename
KelvinType
>
inline
void
KelvinIntegrator
<
1
>::
integrate
(
OutType
&&
out
,
InType
&&
in
,
KelvinType
&&
kelvin
,
VectorProxy
<
const
Real
,
dim_q
>&
q_
,
Real
dij
,
Real
dl
)
{
const
auto
q
=
q_
.
l2norm
();
const
auto
alpha
=
q
*
dl
,
beta
=
q
*
dij
;
constexpr
bool
domain
=
(
shape_num
==
0
)
?
yj_xi
>=
0
:
yj_xi
>
0
;
auto
g0_mul
=
kelvin
.
template
applyU0
<
domain
>
(
q_
,
in
);
auto
g1_mul
=
kelvin
.
template
applyU1
<
domain
>
(
q_
,
in
);
decltype
(
g0_mul
)
res
;
res
=
detail
::
LinearShape
<
shape_num
>::
template
integrateG0
<
domain
>
(
g0_mul
,
alpha
,
beta
);
res
+=
detail
::
LinearShape
<
shape_num
>::
template
integrateG1
<
domain
>
(
g1_mul
,
alpha
,
beta
);
res
*=
dl
;
out
+=
res
;
}
/* -------------------------------------------------------------------------- */
/// Class for the Boussinesq tensor
template
<
UInt
dim
,
UInt
derivative_order
>
class
Boussinesq
;
template
<>
class
Boussinesq
<
3
,
0
>
{
protected
:
static
constexpr
UInt
dim
=
3
;
public
:
/// Constructor
Boussinesq
(
Real
mu
,
Real
nu
)
:
mu
(
mu
),
nu
(
nu
),
lambda
(
2
*
mu
*
nu
/
(
1
-
2
*
nu
))
{}
template
<
typename
ST
>
inline
Vector
<
Complex
,
dim
>
applyU0
(
StaticVector
<
Complex
,
ST
,
dim
>&
t
,
VectorProxy
<
const
Real
,
dim
-
1
>&
q
)
const
;
template
<
typename
ST
>
inline
Vector
<
Complex
,
dim
>
applyU1
(
StaticVector
<
Complex
,
ST
,
dim
>&
t
,
VectorProxy
<
const
Real
,
dim
-
1
>&
q
)
const
;
protected
:
const
Real
mu
,
nu
,
lambda
;
};
template
<
typename
ST
>
inline
Vector
<
Complex
,
Boussinesq
<
3
,
0
>::
dim
>
Boussinesq
<
3
,
0
>::
applyU0
(
StaticVector
<
Complex
,
ST
,
dim
>&
t
,
VectorProxy
<
const
Real
,
dim
-
1
>&
q
)
const
{
// anticipating K(-q) convolution (no change for d2 cause d2xd2 not sensitive)
auto
d
=
computeD
<
true
>
(
q
),
dbar
=
computeD
<
false
>
(
q
),
d2
=
computeD2
(
q
);
auto
res
=
2.
*
t
;
res
+=
((
1.
-
2.
*
nu
)
*
dbar
.
dot
(
t
))
*
d
;
res
+=
d2
.
dot
(
t
)
*
d2
;
res
(
2
)
-=
t
(
2
);
res
*=
(
1.
/
(
2.
*
mu
*
q
.
l2norm
()));
return
res
;
}
template
<
typename
ST
>
inline
Vector
<
Complex
,
Boussinesq
<
3
,
0
>::
dim
>
Boussinesq
<
3
,
0
>::
applyU1
(
StaticVector
<
Complex
,
ST
,
dim
>&
t
,
VectorProxy
<
const
Real
,
dim
-
1
>&
q
)
const
{
// anticipating K(-q) convolution (no change for d2 cause d2xd2 not sensitive)
auto
dbar
=
computeD
<
false
>
(
q
);
auto
res
=
(
dbar
.
dot
(
t
)
/
(
2.
*
mu
*
q
.
l2norm
()))
*
dbar
;
return
res
;
}
/* -------------------------------------------------------------------------- */
template
<>
class
Boussinesq
<
3
,
1
>
:
protected
Boussinesq
<
3
,
0
>
{
protected
:
using
parent
=
Boussinesq
<
3
,
0
>
;
static
const
UInt
dim
=
parent
::
dim
;
public
:
using
parent
::
parent
;
template
<
typename
ST
>
inline
Matrix
<
Complex
,
dim
,
dim
>
applyU0
(
StaticVector
<
Complex
,
ST
,
dim
>&
t
,
VectorProxy
<
const
Real
,
dim
-
1
>&
q
)
const
;
template
<
typename
ST
>
inline
Matrix
<
Complex
,
dim
,
dim
>
applyU1
(
StaticVector
<
Complex
,
ST
,
dim
>&
t
,
VectorProxy
<
const
Real
,
dim
-
1
>&
q
)
const
;
};
template
<
typename
ST
>
inline
Matrix
<
Complex
,
Boussinesq
<
3
,
1
>::
dim
,
Boussinesq
<
3
,
1
>::
dim
>
Boussinesq
<
3
,
1
>::
applyU0
(
StaticVector
<
Complex
,
ST
,
dim
>&
t
,
VectorProxy
<
const
Real
,
dim
-
1
>&
q
)
const
{
// taking into account K(-q) convolution
auto
dbar
=
computeD
<
false
>
(
q
);
Vector
<
Real
,
dim
>
e3
{{{
0
,
0
,
1
}}};
auto
res
=
outer
(
dbar
,
Boussinesq
<
3
,
0
>::
applyU0
(
t
,
q
));
res
*=
-
1.
;
res
+=
outer
(
e3
,
Boussinesq
<
3
,
0
>::
applyU1
(
t
,
q
));
res
*=
q
.
l2norm
();
return
res
;
}
template
<
typename
ST
>
inline
Matrix
<
Complex
,
Boussinesq
<
3
,
1
>::
dim
,
Boussinesq
<
3
,
1
>::
dim
>
Boussinesq
<
3
,
1
>::
applyU1
(
StaticVector
<
Complex
,
ST
,
dim
>&
t
,
VectorProxy
<
const
Real
,
dim
-
1
>&
q
)
const
{
// taking into account K(-q) convolution
auto
dbar
=
computeD
<
false
>
(
q
);
auto
res
=
outer
(
dbar
,
Boussinesq
<
3
,
0
>::
applyU1
(
t
,
q
));
res
*=
-
q
.
l2norm
();
return
res
;
}
}
// namespace influence
/* -------------------------------------------------------------------------- */
__END_TAMAAS__
#endif
// INFLUENCE_HH
Event Timeline
Log In to Comment