v0.14.0
Loading...
Searching...
No Matches
MultifieldPlasticity.cpp
Go to the documentation of this file.
1/** \file MultifieldPlasticity.cpp
2 */
3
4/* This file is part of MoFEM.
5 * MoFEM is free software: you can redistribute it and/or modify it under
6 * the terms of the GNU Lesser General Public License as published by the
7 * Free Software Foundation, either version 3 of the License, or (at your
8 * option) any later version.
9 *
10 * MoFEM is distributed in the hope that it will be useful, but WITHOUT
11 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
12 * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
13 * License for more details.
14 *
15 * You should have received a copy of the GNU Lesser General Public
16 * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>. */
17
18#include <MoFEM.hpp>
19using namespace MoFEM;
20using namespace FTensor;
21
24#include <RigidBodies.hpp>
25#include <BlockMatData.hpp>
27#include <BasicFeTools.hpp>
28
29
30namespace MoFEM {
31
34 SmartPetscObj<DM> &dm_plastic) {
36 auto simple = m_field.getInterface<Simple>();
37
38 dm_plastic = createSmartDM(m_field.get_comm(), "DMMOFEM");
39 CHKERR DMMoFEMCreateSubDM(dm_plastic, dm, "PLASTIC_EP_PROBLEM");
40 CHKERR DMAppendOptionsPrefix(dm_plastic, "ep_");
41 CHKERR DMSetFromOptions(dm_plastic);
42
43 // CHKERR m_field.getInterface<ProblemsManager>()->addFieldToEmptyFieldBlocks(
44 // "PLASTIC_EP_PROBLEM", "U", "TAU");
45 CHKERR DMMoFEMSetDestroyProblem(dm_plastic, PETSC_TRUE);
46 CHKERR DMMoFEMAddSubFieldRow(dm_plastic, "U");
47 CHKERR DMMoFEMAddSubFieldRow(dm_plastic, "TAU");
48 //FIXME:
49 if (m_field.check_field("SIGMA"))
50 CHKERR DMMoFEMAddSubFieldRow(dm_plastic, "SIGMA");
51 CHKERR DMMoFEMAddElement(dm_plastic, simple->getDomainFEName().c_str());
52 CHKERR DMMoFEMAddElement(dm_plastic, simple->getBoundaryFEName().c_str());
53 // CHKERR DMMoFEMAddElement(dm_plastic, simple->getSkeletonFEName().c_str());
54 CHKERR DMMoFEMSetSquareProblem(dm_plastic, PETSC_TRUE);
55 CHKERR DMMoFEMSetIsPartitioned(dm_plastic, PETSC_TRUE);
56 CHKERR DMSetUp(dm_plastic);
57
58 if (0) {
61 ->createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>("PLASTIC_EP_PROBLEM",
62 m);
63 // MatView(m, PETSC_VIEWER_STDOUT_WORLD);
64 // MatView(m, PETSC_VIEWER_DRAW_WORLD);
65 // MatView(m, PETSC_VIEWER_MATLAB_WORLD);
66 std::string wait;
67 // std::cin >> wait;
69 ->checkMPIAIJWithArraysMatrixFillIn<PetscGlobalIdx_mi_tag>(
70 "PLASTIC_EP_PROBLEM", 1, 1, 1);
71 }
72
74}
75
78 SmartPetscObj<DM> &dm_plastic_tau) {
80 auto simple = m_field.getInterface<Simple>();
81
82 dm_plastic_tau = createSmartDM(m_field.get_comm(), "DMMOFEM");
83 CHKERR DMMoFEMCreateSubDM(dm_plastic_tau, dm, "PLASTIC_TAU_PROBLEM");
84 CHKERR DMAppendOptionsPrefix(dm_plastic_tau, "tau_");
85 CHKERR DMSetFromOptions(dm_plastic_tau);
86 CHKERR DMMoFEMSetDestroyProblem(dm_plastic_tau, PETSC_TRUE);
87 CHKERR DMMoFEMAddSubFieldRow(dm_plastic_tau, "U");
88 //FIXME:
89 if (m_field.check_field("SIGMA"))
90 CHKERR DMMoFEMAddSubFieldRow(dm_plastic_tau, "SIGMA");
91 CHKERR DMMoFEMAddElement(dm_plastic_tau, simple->getDomainFEName().c_str());
92 CHKERR DMMoFEMAddElement(dm_plastic_tau, simple->getBoundaryFEName().c_str());
93 // CHKERR DMMoFEMAddElement(dm_plastic_tau,
94 // simple->getSkeletonFEName().c_str());
95 CHKERR DMMoFEMSetSquareProblem(dm_plastic_tau, PETSC_TRUE);
96 CHKERR DMMoFEMSetIsPartitioned(dm_plastic_tau, PETSC_TRUE);
97 CHKERR DMSetUp(dm_plastic_tau);
98
100}
101
103 SmartPetscObj<DM> &dm_contact) {
105 auto simple = m_field.getInterface<Simple>();
106
107 dm_contact = createSmartDM(m_field.get_comm(), "DMMOFEM");
108 CHKERR DMMoFEMCreateSubDM(dm_contact, dm, "CONTACT_PROBLEM");
109 CHKERR DMAppendOptionsPrefix(dm_contact, "contact_");
110 CHKERR DMSetFromOptions(dm_contact);
111 CHKERR DMMoFEMSetDestroyProblem(dm_contact, PETSC_TRUE);
112 CHKERR DMMoFEMAddSubFieldRow(dm_contact, "U");
113
114 if (m_field.check_field("TAU"))
115 CHKERR DMMoFEMAddSubFieldRow(dm_contact, "TAU");
116 if (m_field.check_field("EP"))
117 CHKERR DMMoFEMAddSubFieldRow(dm_contact, "EP");
118
119 CHKERR DMMoFEMAddElement(dm_contact, simple->getDomainFEName().c_str());
120 CHKERR DMMoFEMAddElement(dm_contact, simple->getBoundaryFEName().c_str());
121 // CHKERR DMMoFEMAddElement(dm_contact, simple->getSkeletonFEName().c_str());
122 CHKERR DMMoFEMSetSquareProblem(dm_contact, PETSC_TRUE);
123 CHKERR DMMoFEMSetIsPartitioned(dm_contact, PETSC_TRUE);
124 CHKERR DMSetUp(dm_contact);
125
127}
128
129} // namespace MoFEM
void simple(double P1[], double P2[], double P3[], double c[], const int N)
Definition acoustic.cpp:69
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
FTensor::Index< 'm', SPACE_DIM > m
PetscErrorCode DMMoFEMSetIsPartitioned(DM dm, PetscBool is_partitioned)
Definition DMMoFEM.cpp:1109
PetscErrorCode DMMoFEMCreateSubDM(DM subdm, DM dm, const char problem_name[])
Must be called by user to set Sub DM MoFEM data structures.
Definition DMMoFEM.cpp:219
PetscErrorCode DMMoFEMAddElement(DM dm, std::string fe_name)
add element to dm
Definition DMMoFEM.cpp:483
PetscErrorCode DMMoFEMSetSquareProblem(DM dm, PetscBool square_problem)
set squared problem
Definition DMMoFEM.cpp:442
PetscErrorCode DMMoFEMAddSubFieldRow(DM dm, const char field_name[])
Definition DMMoFEM.cpp:242
virtual bool check_field(const std::string &name) const =0
check if field is in database
Tensors class implemented by Walter Landry.
Definition FTensor.hpp:51
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
implementation of Data Operators for Forces and Sources
Definition Common.hpp:10
PetscErrorCode DMMoFEMSetDestroyProblem(DM dm, PetscBool destroy_problem)
Definition DMMoFEM.cpp:424
DEPRECATED auto createSmartDM(MPI_Comm comm, const std::string dm_type_name)
MoFEMErrorCode set_contact_dm(MoFEM::Interface &m_field, SmartPetscObj< DM > &dm, SmartPetscObj< DM > &dm_contact)
MoFEMErrorCode set_plastic_tau_dm(MoFEM::Interface &m_field, SmartPetscObj< DM > &dm, SmartPetscObj< DM > &dm_plastic_tau)
MoFEMErrorCode set_plastic_ep_dm(MoFEM::Interface &m_field, SmartPetscObj< DM > &dm, SmartPetscObj< DM > &dm_plastic)
virtual MPI_Comm & get_comm() const =0
Deprecated interface functions.
Matrix manager is used to build and partition problems.
Simple interface for fast problem set-up.
Definition Simple.hpp:27
intrusive_ptr for managing petsc objects
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.