Public Member Functions | Public Attributes | Private Attributes | List of all members
ScaleCoupledBeam Class Reference
+ Inheritance diagram for ScaleCoupledBeam:

Public Member Functions

 ScaleCoupledBeam ()
 
void setupOomph ()
 
void setupMercury ()
 
void actionsBeforeSolve () override
 Write header of output file. More...
 
void actionsAfterSolve () override
 Write header of output file. More...
 
void actionsBeforeOomphTimeStep () override
 Each time step, compute deflection, elastic, kinetic and gravitational energy, and write to output file. More...
 
double getBeamDeflection () const
 Computes beam deflection. More...
 
 ScaleCoupledBeam ()
 
void setupOomph ()
 
void setupMercury ()
 
void actionsBeforeSolve () override
 Write header of output file. More...
 
void actionsAfterSolve () override
 Write header of output file. More...
 
void actionsBeforeOomphTimeStep () override
 Each time step, compute deflection, elastic, kinetic and gravitational energy, and write to output file. More...
 
double getBeamDeflection () const
 Computes beam deflection. More...
 
- Public Member Functions inherited from ScaleCoupling< M, O >
const std::vector< CoupledElement > & getCoupledElements ()
 get coupled element More...
 
const std::vector< CoupledParticle > & getCoupledParticles ()
 get coupled particle More...
 
void setPenalty (double penalty)
 Sets penalty parameter. More...
 
void setCouplingWeight (const std::function< double(double x, double y, double z)> &couplingWeight)
 Sets weight function (determines which nodes/el_pts are in coupling zone) More...
 
void solveScaleCoupling ()
 Computes nStep, the ratio of FEM and DEM time steps, then calls solveSurfaceCoupling(nStep) More...
 
- Public Member Functions inherited from BaseCoupling< M, O >
 BaseCoupling ()=default
 
void setName (std::string name)
 
std::string getName () const
 
void removeOldFiles () const
 
void writeEneTimeStep (std::ostream &os) const override
 
void writeEneHeader (std::ostream &os) const override
 
void solveOomph (int max_adapt=0)
 
void checkResidual ()
 
void solveMercury (unsigned long nt)
 
void setCGWidth (const double &width)
 
double getCGWidth ()
 
bool useCGMapping ()
 
CGFunctions::LucyXYZ getCGFunction ()
 

Public Attributes

std::ofstream out
 output file stream More...
 
std::ofstream outFile
 output file stream More...
 

Private Attributes

double length = 20.0
 
double distance = 0.4
 
double bulkDensity = 1309
 
double elasticModulus = 1e8
 
double overlapLength = 10.0*distance
 
double penalty = 8e9
 
double velocity = 1e-1
 

Detailed Description

Define a coupled problem

Constructor & Destructor Documentation

◆ ScaleCoupledBeam() [1/2]

ScaleCoupledBeam::ScaleCoupledBeam ( )
inline
26  {
27  //set name
28  setName("ScaleCoupledBeam");
29  //remove existing output files
30  removeOldFiles();
31 
32  // setup steps
33  setupOomph();
34  setupMercury(); //sets timeMax
35 
36  // set weight function
37  auto couplingWeight = [this] (double x,double y,double z) {
38  if (x<0.5*(length-overlapLength)) {
39  // DEM region
40  return 1.0;
41  } else if (x>0.5*(length+overlapLength)) {
42  // FEM region
43  return 0.0;
44  } else {
45  // coupled region
46  //return (0.5*(length+overlapLength)-x)/overlapLength;
47  double unit = (0.5*(length+overlapLength)-x)/overlapLength;
48  return 0.5+0.4*std::sin(constants::pi*(unit-0.5));
49  }
50  };
51  setCouplingWeight(couplingWeight);
52  setPenalty(penalty);
53 
54  // setup time
55  double waveSpeed = sqrt(elasticModulus/bulkDensity);
56  double propagationTime = length/waveSpeed;
57  setTimeMax(1.25*propagationTime);
58  setOomphTimeStep(getTimeStep());
59  //setOomphTimeStep(propagationTime/100.0);
60  setSaveCount(getTimeMax()/getTimeStep()/100.0);
61  logger(INFO,"TimeMax: %, nTimeSteps %", getTimeMax(), getTimeMax()/getTimeStep());
62 
63  // solve
64  writeToVTK();
65  solveScaleCoupling();
66  saveSolidMesh();
67  }
sin
AnnoyingScalar sin(const AnnoyingScalar &x)
Definition: AnnoyingScalar.h:137
sqrt
AnnoyingScalar sqrt(const AnnoyingScalar &x)
Definition: AnnoyingScalar.h:134
logger
Logger< MERCURYDPM_LOGLEVEL > logger("MercuryKernel")
Definition of different loggers with certain modules. A user can define its own custom logger here.
Log::INFO
@ INFO
BaseCoupling::removeOldFiles
void removeOldFiles() const
Definition: BaseCoupling.h:47
BaseCoupling::setName
void setName(std::string name)
Definition: BaseCoupling.h:38
ScaleCoupledBeam::elasticModulus
double elasticModulus
Definition: ScaleCoupledBeam.cpp:19
ScaleCoupledBeam::bulkDensity
double bulkDensity
Definition: ScaleCoupledBeam.cpp:18
ScaleCoupledBeam::overlapLength
double overlapLength
Definition: ScaleCoupledBeam.cpp:20
ScaleCoupledBeam::penalty
double penalty
Definition: ScaleCoupledBeam.cpp:21
ScaleCoupledBeam::setupOomph
void setupOomph()
Definition: ScaleCoupledBeam.cpp:69
ScaleCoupledBeam::setupMercury
void setupMercury()
Definition: ScaleCoupledBeam.cpp:81
ScaleCoupledBeam::length
double length
Definition: ScaleCoupledBeam.cpp:16
ScaleCoupling::setCouplingWeight
void setCouplingWeight(const std::function< double(double x, double y, double z)> &couplingWeight)
Sets weight function (determines which nodes/el_pts are in coupling zone)
Definition: ScaleCoupling.h:100
ScaleCoupling::solveScaleCoupling
void solveScaleCoupling()
Computes nStep, the ratio of FEM and DEM time steps, then calls solveSurfaceCoupling(nStep)
Definition: ScaleCoupling.h:105
ScaleCoupling::setPenalty
void setPenalty(double penalty)
Sets penalty parameter.
Definition: ScaleCoupling.h:95
y
Scalar * y
Definition: level1_cplx_impl.h:128
constants::pi
const Mdouble pi
Definition: ExtendedMath.h:23
plotDoE.x
list x
Definition: plotDoE.py:28

References bulkDensity, elasticModulus, INFO, length, logger, overlapLength, penalty, constants::pi, BaseCoupling< M, O >::removeOldFiles(), ScaleCoupling< M, O >::setCouplingWeight(), BaseCoupling< M, O >::setName(), ScaleCoupling< M, O >::setPenalty(), setupMercury(), setupOomph(), sin(), ScaleCoupling< M, O >::solveScaleCoupling(), sqrt(), plotDoE::x, and y.

◆ ScaleCoupledBeam() [2/2]

ScaleCoupledBeam::ScaleCoupledBeam ( )
inline
26  {
27  //set name
28  setName("ScaleCoupledBeamUnitTest");
29  //remove existing output files
30  removeOldFiles();
31 
32  // setup steps
33  setupOomph();
34  setupMercury(); //sets timeMax
35 
36  // set weight function
37  auto couplingWeight = [this] (double x,double y,double z) {
38  if (x<0.5*(length-overlapLength)) {
39  // DEM region
40  return 1.0;
41  } else if (x>0.5*(length+overlapLength)) {
42  // FEM region
43  return 0.0;
44  } else {
45  // coupled region
46  //return (0.5*(length+overlapLength)-x)/overlapLength;
47  double unit = (0.5*(length+overlapLength)-x)/overlapLength;
48  return 0.5+0.4*std::sin(constants::pi*(unit-0.5));
49  }
50  };
51  setCouplingWeight(couplingWeight);
52  setPenalty(penalty);
53 
54  // setup time
55  setTimeMax(2.0*getTimeStep()); // different to ScaleCoupledBeam
56  setOomphTimeStep(getTimeStep());
57  setSaveCount(getTimeMax()/getTimeStep()/100.0);
58  logger(INFO,"TimeMax: %, nTimeSteps %", getTimeMax(), getTimeMax()/getTimeStep());
59 
60  // solve
61  writeToVTK();
62  solveScaleCoupling();
63  saveSolidMesh();
64  }

References INFO, length, logger, overlapLength, penalty, constants::pi, BaseCoupling< M, O >::removeOldFiles(), ScaleCoupling< M, O >::setCouplingWeight(), BaseCoupling< M, O >::setName(), ScaleCoupling< M, O >::setPenalty(), setupMercury(), setupOomph(), sin(), ScaleCoupling< M, O >::solveScaleCoupling(), plotDoE::x, and y.

Member Function Documentation

◆ actionsAfterSolve() [1/2]

void ScaleCoupledBeam::actionsAfterSolve ( )
inlineoverride

Write header of output file.

129  {
130  out.close();
131  }
ScaleCoupledBeam::out
std::ofstream out
output file stream
Definition: ScaleCoupledBeam.cpp:157

References out.

◆ actionsAfterSolve() [2/2]

void ScaleCoupledBeam::actionsAfterSolve ( )
inlineoverride

Write header of output file.

121  {
122  outFile.close();
123  // checking errors
124  Vec3D couplingForce = getCoupledParticles().back().couplingForce;
125  helpers::check(couplingForce.X,-15917,0.1,"Coupling force");
126  }
ScaleCoupledBeam::outFile
std::ofstream outFile
output file stream
Definition: ScaleCoupledBeamUnitTest.cpp:152
ScaleCoupling::getCoupledParticles
const std::vector< CoupledParticle > & getCoupledParticles()
get coupled particle
Definition: ScaleCoupling.h:90
Vec3D
Definition: Kernel/Math/Vector.h:30
Vec3D::X
Mdouble X
the vector components
Definition: Kernel/Math/Vector.h:45
helpers::check
void check(double real, double ideal, double error, std::string errorMessage)
Definition: TestHelpers.cc:16

References helpers::check(), ScaleCoupling< M, O >::getCoupledParticles(), outFile, and Vec3D::X.

◆ actionsBeforeOomphTimeStep() [1/2]

void ScaleCoupledBeam::actionsBeforeOomphTimeStep ( )
inlineoverride

Each time step, compute deflection, elastic, kinetic and gravitational energy, and write to output file.

134  {
135  double mass, elasticEnergy, kineticEnergy;
136  Vector<double> com(3), linearMomentum(3), angularMomentum(3);
137  getMassMomentumEnergy(mass, com, linearMomentum, angularMomentum, elasticEnergy, kineticEnergy);
138  static double comZ0 = com[2];
139  double gravEnergy = 9.8*mass*(com[2]-comZ0);
140  out << getOomphTime() << ' ' << getBeamDeflection() << ' ' << elasticEnergy << ' ' << kineticEnergy << ' ' << gravEnergy << std::endl;
141  std::cout << getOomphTime() << ' ' << getBeamDeflection() << ' ' << elasticEnergy << ' ' << kineticEnergy << ' ' << gravEnergy << '\n';
142  }
ScaleCoupledBeam::getBeamDeflection
double getBeamDeflection() const
Computes beam deflection.
Definition: ScaleCoupledBeam.cpp:145

References getBeamDeflection(), and out.

◆ actionsBeforeOomphTimeStep() [2/2]

void ScaleCoupledBeam::actionsBeforeOomphTimeStep ( )
inlineoverride

Each time step, compute deflection, elastic, kinetic and gravitational energy, and write to output file.

129  {
130  double mass, elasticEnergy, kineticEnergy;
131  Vector<double> com(3), linearMomentum(3), angularMomentum(3);
132  getMassMomentumEnergy(mass, com, linearMomentum, angularMomentum, elasticEnergy, kineticEnergy);
133  static double comZ0 = com[2];
134  double gravEnergy = 9.8*mass*(com[2]-comZ0);
135  outFile << getOomphTime() << ' ' << getBeamDeflection() << ' ' << elasticEnergy << ' ' << kineticEnergy << ' ' << gravEnergy << std::endl;
136  std::cout << getOomphTime() << ' ' << getBeamDeflection() << ' ' << elasticEnergy << ' ' << kineticEnergy << ' ' << gravEnergy << '\n';
137  }

References getBeamDeflection(), and outFile.

◆ actionsBeforeSolve() [1/2]

void ScaleCoupledBeam::actionsBeforeSolve ( )
inlineoverride

Write header of output file.

121  {
122  helpers::writeToFile(getName()+".gnu", "set key autotitle columnheader\n"
123  "p 'SolidBeamUnsteady.out' u 1:3 w l, '' u 1:4 w l, '' u 1:5 w l, '' u 1:($3+$4+$5) t 'totalEnergy' w l");
124  out.open(getName()+".out");
125  out << "time deflection elasticEnergy kineticEnergy gravEnergy\n";
126  }
BaseCoupling::getName
std::string getName() const
Definition: BaseCoupling.h:43
helpers::writeToFile
bool writeToFile(const std::string &filename, const std::string &filecontent)
Writes a string to a file.
Definition: FileIOHelpers.cc:29

References BaseCoupling< M, O >::getName(), out, and helpers::writeToFile().

◆ actionsBeforeSolve() [2/2]

void ScaleCoupledBeam::actionsBeforeSolve ( )
inlineoverride

Write header of output file.

113  {
114  helpers::writeToFile(getName()+".gnu", "set key autotitle columnheader\n"
115  "p 'SolidBeamUnsteady.out' u 1:3 w l, '' u 1:4 w l, '' u 1:5 w l, '' u 1:($3+$4+$5) t 'totalEnergy' w l");
116  outFile.open(getName()+".out");
117  outFile << "time deflection elasticEnergy kineticEnergy gravEnergy\n";
118  }

References BaseCoupling< M, O >::getName(), outFile, and helpers::writeToFile().

◆ getBeamDeflection() [1/2]

double ScaleCoupledBeam::getBeamDeflection ( ) const
inline

Computes beam deflection.

145  {
146  std::array<double, 3> min, max;
147  getDomainSize(min, max);
148 
149  Vector<double> xi(3);
150  xi[0] = max[0];
151  xi[1] = 0.5 * (max[1] + min[1]);
152  xi[2] = 0.5 * (max[2] + min[2]);
153  return getDeflection(xi, 2);
154  }
min
#define min(a, b)
Definition: datatypes.h:22
max
#define max(a, b)
Definition: datatypes.h:23

References max, and min.

Referenced by actionsBeforeOomphTimeStep().

◆ getBeamDeflection() [2/2]

double ScaleCoupledBeam::getBeamDeflection ( ) const
inline

Computes beam deflection.

140  {
141  std::array<double, 3> min, max;
142  getDomainSize(min, max);
143 
144  Vector<double> xi(3);
145  xi[0] = max[0];
146  xi[1] = 0.5 * (max[1] + min[1]);
147  xi[2] = 0.5 * (max[2] + min[2]);
148  return getDeflection(xi, 2);
149  }

References max, and min.

◆ setupMercury() [1/2]

void ScaleCoupledBeam::setupMercury ( )
inline
82  {
83  // positive overlap means adhesive forces
84  double overlap = 0.0*distance; //overlap between particles
85  double radius = 0.5*(distance+overlap); // particle radius
86  double particleDensity = bulkDensity*mathsFunc::cubic(distance) / (constants::pi/6.0*mathsFunc::cubic(2.0*radius));
87  double stiffness = elasticModulus*distance; // particle stiffness
88 
89  setDomain(Vec3D(-radius,-radius,0),Vec3D(radius,radius,0.5*(length+overlapLength)));
90  setParticlesWriteVTK(true);
91 
92  auto species = speciesHandler.copyAndAddObject(LinearViscoelasticReversibleAdhesiveSpecies());
93  species->setDensity(particleDensity);
94  species->setStiffness(stiffness);
95  species->setAdhesionStiffness(stiffness);
96  species->setAdhesionForceMax(stiffness*overlap);
97 
98  SphericalParticle p(species);
99  p.setRadius(radius);
100  //p.setVelocity(Vec3D(1,0,0));
101  double dx = distance;
102  auto n = (unsigned)(0.5*(length+overlapLength)/dx);
103  for (int i = 0; i<n; ++i) {
104  for (int j = 0; j<2; ++j) {
105  for (int k = 0; k<2; ++k) {
106  p.setPosition(Vec3D(1+2*i, j==0?-1:1, k==0?-1:1)*0.5*distance);
107  particleHandler.copyAndAddObject(p);
108  if (i==0) {
109  particleHandler.getLastObject()->setVelocity(Vec3D(velocity,0,0));
110  //particleHandler.getLastObject()->fixParticle();
111  //particleHandler.getLastObject()->setPrescribedVelocity(velocityAtBoundary);
112  }
113  }
114  }
115  }
116 
117  setTimeStep(0.5*species->computeTimeStep(p.getMass()));
118  }
i
int i
Definition: BiCGSTAB_step_by_step.cpp:9
n
const unsigned n
Definition: CG3DPackingUnitTest.cpp:11
LinearViscoelasticReversibleAdhesiveSpecies
Species< LinearViscoelasticNormalSpecies, EmptyFrictionSpecies, ReversibleAdhesiveSpecies > LinearViscoelasticReversibleAdhesiveSpecies
Definition: LinearViscoelasticReversibleAdhesiveSpecies.h:13
p
float * p
Definition: Tutorial_Map_using.cpp:9
ScaleCoupledBeam::distance
double distance
Definition: ScaleCoupledBeam.cpp:17
ScaleCoupledBeam::velocity
double velocity
Definition: ScaleCoupledBeam.cpp:22
SphericalParticle
A spherical particle is the most simple particle used in MercuryDPM.
Definition: SphericalParticle.h:16
k
char char char int int * k
Definition: level2_impl.h:374
UniformPSDSelfTest.radius
radius
Definition: UniformPSDSelfTest.py:15
mathsFunc::cubic
T cubic(const T val)
calculates the cube of a number
Definition: ExtendedMath.h:95
j
std::ptrdiff_t j
Definition: tut_arithmetic_redux_minmax.cpp:2

References bulkDensity, mathsFunc::cubic(), distance, elasticModulus, i, j, k, length, n, overlapLength, p, constants::pi, UniformPSDSelfTest::radius, and velocity.

Referenced by ScaleCoupledBeam().

◆ setupMercury() [2/2]

void ScaleCoupledBeam::setupMercury ( )
inline
79  {
80  // positive overlap means adhesive forces
81  double overlap = 0.0*distance; //overlap between particles
82  double radius = 0.5*(distance+overlap); // particle radius
83  double particleDensity = bulkDensity*mathsFunc::cubic(distance) / (constants::pi/6.0*mathsFunc::cubic(2.0*radius));
84  double stiffness = elasticModulus*distance; // particle stiffness
85 
86  setDomain(Vec3D(-radius,-radius,0),Vec3D(radius,radius,0.5*(length+overlapLength)));
87  setParticlesWriteVTK(true);
88 
89  auto species = speciesHandler.copyAndAddObject(LinearViscoelasticReversibleAdhesiveSpecies());
90  species->setDensity(particleDensity);
91  species->setStiffness(stiffness);
92  species->setAdhesionStiffness(stiffness);
93  species->setAdhesionForceMax(stiffness*overlap);
94 
95  SphericalParticle p(species);
96  p.setRadius(radius);
97  p.setVelocity(Vec3D(1,0,0)); // different to ScaleCoupledBeam
98  double dx = distance;
99  auto n = (unsigned)(0.5*(length+overlapLength)/dx);
100  for (int i = 0; i<n; ++i) {
101  for (int j = 0; j<2; ++j) {
102  for (int k = 0; k<2; ++k) {
103  p.setPosition(Vec3D(1+2*i, j==0?-1:1, k==0?-1:1)*0.5*distance);
104  particleHandler.copyAndAddObject(p); // different to ScaleCoupledBeam
105  }
106  }
107  }
108 
109  setTimeStep(0.5*species->computeTimeStep(p.getMass()));
110  }

References bulkDensity, mathsFunc::cubic(), distance, elasticModulus, i, j, k, length, n, overlapLength, p, constants::pi, and UniformPSDSelfTest::radius.

◆ setupOomph() [1/2]

void ScaleCoupledBeam::setupOomph ( )
inline
69  {
70  setElasticModulus(elasticModulus);
71  setDensity(bulkDensity);
72  setSolidCubicMesh(round(0.5*(length+overlapLength)/distance), 2, 2,
73  0.5*(length-overlapLength), length, -distance, distance, -distance, distance);
74  //pinBoundaries({Beam::Boundary::X_MIN});
75  setNewtonSolverTolerance(1e-2);
76  prepareForSolve();
77  linear_solver_pt()->disable_doc_time();
78  //disable_info_in_newton_solve();
79  }
e
Array< double, 1, 3 > e(1./3., 0.5, 2.)
Eigen::bfloat16_impl::round
EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC bfloat16 round(const bfloat16 &a)
Definition: BFloat16.h:646

References bulkDensity, distance, e(), elasticModulus, length, overlapLength, and Eigen::bfloat16_impl::round().

Referenced by ScaleCoupledBeam().

◆ setupOomph() [2/2]

void ScaleCoupledBeam::setupOomph ( )
inline
66  {
67  setElasticModulus(elasticModulus);
68  setDensity(bulkDensity);
69  setSolidCubicMesh(round(0.5*(length+overlapLength)/distance), 2, 2,
70  0.5*(length-overlapLength), length, -distance, distance, -distance, distance);
71  //pinBoundaries({Beam::Boundary::X_MIN});
72  setNewtonSolverTolerance(1e-2);
73  prepareForSolve();
74  linear_solver_pt()->disable_doc_time();
75  //disable_info_in_newton_solve();
76  }

References bulkDensity, distance, e(), elasticModulus, length, overlapLength, and Eigen::bfloat16_impl::round().

Member Data Documentation

◆ bulkDensity

double ScaleCoupledBeam::bulkDensity = 1309
private

Referenced by ScaleCoupledBeam(), setupMercury(), and setupOomph().

◆ distance

double ScaleCoupledBeam::distance = 0.4
private

Referenced by setupMercury(), and setupOomph().

◆ elasticModulus

double ScaleCoupledBeam::elasticModulus = 1e8
private

Referenced by ScaleCoupledBeam(), setupMercury(), and setupOomph().

◆ length

double ScaleCoupledBeam::length = 20.0
private

Referenced by ScaleCoupledBeam(), setupMercury(), and setupOomph().

◆ out

std::ofstream ScaleCoupledBeam::out

output file stream

Referenced by actionsAfterSolve(), actionsBeforeOomphTimeStep(), and actionsBeforeSolve().

◆ outFile

std::ofstream ScaleCoupledBeam::outFile

output file stream

Referenced by actionsAfterSolve(), actionsBeforeOomphTimeStep(), and actionsBeforeSolve().

◆ overlapLength

double ScaleCoupledBeam::overlapLength = 10.0*distance
private

Referenced by ScaleCoupledBeam(), setupMercury(), and setupOomph().

◆ penalty

double ScaleCoupledBeam::penalty = 8e9
private

Referenced by ScaleCoupledBeam().

◆ velocity

double ScaleCoupledBeam::velocity = 1e-1
private

Referenced by setupMercury().

The documentation for this class was generated from the following files: