-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathstiffness_tensors.py
235 lines (202 loc) · 8.64 KB
/
stiffness_tensors.py
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
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
"""File containing classes for a fourth order stiffness tensor.
Analogous to the tensors.py module within PorePy, but this one provides a tensor
corresponding to an anisotropic medium."""
import numpy as np
class StiffnessTensorInnerVTI(object):
"""Cell-wise representation of fourth order tensor.
The tensor is represented by a (3^2, 3^2 ,Nc)-array, where Nc denotes the number of cells, i.e. the tensor values are stored discretely.
See pp.FourthOrderTensor in PorePy for more extensive documentation.
The only constructor available (here, in this file) for the moment is based on five
independent material parameters. It should also be noted that the construction of
the stiffness tensor has one intended use as of now: We create a tensor that
represents an outer isotropic domain and an inner vertically transverse isotropic
inner domain. Therefore you will see 7 material parameter attributes instead of only
5. 2 are for the isotropic domain, 5 are for the anisotropic one.
Attributes:
values (numpy.ndarray): dimensions (3^2, 3^2, nc), cell-wise
representation of the stiffness matrix.
lmbda (np.ndarray): Nc array of first Lame parameter.
mu (np.ndarray): Nc array of second Lame parameter.
mu_parallel (np.ndarray): Nc array of transverse shear parameter.
mu_orthogonal (np.ndarray): Nc array of transverse-to-perpendicular shear
parameter.
lambda_parallell (np.ndarray): Transverse compressive stress parameter.
lambda_orthogonal (np.ndarray): Perpendicular compressive stress parameter.
volumetric_compr_lambda (np.ndarray): Volumetric compressive stress parameter.
"""
def __init__(
self,
mu: np.ndarray,
lmbda: np.ndarray,
mu_parallel: np.ndarray,
mu_orthogonal: np.ndarray,
lambda_parallel: np.ndarray,
lambda_orthogonal: np.ndarray,
volumetric_compr_lambda: np.ndarray,
):
"""Constructor for the fourth order tensor.
Parameters:
lmbda: Nc array of first Lame parameter.
mu: Nc array of second Lame parameter.
mu_parallel: Nc array of transverse shear parameter.
mu_orthogonal: Nc array of transverse-to-perpendicular shear parameter.
lambda_parallell: Transverse compressive stress parameter.
lambda_orthogonal: Perpendicular compressive stress parameter.
volumetric_compr_lambda: Volumetric compressive stress parameter.
"""
# Creating attributes of the values, as this might come in handy later
self.lmbda = lmbda
self.mu = mu
self.mu_parallel = mu_parallel
self.mu_orthogonal = mu_orthogonal
self.volumetric_compr_lambda = volumetric_compr_lambda
self.lambda_parallel = lambda_parallel
self.lambda_orthogonal = lambda_orthogonal
# Basis for mu and lambda contribution for tensor build
(
l_mat,
l_par_mat,
l_ort_mat,
m_par_mat,
m_ort_mat,
) = self.hardcoded_mu_lam_basis()
# Matrices for the isotropic tensor build
lmbda_mat = l_mat
mu_mat = m_par_mat + m_ort_mat
# Adding axis to isotropic related matrices
mu_mat = mu_mat[:, :, np.newaxis]
lmbda_mat = lmbda_mat[:, :, np.newaxis]
# Adding axis to anisotropy related matrices
l_mat = l_mat[:, :, np.newaxis]
l_par_mat = l_par_mat[:, :, np.newaxis]
l_ort_mat = l_ort_mat[:, :, np.newaxis]
m_par_mat = m_par_mat[:, :, np.newaxis]
m_ort_mat = m_ort_mat[:, :, np.newaxis]
c = (
# Isotropic outer domain
mu_mat * mu
+ lmbda_mat * lmbda
# VTI inner domain
+ l_mat * volumetric_compr_lambda
+ l_par_mat * lambda_parallel
+ l_ort_mat * lambda_orthogonal
+ m_par_mat * mu_parallel
+ m_ort_mat * mu_orthogonal
)
self.values = c
def copy(self) -> "StiffnessTensorInnerVTI":
"""
Define a deep copy of the tensor.
Returns:
StiffnessTensorInnerVTI: New tensor with identical fields, but separate
arrays (in the memory sense).
"""
C = StiffnessTensorInnerVTI(
mu=self.mu,
lmbda=self.lmbda,
mu_parallel=self.mu_parallel,
mu_orthogonal=self.mu_orthogonal,
lambda_parallel=self.lambda_parallel,
lambda_orthogonal=self.lambda_orthogonal,
volumetric_compr_lambda=self.volumetric_compr_lambda,
)
C.values = self.values.copy()
return C
def hardcoded_mu_lam_basis(self) -> tuple:
"""Basis for the contributions of mu and lambda in a VTI stiffness tensor.
This method contains the five matrices needed for creating a vertical
transversely isotropic stiffness tensor.
Returns:
A tuple of all the matrices in the following order:
lambda_mat,
lambda_parallel_mat,
lambda_orthogonal_mat,
mu_parallel_mat,
mu_orthogonal_mat.
Additional notes:
lambda_mat is the basis related to the volumetric compressive stress
parameter.
The two mu-matrices, namely mu_parallel_mat and mu_orthogonal_mat, are the
bases related to the transverse and transverse-to-perpendicular shear
parameter, respectively.
The two remaining matrices, lambda_parallel_mat and lambda_orthogonal_mat,
are the bases related to the transverse compressive stress parameter and the
perpendicular compressive stress parameter, respectively.
For an isotropic media, lambda_mat is kept as is, and the two mu-matrices
are summed up (and provided the same values for orthogonal and parallel mu).
The remaining lambda matrices contribute nothing, as the transverse
compressive stress parameter and perpendicular compressive stress parameter
are then set to zero.
"""
lambda_mat = np.array(
[
[1, 0, 0, 0, 1, 0, 0, 0, 1],
[0, 0, 0, 0, 0, 0, 0, 0, 0],
[0, 0, 0, 0, 0, 0, 0, 0, 0],
[0, 0, 0, 0, 0, 0, 0, 0, 0],
[1, 0, 0, 0, 1, 0, 0, 0, 1],
[0, 0, 0, 0, 0, 0, 0, 0, 0],
[0, 0, 0, 0, 0, 0, 0, 0, 0],
[0, 0, 0, 0, 0, 0, 0, 0, 0],
[1, 0, 0, 0, 1, 0, 0, 0, 1],
]
)
lambda_parallel_mat = np.array(
[
[1, 0, 0, 0, 1, 0, 0, 0, 0],
[0, 0, 0, 0, 0, 0, 0, 0, 0],
[0, 0, 0, 0, 0, 0, 0, 0, 0],
[0, 0, 0, 0, 0, 0, 0, 0, 0],
[1, 0, 0, 0, 1, 0, 0, 0, 0],
[0, 0, 0, 0, 0, 0, 0, 0, 0],
[0, 0, 0, 0, 0, 0, 0, 0, 0],
[0, 0, 0, 0, 0, 0, 0, 0, 0],
[0, 0, 0, 0, 0, 0, 0, 0, 0],
]
)
lambda_orthogonal_mat = np.array(
[
[0, 0, 0, 0, 0, 0, 0, 0, 0],
[0, 0, 0, 0, 0, 0, 0, 0, 0],
[0, 0, 0, 0, 0, 0, 0, 0, 0],
[0, 0, 0, 0, 0, 0, 0, 0, 0],
[0, 0, 0, 0, 0, 0, 0, 0, 0],
[0, 0, 0, 0, 0, 0, 0, 0, 0],
[0, 0, 0, 0, 0, 0, 0, 0, 0],
[0, 0, 0, 0, 0, 0, 0, 0, 0],
[0, 0, 0, 0, 0, 0, 0, 0, 1],
]
)
mu_parallel_mat = np.array(
[
[2, 0, 0, 0, 0, 0, 0, 0, 0],
[0, 1, 0, 1, 0, 0, 0, 0, 0],
[0, 0, 0, 0, 0, 0, 0, 0, 0],
[0, 1, 0, 1, 0, 0, 0, 0, 0],
[0, 0, 0, 0, 2, 0, 0, 0, 0],
[0, 0, 0, 0, 0, 0, 0, 0, 0],
[0, 0, 0, 0, 0, 0, 0, 0, 0],
[0, 0, 0, 0, 0, 0, 0, 0, 0],
[0, 0, 0, 0, 0, 0, 0, 0, 0],
]
)
mu_orthogonal_mat = np.array(
[
[0, 0, 0, 0, 0, 0, 0, 0, 0],
[0, 0, 0, 0, 0, 0, 0, 0, 0],
[0, 0, 1, 0, 0, 0, 1, 0, 0],
[0, 0, 0, 0, 0, 0, 0, 0, 0],
[0, 0, 0, 0, 0, 0, 0, 0, 0],
[0, 0, 0, 0, 0, 1, 0, 1, 0],
[0, 0, 1, 0, 0, 0, 1, 0, 0],
[0, 0, 0, 0, 0, 1, 0, 1, 0],
[0, 0, 0, 0, 0, 0, 0, 0, 2],
]
)
return (
lambda_mat,
lambda_parallel_mat,
lambda_orthogonal_mat,
mu_parallel_mat,
mu_orthogonal_mat,
)