Page Menu
Home
c4science
Search
Configure Global Search
Log In
Files
F77507299
bridging_atomic_continuum.cc
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
Wed, Aug 14, 21:31
Size
12 KB
Mime Type
text/x-c++
Expires
Fri, Aug 16, 21:31 (2 d)
Engine
blob
Format
Raw Data
Handle
19880484
Attached To
rLIBMULTISCALE LibMultiScale
bridging_atomic_continuum.cc
View Options
/**
* @file bridging_atomic_continuum.cc
*
* @author Guillaume Anciaux <guillaume.anciaux@epfl.ch>
*
* @date Mon Oct 28 19:23:14 2013
*
* @brief Bridging object between atomistic and finite elements
*
* @section LICENSE
*
* Copyright (©) 2010-2011 EPFL (Ecole Polytechnique Fédérale de Lausanne)
* Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
*
* LibMultiScale 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.
*
* LibMultiScale 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 LibMultiScale. If not, see <http://www.gnu.org/licenses/>.
*
*/
/* -------------------------------------------------------------------------- */
//#define TIMER
/* -------------------------------------------------------------------------- */
#include "bridging_atomic_continuum.hh"
#include "compute_extract.hh"
#include "lib_bridging.hh"
#include "lm_common.hh"
#include "reference_manager.hh"
#include "trace_atom.hh"
/* -------------------------------------------------------------------------- */
__BEGIN_LIBMULTISCALE__
template
<
typename
DomainA
,
typename
DomainC
>
BridgingAtomicContinuum
<
DomainA
,
DomainC
>::
BridgingAtomicContinuum
(
const
std
::
string
&
name
,
DomainA
&
A
,
DomainC
&
C
)
:
LMObject
(
name
),
Bridging
<
typename
DomainA
::
ContainerAtoms
,
typename
DomainC
::
ContainerPoints
>
(
name
,
A
.
getContainer
(),
C
.
getContainer
()),
domA
(
A
),
domC
(
C
),
comm_group_continuum
(
domC
.
getCommGroup
()),
comm_group_atomic
(
domA
.
getCommGroup
())
{
if
(
domA
.
getContainer
().
hasRefManager
())
domA
.
getContainer
().
getRefManager
()
->
addSubSet
(
this
->
pointList
);
}
/* -------------------------------------------------------------------------- */
template
<
typename
DomainA
,
typename
DomainC
>
BridgingAtomicContinuum
<
DomainA
,
DomainC
>::~
BridgingAtomicContinuum
()
{}
/* -------------------------------------------------------------------------- */
template
<
typename
DomainA
,
typename
DomainC
>
void
BridgingAtomicContinuum
<
DomainA
,
DomainC
>::
clearAll
()
{
Bridging
<
typename
DomainA
::
ContainerPoints
,
typename
DomainC
::
ContainerPoints
>::
clearAll
();
}
/* -------------------------------------------------------------------------- */
template
<
typename
DomainA
,
typename
DomainC
>
void
BridgingAtomicContinuum
<
DomainA
,
DomainC
>::
updateForMigration
()
{
Bridging
<
typename
DomainA
::
ContainerPoints
,
typename
DomainC
::
ContainerPoints
>::
updateForMigration
();
if
(
this
->
check_coherency
)
{
checkPositionCoherency
();
}
}
/* -------------------------------------------------------------------------- */
template
<
typename
DomainA
,
typename
DomainC
>
void
BridgingAtomicContinuum
<
DomainA
,
DomainC
>::
init
()
{
Bridging
<
typename
DomainA
::
ContainerPoints
,
typename
DomainC
::
ContainerPoints
>::
init
();
DUMP
(
this
->
getID
()
<<
" : attach duo vector"
,
DBG_MESSAGE
);
if
(
this
->
in_group_A
)
{
try
{
this
->
contA
.
getRefManager
()
->
attachObject
(
this
->
getDuoVector
(),
this
->
pointList
);
DUMP
(
this
->
getID
()
<<
" : attaching Parent::duo_vector"
,
DBG_INFO
);
}
catch
(
LibMultiScaleException
&
e
)
{
}
}
MPI_Barrier
(
MPI_COMM_WORLD
);
}
/* -------------------------------------------------------------------------- */
template
<
typename
DomainA
,
typename
DomainC
>
void
BridgingAtomicContinuum
<
DomainA
,
DomainC
>::
averageDOFsOnElements
(
ContainerArray
<
Real
>
&
data_mesh
,
FieldType
field
)
{
LM_FATAL
(
"to adapt"
);
if
(
comm_group_atomic
==
comm_group_continuum
)
{
ComputeExtract
extracted_field
(
"ComputeExtract:"
+
this
->
getID
());
extracted_field
.
setParam
(
"FIELD"
,
field
);
extracted_field
.
build
(
*
this
->
pointList
);
this
->
averageOnElements
(
extracted_field
.
getArray
(),
data_mesh
);
return
;
}
ComputeExtract
extracted_field
(
"ComputeExtract:"
+
this
->
getID
());
extracted_field
.
setParam
(
"FIELD"
,
field
);
if
(
this
->
is_in_atomic
)
extracted_field
.
build
(
*
this
->
pointList
);
LM_TOIMPLEMENT
;
// extracted_field.resize(Parent::local_points * Dim);
/* do communication */
this
->
communicateBufferAtom2Continuum
(
extracted_field
.
getArray
());
if
(
this
->
is_in_continuum
)
this
->
averageOnElements
(
extracted_field
.
getArray
(),
data_mesh
);
}
/* -------------------------------------------------------------------------- */
// template <typename DomainA, typename DomainC,UInt Dim>
// template <FieldType field>
// void BridgingAtomicContinuum<DomainA,DomainC,Dim>::
// buildAtomicDOFsVector(std::vector<Real> & v_out){
// IteratorAtomsSubset it = this->pointList->begin();
// UInt nmax = this->pointList->size();
// UInt i =0;
// Real * v_out_ptr = &v_out[0];
// for (RefAtom at = it.getFirst(); !it.end() ;
// at = it.getNext() , ++i, v_out_ptr+= Dim) {
// LM_ASSERT(i < nmax,"overflow detected");
// at.template getField<field>(v_out_ptr);
// }
// }
// /* --------------------------------------------------------------------------
// */
// template <typename DomainA, typename DomainC,UInt Dim>
// void BridgingAtomicContinuum<DomainA,DomainC,Dim>::
// buildAtomicDOFsVector(std::vector<Real> & v_out,
// FieldType field){
// switch(field){
// case _position0: buildAtomicDOFsVector<_position0 >(v_out); break;
// case _position: buildAtomicDOFsVector<_position >(v_out); break;
// case _displacement: buildAtomicDOFsVector<_displacement>(v_out); break;
// case _velocity: buildAtomicDOFsVector<_velocity >(v_out); break;
// case _force: buildAtomicDOFsVector<_force >(v_out);
// break;
// default: LM_FATAL("not handled field");
// }
// }
/* -------------------------------------------------------------------------- */
template
<
typename
DomainA
,
typename
DomainC
>
void
BridgingAtomicContinuum
<
DomainA
,
DomainC
>::
leastSquareAtomicData
(
ContainerArray
<
Real
>
&
data_mesh
,
ContainerArray
<
Real
>
&
data_atom
)
{
if
(
comm_group_atomic
==
comm_group_continuum
)
{
this
->
solveLeastSquare
(
data_mesh
,
data_atom
);
return
;
}
/* do communication */
this
->
communicateBufferAtom2Continuum
(
data_atom
);
if
(
this
->
is_in_continuum
)
{
this
->
solveLeastSquare
(
data_mesh
,
data_atom
);
}
}
/* -------------------------------------------------------------------------- */
template
<
typename
DomainA
,
typename
DomainC
>
void
BridgingAtomicContinuum
<
DomainA
,
DomainC
>::
leastSquareAtomicDOFs
(
ContainerArray
<
Real
>
&
data_mesh
,
FieldType
field
)
{
LM_FATAL
(
"to adapt"
);
if
(
comm_group_atomic
==
comm_group_continuum
)
{
ComputeExtract
extracted_field
(
"ComputeExtract:"
+
this
->
getID
());
extracted_field
.
setParam
(
"FIELD"
,
field
);
extracted_field
.
build
<
typename
DomainA
::
ContainerSubset
>
(
this
->
pointList
);
this
->
solveLeastSquare
(
data_mesh
,
extracted_field
.
getArray
());
return
;
}
ComputeExtract
extracted_field
(
"ComputeExtract:"
+
this
->
getID
());
extracted_field
.
setParam
(
"FIELD"
,
field
);
if
(
this
->
is_in_atomic
)
extracted_field
.
build
<
typename
DomainA
::
ContainerSubset
>
(
this
->
pointList
);
LM_TOIMPLEMENT
;
// extracted_field.resize(Parent::local_points * Dim);
this
->
leastSquareAtomicData
(
data_mesh
,
extracted_field
.
getArray
());
}
/* -------------------------------------------------------------------------- */
template
<
typename
DomainA
,
typename
DomainC
>
void
BridgingAtomicContinuum
<
DomainA
,
DomainC
>::
buildAtomicDOFsNoiseVector
(
std
::
vector
<
Real
>
&
,
FieldType
field
)
{
ComputeExtract
extracted_field
(
"ComputeExtract:"
+
this
->
getID
());
extracted_field
.
setParam
(
"FIELD"
,
field
);
if
(
this
->
is_in_atomic
)
extracted_field
.
build
(
*
this
->
pointList
);
if
(
this
->
is_in_continuum
)
{
this
->
interpolateAtomicDOFs
(
field
,
extracted_field
.
getArray
());
for
(
UInt
i
=
0
;
i
<
this
->
local_points
*
DomainA
::
ContainerPoints
::
Dim
;
++
i
)
LM_TOIMPLEMENT
;
// extracted_field[i] *= -1;
}
/* do communication */
this
->
synchSumBuffer
(
extracted_field
.
getArray
());
}
/* -------------------------------------------------------------------------- */
template
<
typename
DomainA
,
typename
DomainC
>
void
BridgingAtomicContinuum
<
DomainA
,
DomainC
>::
checkPositionCoherency
()
{
constexpr
UInt
Dim
=
DomainA
::
ContainerPoints
::
Dim
;
if
(
this
->
in_group_B
)
{
// interpolate atomic sites to check everything is ok from this
auto
&&
field
=
domC
.
getField
(
_position0
);
UInt
*
assoc_tmp
=
&
this
->
pointToElement
[
0
];
for
(
UInt
i
=
0
;
i
<
this
->
local_points
;
++
i
)
{
Vector
<
Dim
>
X
;
Vector
<
Dim
>
pos
=
VectorView
<
Dim
>
(
&
this
->
positions
[
Dim
*
i
]);
X
=
this
->
getShapeMatrix
().
interpolate
(
i
,
field
,
assoc_tmp
[
i
]);
for
(
UInt
k
=
0
;
k
<
Dim
;
++
k
)
{
if
(
X
[
k
]
!=
pos
[
k
])
if
(
fabs
(
pos
[
k
]
-
X
[
k
])
>
1e-5
)
{
searchAtomInLocalPool
(
static_cast
<
typename
DomainA
::
ContainerSubset
&>
(
this
->
pointList
),
pos
);
DUMP
(
"mismatch in positions (X), atom "
<<
i
<<
" dist is "
<<
pos
[
k
]
-
X
[
k
]
<<
" : probably a redistribution trouble check "
"migrations
\n
"
<<
"computed position is "
<<
X
[
0
]
<<
" "
<<
X
[
1
]
<<
" "
<<
X
[
2
]
<<
" position was "
<<
pos
[
0
]
<<
" "
<<
pos
[
1
]
<<
" "
<<
pos
[
2
],
DBG_MESSAGE
);
}
}
}
}
this
->
distributeVectorB2A
(
"temp_positions"
,
this
->
positions
);
if
(
this
->
in_group_A
)
{
UInt
i
=
0
;
UInt
global_flag
=
0
;
for
(
auto
&&
at
:
this
->
pointList
)
{
UInt
local_flag
=
0
;
Vector
<
Dim
>
pos
=
VectorView
<
Dim
>
(
&
this
->
positions
[
Dim
*
i
]);
Vector
<
Dim
>
at_pos
=
at
.
position0
();
for
(
UInt
k
=
0
;
k
<
Dim
;
++
k
)
{
if
(
pos
[
k
]
!=
at_pos
[
k
])
if
(
!
local_flag
&&
fabs
(
pos
[
k
]
-
at_pos
[
k
])
>
1e-5
)
{
searchAtomInLocalPool
(
this
->
pointList
,
pos
);
DUMP
(
"mismatch in positions (X), atom "
<<
i
<<
" dist is "
<<
this
->
positions
[
Dim
*
i
]
-
at
.
position0
()[
0
]
<<
" : probably a redistribution trouble check migrations"
<<
std
::
endl
<<
"ref is "
<<
at
<<
"
\n
and send position is "
<<
pos
[
0
]
<<
" "
<<
pos
[
1
]
<<
" "
<<
pos
[
2
],
DBG_MESSAGE
);
local_flag
=
1
;
global_flag
=
1
;
}
}
++
i
;
}
if
(
global_flag
)
{
this
->
getDuoVector
().
print
(
"buggy-redistrib"
);
LM_FATAL
(
"abort due to migration problem"
);
}
}
DUMP
(
"position coherency OK"
,
DBG_INFO
);
}
/* -------------------------------------------------------------------------- */
template
<
typename
DomainA
,
typename
DomainC
>
template
<
typename
Cont
,
typename
Vec
>
void
BridgingAtomicContinuum
<
DomainA
,
DomainC
>::
searchAtomInLocalPool
(
Cont
&
container
,
Vec
&
x
)
{
UInt
i
=
0
;
for
(
auto
&&
at
:
container
)
{
#ifdef TRACE_ATOM
auto
pos
=
at
.
position0
();
#endif
IF_COORD_EQUALS
(
pos
,
x
)
{
DUMP
(
"actually searched atom is at ref "
<<
at
<<
" bridge index "
<<
i
,
DBG_MESSAGE
);
break
;
}
++
i
;
}
if
(
i
==
container
.
size
())
DUMP
(
"atom "
<<
x
[
0
]
<<
" "
<<
x
[
1
]
<<
" "
<<
x
[
2
]
<<
" "
<<
" not found locally : probably has migrated"
<<
" or is simply lost !"
<<
" NOT GOOD"
,
DBG_MESSAGE
);
i
=
0
;
for
(
auto
&&
at
:
domA
.
getContainer
())
{
#ifdef TRACE_ATOM
auto
pos
=
at
.
position0
();
#endif
IF_COORD_EQUALS
(
pos
,
x
)
{
DUMP
(
"actually searched atom is at ref "
<<
at
<<
" bridge index "
<<
i
,
DBG_MESSAGE
);
break
;
}
++
i
;
}
if
(
i
==
domA
.
getContainer
().
size
())
DUMP
(
"atom not found locally : probably has migrated"
<<
"or is simply lost !"
<<
" NOT GOOD"
,
DBG_MESSAGE
);
}
/* -------------------------------------------------------------------------- */
DECLARE_ATOMIC_CONTINUUM_TEMPLATE
(
BridgingAtomicContinuum
)
__END_LIBMULTISCALE__
Event Timeline
Log In to Comment