forked from gismo/gsElasticity
-
Notifications
You must be signed in to change notification settings - Fork 0
/
gsElTimeIntegrator.h
152 lines (123 loc) · 4.99 KB
/
gsElTimeIntegrator.h
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
/** @file gsElTimeIntegrator.h
@brief Provides time integration for dynamical elasticity.
This file is part of the G+Smo library.
This Source Code Form is subject to the terms of the Mozilla Public
License, v. 2.0. If a copy of the MPL was not distributed with this
file, You can obtain one at http://mozilla.org/MPL/2.0/.
Author(s):
A.Shamanskiy (2016 - ...., TU Kaiserslautern)
*/
#pragma once
#include <gsElasticity/gsBaseAssembler.h>
#include <gsElasticity/gsBaseUtils.h>
namespace gismo
{
template <class T>
class gsElasticityAssembler;
template <class T>
class gsMassAssembler;
/** @brief Time integation for equations of dynamic elasticity with implicit schemes
*/
template <class T>
class gsElTimeIntegrator : public gsBaseAssembler<T>
{
public:
typedef gsBaseAssembler<T> Base;
/// constructor method. requires a gsElasticityAssembler for construction of the static linear system
/// and a gsMassAssembler for the mass matrix
gsElTimeIntegrator(gsElasticityAssembler<T> & stiffAssembler_,
gsMassAssembler<T> & massAssembler_);
/// @brief Returns the list of default options for assembly
static gsOptionList defaultOptions();
/// set intial conditions
void setDisplacementVector(const gsMatrix<T> & displacementVector)
{
GISMO_ENSURE(displacementVector.rows() == stiffAssembler.numDofs(),
"Wrong size of the displacement vector: " + util::to_string(displacementVector.rows()) +
". Must be: " + util::to_string(stiffAssembler.numDofs()));
dispVector = displacementVector;
initialized = false;
}
void setVelocityVector(const gsMatrix<T> & velocityVector)
{
GISMO_ENSURE(velocityVector.rows() == stiffAssembler.numDofs(),
"Wrong size of the velocity vector: " + util::to_string(velocityVector.rows()) +
". Must be: " + util::to_string(stiffAssembler.numDofs()));
velVector = velocityVector;
initialized = false;
}
/// make a time step according to a chosen scheme
void makeTimeStep(T timeStep);
/// assemble the linear system for the nonlinear solver
virtual bool assemble(const gsMatrix<T> & solutionVector,
const std::vector<gsMatrix<T> > & fixedDoFs);
/// return the number of free degrees of freedom
virtual int numDofs() const { return stiffAssembler.numDofs(); }
/// returns vector of displacement DoFs
const gsMatrix<T> & displacementVector() const
{
GISMO_ENSURE(dispVector.rows() == stiffAssembler.numDofs(),
"No initial conditions provided!");
return dispVector;
}
/// returns vector of velocity DoFs
const gsMatrix<T> & velocityVector() const
{
GISMO_ENSURE(velVector.rows() == stiffAssembler.numDofs(),
"No initial conditions provided!");
return velVector;
}
/// save solver state
void saveState();
/// recover solver state from saved state
void recoverState();
/// number of iterations Newton's method required at the last time step
index_t numberIterations() const { return numIters;}
/// construct solution using the stiffness assembler
void constructSolution(gsMultiPatch<T> & solution) const;
/// assemblers' accessors
gsBaseAssembler<T> & mAssembler() { return massAssembler; }
gsBaseAssembler<T> & assembler() { return stiffAssembler; }
protected:
void initialize();
/// time integraton schemes
gsMatrix<T> implicitLinear();
gsMatrix<T> implicitNonlinear();
/// time integration scheme coefficients
T alpha1() {return 1./m_options.getReal("Beta")/pow(tStep,2); }
T alpha2() {return 1./m_options.getReal("Beta")/tStep; }
T alpha3() {return (1-2*m_options.getReal("Beta"))/2/m_options.getReal("Beta"); }
T alpha4() {return m_options.getReal("Gamma")/m_options.getReal("Beta")/tStep; }
T alpha5() {return 1 - m_options.getReal("Gamma")/m_options.getReal("Beta"); }
T alpha6() {return (1-m_options.getReal("Gamma")/m_options.getReal("Beta")/2)*tStep; }
protected:
/// assembler object that generates the static system
gsElasticityAssembler<T> & stiffAssembler;
/// assembler object that generates the mass matrix
gsMassAssembler<T> & massAssembler;
/// initialization flag
bool initialized;
/// time step length
T tStep;
/// vector of displacement DoFs
gsMatrix<T> dispVector;
/// vector of velocity DoFs
gsMatrix<T> velVector;
/// vector of acceleration DoFs
gsMatrix<T> accVector;
using Base::m_system;
using Base::m_options;
using Base::m_ddof;
/// number of iterations Newton's method took to converge at the last time step
index_t numIters;
/// saved state
bool hasSavedState;
gsMatrix<T> dispVecSaved;
gsMatrix<T> velVecSaved;
gsMatrix<T> accVecSaved;
std::vector<gsMatrix<T> > ddofsSaved;
};
}
#ifndef GISMO_BUILD_LIB
#include GISMO_HPP_HEADER(gsElTimeIntegrator.hpp)
#endif