Skip to content

Commit 813f3a6

Browse files
mohd-afeef-badrimohd-afeef-badri
committed
refactor Dirichlet/Source-Term/Traction
1 parent f4a24ec commit 813f3a6

23 files changed

+726
-550
lines changed

‎femutils/ArcaneFemFunctions.h‎

Lines changed: 80 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -1136,7 +1136,7 @@ class ArcaneFemFunctions
11361136
{
11371137
public:
11381138

1139-
static inline void applyDirichletToNodeGroupRhsOnly(const Int32 dof_index, Real value, const IndexedNodeDoFConnectivityView& node_dof, DoFLinearSystem& linear_system, VariableDoFReal& rhs_values, NodeGroup& node_group)
1139+
static inline void applyDirichletToNodeGroupRhsOnly(const Int32 dof_index, Real value, const IndexedNodeDoFConnectivityView& node_dof, VariableDoFReal& rhs_values, NodeGroup& node_group)
11401140
{
11411141
ENUMERATE_ (Node, inode, node_group) {
11421142
Node node = *inode;
@@ -1272,6 +1272,85 @@ class ArcaneFemFunctions
12721272
}
12731273
}
12741274
}
1275+
1276+
/*---------------------------------------------------------------------------*/
1277+
/**
1278+
* @brief Applies Dirichlet boundary conditions to RHS.
1279+
*
1280+
* Updates the RHS vector to enforce Dirichlet conditions.
1281+
*
1282+
* - For RHS vector `𝐛`, the Dirichlet DOF term is scaled by `𝑃`.
1283+
*
1284+
* @param [IN] bs : Boundary condition values.
1285+
* @param [IN] node_dof : DOF connectivity view.
1286+
* @param [IN] node_coord : Node coordinates.
1287+
* @param [OUT] rhs_values RHS : RHS values to update.
1288+
*/
1289+
/*---------------------------------------------------------------------------*/
1290+
static inline void applyDirichletToRhs(BC::IDirichletBoundaryCondition* bs, const IndexedNodeDoFConnectivityView& node_dof, VariableDoFReal& rhs_values)
1291+
{
1292+
FaceGroup face_group = bs->getSurface();
1293+
NodeGroup node_group = face_group.nodeGroup();
1294+
const StringConstArrayView u_dirichlet_string = bs->getValue();
1295+
for (Int32 dof_index = 0; dof_index < u_dirichlet_string.size(); ++dof_index) {
1296+
if (u_dirichlet_string[dof_index] != "NULL") {
1297+
Real value = std::stod(u_dirichlet_string[dof_index].localstr());
1298+
if (bs->getEnforceDirichletMethod() == "Penalty") {
1299+
Real penalty = bs->getPenalty();
1300+
value = value * penalty;
1301+
ArcaneFemFunctions::BoundaryConditionsHelpers::applyDirichletToNodeGroupRhsOnly(dof_index, value, node_dof, rhs_values, node_group);
1302+
}
1303+
else if (bs->getEnforceDirichletMethod() == "RowElimination") {
1304+
ArcaneFemFunctions::BoundaryConditionsHelpers::applyDirichletToNodeGroupRhsOnly(dof_index, value, node_dof, rhs_values, node_group);
1305+
}
1306+
else if (bs->getEnforceDirichletMethod() == "RowColumnElimination") {
1307+
ArcaneFemFunctions::BoundaryConditionsHelpers::applyDirichletToNodeGroupRhsOnly(dof_index, value, node_dof, rhs_values, node_group);
1308+
}
1309+
else {
1310+
ARCANE_FATAL("Unknown Dirichlet method");
1311+
}
1312+
}
1313+
}
1314+
}
1315+
1316+
/*---------------------------------------------------------------------------*/
1317+
/**
1318+
* @brief Applies Point Dirichlet boundary conditions to RHS.
1319+
*
1320+
* Updates the RHS vector to enforce the Dirichlet.
1321+
*
1322+
* - For RHS vector `𝐛`, the Dirichlet DOF term is scaled by `𝑃`.
1323+
*
1324+
* @param [IN] bs : Boundary condition values.
1325+
* @param [IN] node_dof : DOF connectivity view.
1326+
* @param [IN] node_coord : Node coordinates.
1327+
* @param [OUT] rhs_values RHS : RHS values to update.
1328+
*/
1329+
/*---------------------------------------------------------------------------*/
1330+
static inline void applyPointDirichletToRhs(BC::IDirichletPointCondition* bs, const IndexedNodeDoFConnectivityView& node_dof, VariableDoFReal& rhs_values)
1331+
{
1332+
NodeGroup node_group = bs->getNode();
1333+
const StringConstArrayView u_dirichlet_string = bs->getValue();
1334+
for (Int32 dof_index = 0; dof_index < u_dirichlet_string.size(); ++dof_index) {
1335+
if (u_dirichlet_string[dof_index] != "NULL") {
1336+
Real value = std::stod(u_dirichlet_string[dof_index].localstr());
1337+
if (bs->getEnforceDirichletMethod() == "Penalty") {
1338+
Real penalty = bs->getPenalty();
1339+
value = value * penalty;
1340+
ArcaneFemFunctions::BoundaryConditionsHelpers::applyDirichletToNodeGroupRhsOnly(dof_index, value, node_dof, rhs_values, node_group);
1341+
}
1342+
else if (bs->getEnforceDirichletMethod() == "RowElimination") {
1343+
ArcaneFemFunctions::BoundaryConditionsHelpers::applyDirichletToNodeGroupRhsOnly(dof_index, value, node_dof, rhs_values, node_group);
1344+
}
1345+
else if (bs->getEnforceDirichletMethod() == "RowColumnElimination") {
1346+
ArcaneFemFunctions::BoundaryConditionsHelpers::applyDirichletToNodeGroupRhsOnly(dof_index, value, node_dof, rhs_values, node_group);
1347+
}
1348+
else {
1349+
ARCANE_FATAL("Unknown Dirichlet method");
1350+
}
1351+
}
1352+
}
1353+
}
12751354
};
12761355

12771356
/*---------------------------------------------------------------------------*/

‎femutils/ArcaneFemFunctionsGpu.h‎

Lines changed: 102 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -430,6 +430,31 @@ namespace BoundaryConditionsHelpers
430430
/*---------------------------------------------------------------------------*/
431431
/*---------------------------------------------------------------------------*/
432432

433+
inline void applyDirichletToNodeGroupRhsOnly(const Int32 dof_index, Real value, Accelerator::RunQueue* queue, IMesh* mesh, DoFLinearSystem& linear_system, const FemDoFsOnNodes& dofs_on_nodes, NodeGroup& node_group)
434+
{
435+
ARCANE_CHECK_PTR(queue);
436+
ARCANE_CHECK_PTR(mesh);
437+
438+
NodeInfoListView nodes_infos(mesh->nodeFamily());
439+
auto node_dof(dofs_on_nodes.nodeDoFConnectivityView());
440+
441+
auto command = Accelerator::makeCommand(queue);
442+
auto in_out_forced_info = Accelerator::viewInOut(command, linear_system.getForcedInfo());
443+
auto in_out_rhs_variable = Accelerator::viewInOut(command, linear_system.rhsVariable());
444+
445+
command << RUNCOMMAND_ENUMERATE(NodeLocalId, node_lid, node_group)
446+
{
447+
if (nodes_infos.isOwn(node_lid)) {
448+
DoFLocalId dof_id = node_dof.dofId(node_lid, dof_index);
449+
in_out_forced_info[dof_id] = true;
450+
in_out_rhs_variable[dof_id] = value;
451+
}
452+
};
453+
}
454+
455+
/*---------------------------------------------------------------------------*/
456+
/*---------------------------------------------------------------------------*/
457+
433458
inline void applyDirichletToNodeGroupViaPenalty(const Int32 dof_index, Real value, Real penalty, Accelerator::RunQueue* queue, IMesh* mesh, DoFLinearSystem& linear_system, const FemDoFsOnNodes& dofs_on_nodes, NodeGroup& node_group)
434459
{
435460
ARCANE_CHECK_PTR(queue);
@@ -564,6 +589,83 @@ namespace BoundaryConditions
564589
}
565590
}
566591

592+
/*---------------------------------------------------------------------------*/
593+
/**
594+
* @brief Applies Dirichlet boundary conditions to RHS.
595+
*
596+
* - For RHS vector `𝐛`, the Dirichlet DOF term is scaled by `𝑃`.
597+
*/
598+
/*---------------------------------------------------------------------------*/
599+
600+
inline void applyDirichletToRhs(BC::IDirichletBoundaryCondition* bs, const FemDoFsOnNodes& dofs_on_nodes, DoFLinearSystem& linear_system, IMesh* mesh, Accelerator::RunQueue* queue)
601+
{
602+
ARCANE_CHECK_PTR(bs);
603+
604+
FaceGroup face_group = bs->getSurface();
605+
NodeGroup node_group = face_group.nodeGroup();
606+
607+
const StringConstArrayView u_dirichlet_string = bs->getValue();
608+
609+
for (Int32 dof_index = 0; dof_index < u_dirichlet_string.size(); ++dof_index) {
610+
if (u_dirichlet_string[dof_index] != "NULL") {
611+
612+
Real value = std::stod(u_dirichlet_string[dof_index].localstr());
613+
614+
if (bs->getEnforceDirichletMethod() == "Penalty") {
615+
Real penalty = bs->getPenalty();
616+
value = value * penalty;
617+
BoundaryConditionsHelpers::applyDirichletToNodeGroupRhsOnly(dof_index, value, queue, mesh, linear_system, dofs_on_nodes, node_group);
618+
}
619+
else if (bs->getEnforceDirichletMethod() == "RowElimination") {
620+
BoundaryConditionsHelpers::applyDirichletToNodeGroupRhsOnly(dof_index, value, queue, mesh, linear_system, dofs_on_nodes, node_group);
621+
}
622+
else if (bs->getEnforceDirichletMethod() == "RowColumnElimination") {
623+
ARCANE_THROW(Arccore::NotImplementedException, "RowColumnElimination is not supported.");
624+
}
625+
else {
626+
ARCANE_FATAL("Unknown method to enforce Dirichlet BC: '{0}'", bs->getEnforceDirichletMethod());
627+
}
628+
}
629+
}
630+
}
631+
632+
/*---------------------------------------------------------------------------*/
633+
/**
634+
* @brief Applies Point Dirichlet boundary conditions to RHS.
635+
*
636+
* - For RHS vector `𝐛`, the Dirichlet DOF term is scaled by `𝑃`.
637+
*/
638+
/*---------------------------------------------------------------------------*/
639+
640+
inline void applyPointDirichletToRhs(BC::IDirichletPointCondition* bs, const FemDoFsOnNodes& dofs_on_nodes, DoFLinearSystem& linear_system, IMesh* mesh, Accelerator::RunQueue* queue)
641+
{
642+
ARCANE_CHECK_PTR(bs);
643+
NodeGroup node_group = bs->getNode();
644+
645+
const StringConstArrayView u_dirichlet_str = bs->getValue();
646+
647+
for (Int32 dof_index = 0; dof_index < u_dirichlet_str.size(); ++dof_index) {
648+
if (u_dirichlet_str[dof_index] != "NULL") {
649+
Real value = std::stod(u_dirichlet_str[dof_index].localstr());
650+
651+
if (bs->getEnforceDirichletMethod() == "Penalty"){
652+
Real penalty = bs->getPenalty();
653+
value = value * penalty;
654+
BoundaryConditionsHelpers::applyDirichletToNodeGroupRhsOnly(dof_index, value, queue, mesh, linear_system, dofs_on_nodes, node_group);
655+
}
656+
else if (bs->getEnforceDirichletMethod() == "RowElimination") {
657+
BoundaryConditionsHelpers::applyDirichletToNodeGroupRhsOnly(dof_index, value, queue, mesh, linear_system, dofs_on_nodes, node_group);
658+
}
659+
else if (bs->getEnforceDirichletMethod() == "RowColumnElimination") {
660+
ARCANE_THROW(Arccore::NotImplementedException, "RowColumnElimination is not supported.");
661+
}
662+
else {
663+
ARCANE_FATAL("Unknown method to enforce Dirichlet BC: '{0}'", bs->getEnforceDirichletMethod());
664+
}
665+
}
666+
}
667+
}
668+
567669
}; // namespace BoundaryConditions
568670

569671
class BoundaryConditions2D

‎modules/elastodynamics/CMakeLists.txt‎

Lines changed: 28 additions & 9 deletions
Original file line numberDiff line numberDiff line change
@@ -40,40 +40,55 @@ enable_testing()
4040

4141
if(FEMUTILS_HAS_SOLVER_BACKEND_PETSC)
4242
add_test(NAME [elastodynamics]2D_bar COMMAND Elastodynamics
43-
-A,//fem/petsc-flags=-ksp_monitor
43+
-A,//fem/petsc-flags=-mat_type\ baij\ -mat_block_size\ 2\ -ksp_view\ -ksp_monitor
4444
inputs/bar.arc)
4545

4646
add_test(NAME [elastodynamics]2D_bar_RowElimination COMMAND Elastodynamics
47-
-A,//fem/enforce-Dirichlet-method=RowElimination
47+
-A,//fem/boundary-conditions/dirichlet[1]/enforce-Dirichlet-method=RowElimination
48+
-A,//fem/petsc-flags=-mat_type\ baij\ -mat_block_size\ 2\ -ksp_view\ -ksp_monitor\ -ksp_type\ gmres\ -pc_type\ bjacobi\ -ksp_rtol\ 1e-10
4849
inputs/bar.arc)
4950

5051
add_test(NAME [elastodynamics]2D_bar_RowColumnElimination COMMAND Elastodynamics
51-
-A,//fem/enforce-Dirichlet-method=RowColumnElimination
52+
-A,//fem/boundary-conditions/dirichlet[1]/enforce-Dirichlet-method=RowColumnElimination
53+
-A,//fem/petsc-flags=-mat_type\ baij\ -mat_block_size\ 2\ -ksp_view\ -ksp_monitor\ -ksp_type\ gmres\ -pc_type\ bjacobi
5254
inputs/bar.arc)
5355

5456
add_test(NAME [elastodynamics]3D_bar COMMAND Elastodynamics
57+
-A,//fem/petsc-flags=-mat_type\ baij\ -mat_block_size\ 3\ -ksp_view\ -ksp_monitor
5558
inputs/bar.3D.arc)
5659

5760
add_test(NAME [elastodynamics]3D_bar_RowColumnElimination_bjacobi COMMAND Elastodynamics
58-
-A,//fem/enforce-Dirichlet-method=RowColumnElimination
59-
-A,//fem/petsc-flags=-mat_type\ baij\ -mat_block_size\ 3\ -ksp_monitor\ -ksp_type\ gmres\ -pc_type\ bjacobi\ -ksp_rtol\ 1e-8
61+
-A,//fem/boundary-conditions/dirichlet[1]/enforce-Dirichlet-method=RowColumnElimination
62+
-A,//fem/petsc-flags=-mat_type\ baij\ -mat_block_size\ 3\ -ksp_monitor\ -ksp_type\ gmres\ -pc_type\ bjacobi\ -ksp_rtol\ 1e-10
6063
inputs/bar.3D.arc)
6164

6265
add_test(NAME [elastodynamics]2D_transient_traction COMMAND Elastodynamics
6366
inputs/bar.transient-traction.arc)
6467

6568
add_test(NAME [elastodynamics]Dirichlet_pointBc_Penalty COMMAND Elastodynamics
66-
-A,//fem/enforce-Dirichlet-method=Penalty
6769
-A,//fem/petsc-flags=-mat_type\ baij\ -mat_block_size\ 2\ -ksp_view\ -ksp_monitor\ -ksp_type\ cg\ -pc_type\ bjacobi\ -ksp_rtol\ 1e-15
6870
inputs/semi-circle.pointBC.arc)
6971

7072
add_test(NAME [elastodynamics]2D_pointBc_RowElimination COMMAND Elastodynamics
71-
-A,//fem/enforce-Dirichlet-method=RowElimination
72-
-A,//fem/petsc-flags=-mat_type\ baij\ -mat_block_size\ 2\ -ksp_monitor\ -ksp_type\ gmres\ -pc_type\ bjacobi\ -ksp_rtol\ 1e-15
73+
-A,//fem/boundary-conditions/dirichlet[1]/enforce-Dirichlet-method=RowElimination
74+
-A,//fem/boundary-conditions/dirichlet-point[1]/enforce-Dirichlet-method=RowElimination
75+
-A,//fem/petsc-flags=-mat_type\ baij\ -mat_block_size\ 2\ -ksp_monitor\ -ksp_type\ gmres\ -pc_type\ bjacobi\ -ksp_rtol\ 1e-10
7376
inputs/semi-circle.pointBC.arc)
7477

7578
add_test(NAME [elastodynamics]3D_truncated-cube_pointBc_Penalty COMMAND Elastodynamics
76-
-A,//fem/enforce-Dirichlet-method=Penalty
79+
-A,//fem/petsc-flags=-mat_type\ baij\ -mat_block_size\ 3\ -ksp_view\ -ksp_monitor\ -ksp_type\ cg\ -pc_type\ bjacobi\ -ksp_rtol\ 1e-15
80+
inputs/truncated-cube.pointBC.arc)
81+
82+
add_test(NAME [elastodynamics]3D_truncated-cube_pointBc_RowColElim COMMAND Elastodynamics
83+
-A,//fem/boundary-conditions/dirichlet[1]/enforce-Dirichlet-method=RowColumnElimination
84+
-A,//fem/boundary-conditions/dirichlet-point[1]/enforce-Dirichlet-method=RowColumnElimination
85+
-A,//fem/petsc-flags=-mat_type\ baij\ -mat_block_size\ 3\ -ksp_view\ -ksp_monitor\ -ksp_type\ cg\ -pc_type\ bjacobi\ -ksp_rtol\ 1e-15
86+
inputs/truncated-cube.pointBC.arc)
87+
88+
add_test(NAME [elastodynamics]3D_truncated-cube_pointBc_RowElim COMMAND Elastodynamics
89+
-A,//fem/boundary-conditions/dirichlet[1]/enforce-Dirichlet-method=RowElimination
90+
-A,//fem/boundary-conditions/dirichlet-point[1]/enforce-Dirichlet-method=RowElimination
91+
-A,//fem/petsc-flags=-mat_type\ baij\ -mat_block_size\ 3\ -ksp_view\ -ksp_monitor\ -ksp_type\ gmres\ -pc_type\ bjacobi\ -ksp_rtol\ 1e-15
7792
inputs/truncated-cube.pointBC.arc)
7893

7994
add_test(NAME [elastodynamics]3D_transient_traction COMMAND Elastodynamics
@@ -93,18 +108,22 @@ if(FEMUTILS_HAS_SOLVER_BACKEND_PETSC)
93108

94109
arcanefem_add_gpu_test(NAME [elastodynamics]2D_bar_petsc_bsr COMMAND ./Elastodynamics ARGS
95110
-A,//fem/matrix-format=BSR
111+
-A,//fem/petsc-flags=-ksp_view\ -ksp_monitor\ -ksp_type\ gmres\ -pc_type\ bjacobi\ -ksp_rtol\ 1e-10
96112
inputs/bar.arc)
97113

98114
arcanefem_add_gpu_test(NAME [elastodynamics]3D_bar_petsc_bsr COMMAND ./Elastodynamics ARGS
99115
-A,//fem/matrix-format=BSR
116+
-A,//fem/petsc-flags=-mat_type\ baij\ -mat_block_size\ 3\ -ksp_view\ -ksp_monitor
100117
inputs/bar.3D.arc)
101118

102119
arcanefem_add_gpu_test(NAME [elastodynamics]2D_bar_petsc_af-bsr COMMAND ./Elastodynamics ARGS
103120
-A,//fem/matrix-format=AF-BSR
121+
-A,//fem/petsc-flags=-mat_type\ baij\ -mat_block_size\ 2\ -ksp_view\ -ksp_monitor
104122
inputs/bar.arc)
105123

106124
arcanefem_add_gpu_test(NAME [elastodynamics]3D_bar_petsc_af-bsr COMMAND ./Elastodynamics ARGS
107125
-A,//fem/matrix-format=AF-BSR
126+
-A,//fem/petsc-flags=-mat_type\ baij\ -mat_block_size\ 3\ -ksp_view\ -ksp_monitor
108127
inputs/bar.3D.arc)
109128

110129
endif()

0 commit comments

Comments
 (0)