|
Alamo
|
Go to the documentation of this file.
11 #ifndef MODEL_SOLID_AFFINE_HEXAGONAL_H_
12 #define MODEL_SOLID_AFFINE_HEXAGONAL_H_
56 virtual void Print(std::ostream &out)
const override
70 ret.
F0 = Set::Matrix::Zero();
90 static Hexagonal Combine(
const std::vector<Hexagonal> &models,
const std::vector<Set::Scalar> &eta,
int order)
94 ret.
F0 = Set::Matrix::Zero();
98 for (
unsigned int n = 0 ; n < models.size(); n++) etasum += eta[n];
99 for (
unsigned int n = 0 ; n < models.size(); n++)
101 ret.
ddw += models[n].ddw * (eta[n] / etasum);
102 ret.
F0 += models[n].F0 * (eta[n] / etasum);
110 #define OP_CLASS Hexagonal
111 #define OP_VARS X(ddw) X(F0)
131 if (i==0)
return name +
".Cxxxx";
139 for (amrex::MFIter mfi(a_dst, amrex::TilingIfNotGPU()); mfi.isValid(); ++mfi)
141 const amrex::Box& bx = mfi.growntilebox(amrex::IntVect(a_nghost));
144 amrex::Array4<const Model::Solid::Affine::Hexagonal>
const & src = ((*this)[a_lev])->array(mfi);
145 amrex::Array4<Set::Scalar>
const & dst = a_dst.array(mfi);
146 amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE(
int i,
int j,
int k)
148 dst(i,j,k,a_dstcomp + 0) = src(i,j,k).ddw(0,0,0,0);
static const KinematicVariable kinvar
Set::Matrix DW(const Set::Matrix &gradu) const override
static Hexagonal Random(Set::Scalar C11, Set::Scalar C12, Set::Scalar C13, Set::Scalar C33, Set::Scalar C44)
Set::Scalar W(const Set::Matrix &gradu) const override
Set::Matrix4< AMREX_SPACEDIM, Set::Sym::MajorMinor > DDW(const Set::Matrix &) const override
static Hexagonal Random()
#define ALAMO_SINGLE_DEFINITION
AMREX_FORCE_INLINE void SetF0(Set::Matrix &a_F0)
void Copy(int, amrex::MultiFab &, int, int) const
static void Parse(Hexagonal &value, IO::ParmParse &pp)
void Define(Set::Scalar C11, Set::Scalar C12, Set::Scalar C13, Set::Scalar C33, Set::Scalar C44, Set::Scalar phi1, Set::Scalar Phi, Set::Scalar phi2, Set::Matrix a_F0=Set::Matrix::Zero())
static void Parse(Hexagonal &value, IO::ParmParse &pp)
void Define(Set::Scalar C11, Set::Scalar C12, Set::Scalar C13, Set::Scalar C33, Set::Scalar C44, Set::Scalar phi1, Set::Scalar Phi, Set::Scalar phi2)
std::string Name(int) const
virtual void Print(std::ostream &out) const override
Hexagonal(Linear::Hexagonal base)
Set::Matrix4< AMREX_SPACEDIM, Set::Sym::MajorMinor > DDW(const Set::Matrix &gradu) const override
Eigen::Matrix< amrex::Real, AMREX_SPACEDIM, AMREX_SPACEDIM > Matrix
static AMREX_FORCE_INLINE Hexagonal Combine(const std::vector< Hexagonal > &models, const std::vector< Set::Scalar > &eta, int order)
void Exception(std::string file, std::string func, int line, Args const &... args)
static Hexagonal Random()
Set::Matrix DW(const Set::Matrix &gradu) const override
Set::Scalar W(const Set::Matrix &gradu) const override
void Define(Set::Scalar C11, Set::Scalar C12, Set::Scalar C13, Set::Scalar C33, Set::Scalar C44, Eigen::Matrix3d R, Set::Matrix a_F0=Set::Matrix::Zero())
Set::Matrix4< AMREX_SPACEDIM, Set::Sym::MajorMinor > ddw