Page Menu
Home
c4science
Search
Configure Global Search
Log In
Files
F62040118
kelvin_helper.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
Fri, May 10, 13:07
Size
6 KB
Mime Type
text/x-c++
Expires
Sun, May 12, 13:07 (2 d)
Engine
blob
Format
Raw Data
Handle
17594213
Attached To
rTAMAAS tamaas
kelvin_helper.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 KELVIN_HELPER_HH
#define KELVIN_HELPER_HH
/* -------------------------------------------------------------------------- */
#include "grid.hh"
#include "grid_hermitian.hh"
#include "influence.hh"
#include "model_type.hh"
/* -------------------------------------------------------------------------- */
__BEGIN_TAMAAS__
namespace
detail
{
/// Trait for kelvin local types
template
<
typename
T
>
struct
KelvinTrait
;
template
<
UInt
dim
>
struct
KelvinTrait
<
influence
::
Kelvin
<
dim
,
0
>>
{
using
source_t
=
VectorProxy
<
Complex
,
dim
>
;
using
out_t
=
VectorProxy
<
Complex
,
dim
>
;
};
template
<
UInt
dim
>
struct
KelvinTrait
<
influence
::
Kelvin
<
dim
,
1
>>
{
using
source_t
=
MatrixProxy
<
Complex
,
dim
,
dim
>
;
using
out_t
=
VectorProxy
<
Complex
,
dim
>
;
};
template
<
UInt
dim
>
struct
KelvinTrait
<
influence
::
Kelvin
<
dim
,
2
>>
{
using
source_t
=
MatrixProxy
<
Complex
,
dim
,
dim
>
;
using
out_t
=
MatrixProxy
<
Complex
,
dim
,
dim
>
;
};
template
<
UInt
dim
>
struct
KelvinTrait
<
influence
::
Boussinesq
<
dim
,
0
>>
{
using
source_t
=
VectorProxy
<
Complex
,
dim
>
;
using
out_t
=
VectorProxy
<
Complex
,
dim
>
;
};
template
<
UInt
dim
>
struct
KelvinTrait
<
influence
::
Boussinesq
<
dim
,
1
>>
{
using
source_t
=
VectorProxy
<
Complex
,
dim
>
;
using
out_t
=
MatrixProxy
<
Complex
,
dim
,
dim
>
;
};
/// Helper to apply integral representation on output
template
<
model_type
type
,
typename
T
>
struct
KelvinHelper
{
using
kelvin_t
=
T
;
using
trait
=
model_type_traits
<
type
>
;
static
constexpr
UInt
dim
=
trait
::
dimension
;
static
constexpr
UInt
bdim
=
trait
::
boundary_dimension
;
KelvinHelper
(
const
T
&
kelvin
,
Real
cutoff
)
:
kelvin
(
kelvin
),
cutoff
(
cutoff
)
{}
template
<
Int
yj_xi
,
UInt
shape_num
>
inline
void
applyIntegral
(
GridHermitian
<
Real
,
bdim
>&
out
,
GridHermitian
<
Real
,
bdim
>&
source
,
Real
dij
,
Real
dl
,
const
Grid
<
Real
,
bdim
>&
wavevectors
)
const
{
Loop
::
stridedLoop
(
[
&
](
typename
KelvinTrait
<
T
>::
out_t
&&
out_local
,
typename
KelvinTrait
<
T
>::
source_t
&&
source_local
,
VectorProxy
<
const
Real
,
bdim
>&&
q
)
{
// Cutoff
if
(
-
q
.
l2norm
()
*
std
::
abs
(
dij
)
<
std
::
log
(
cutoff
))
return
;
influence
::
KelvinIntegrator
<
1
>::
integrate
<
yj_xi
,
shape_num
>
(
out_local
,
source_local
,
kelvin
,
q
,
dij
,
dl
);
},
out
,
source
,
wavevectors
);
}
inline
void
applyFreeTerm
(
GridHermitian
<
Real
,
bdim
>&
out
,
GridHermitian
<
Real
,
bdim
>&
source
,
const
Grid
<
Real
,
bdim
>&
wavevectors
)
const
{}
const
T
&
kelvin
;
const
Real
cutoff
;
};
template
<>
inline
void
KelvinHelper
<
model_type
::
volume_2d
,
influence
::
Kelvin
<
3
,
2
>>::
applyFreeTerm
(
GridHermitian
<
Real
,
bdim
>&
out
,
GridHermitian
<
Real
,
bdim
>&
source
,
const
Grid
<
Real
,
bdim
>&
wavevectors
)
const
{
Loop
::
stridedLoop
(
[
&
](
KelvinTrait
<
kelvin_t
>::
out_t
&&
out_local
,
KelvinTrait
<
kelvin_t
>::
source_t
&&
source_local
,
VectorProxy
<
const
Real
,
bdim
>&&
q
)
{
out_local
+=
kelvin
.
applyDiscontinuityTerm
(
q
,
source_local
);
},
out
,
source
,
wavevectors
);
}
/* -------------------------------------------------------------------------- */
template
<
model_type
type
,
typename
T
>
class
LinearElement
{
using
trait
=
model_type_traits
<
type
>
;
static
constexpr
UInt
bdim
=
trait
::
boundary_dimension
;
public
:
LinearElement
(
const
thrust
::
pair
<
UInt
,
Real
>&
integration_node
,
const
T
&
kelvin
,
std
::
vector
<
GridHermitian
<
Real
,
bdim
>>&
sources
,
GridHermitian
<
Real
,
bdim
>&
out
,
const
Grid
<
Real
,
bdim
>&
wavevectors
,
Real
cutoff
=
1e-2
)
:
integration_node
(
integration_node
),
sources
(
sources
),
out
(
out
),
wavevectors
(
wavevectors
),
helper
(
kelvin
,
cutoff
)
{}
private
:
template
<
UInt
shape_num
>
void
integrateNode
(
UInt
element
,
Real
yj
,
Real
el_size
)
{
static_assert
(
shape_num
==
0
or
shape_num
==
1
,
"Shape number is invalid"
);
const
Real
dij
=
yj
-
integration_node
.
second
;
const
UInt
node
=
element
+
shape_num
;
if
(
node
>
integration_node
.
first
)
helper
.
template
applyIntegral
<
1
,
shape_num
>
(
out
,
sources
[
node
],
dij
,
el_size
,
wavevectors
);
else
if
(
node
==
integration_node
.
first
)
{
helper
.
template
applyIntegral
<
0
,
shape_num
>
(
out
,
sources
[
node
],
dij
,
el_size
,
wavevectors
);
if
(
shape_num
==
0
)
// apply free term only once
helper
.
applyFreeTerm
(
out
,
sources
[
node
],
wavevectors
);
}
else
helper
.
template
applyIntegral
<-
1
,
shape_num
>
(
out
,
sources
[
node
],
dij
,
el_size
,
wavevectors
);
}
public
:
void
integrateElement
(
UInt
e
,
const
thrust
::
pair
<
Real
,
Real
>&
positions
)
{
const
Real
dl
=
std
::
abs
(
positions
.
second
-
positions
.
first
);
integrateNode
<
0
>
(
e
,
positions
.
first
,
dl
);
integrateNode
<
1
>
(
e
,
positions
.
second
,
dl
);
}
private
:
thrust
::
pair
<
UInt
,
Real
>
integration_node
;
std
::
vector
<
GridHermitian
<
Real
,
bdim
>>&
sources
;
GridHermitian
<
Real
,
bdim
>&
out
;
const
Grid
<
Real
,
bdim
>&
wavevectors
;
KelvinHelper
<
type
,
T
>
helper
;
};
}
// namespace detail
__END_TAMAAS__
#endif
// KELVIN_HELPER_HH
Event Timeline
Log In to Comment