OOFEM  2.4
OOFEM.org - Object Oriented Finite Element Solver
lsmastermat.C
Go to the documentation of this file.
1 /*
2  *
3  * ##### ##### ###### ###### ### ###
4  * ## ## ## ## ## ## ## ### ##
5  * ## ## ## ## #### #### ## # ##
6  * ## ## ## ## ## ## ## ##
7  * ## ## ## ## ## ## ## ##
8  * ##### ##### ## ###### ## ##
9  *
10  *
11  * OOFEM : Object Oriented Finite Element Code
12  *
13  * Copyright (C) 1993 - 2013 Borek Patzak
14  *
15  *
16  *
17  * Czech Technical University, Faculty of Civil Engineering,
18  * Department of Structural Mechanics, 166 29 Prague, Czech Republic
19  *
20  * This library is free software; you can redistribute it and/or
21  * modify it under the terms of the GNU Lesser General Public
22  * License as published by the Free Software Foundation; either
23  * version 2.1 of the License, or (at your option) any later version.
24  *
25  * This program is distributed in the hope that it will be useful,
26  * but WITHOUT ANY WARRANTY; without even the implied warranty of
27  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
28  * Lesser General Public License for more details.
29  *
30  * You should have received a copy of the GNU Lesser General Public
31  * License along with this library; if not, write to the Free Software
32  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
33  */
34 
35 #include "lsmastermat.h"
36 #include "Materials/isolinearelasticmaterial.h"
37 #include "gausspoint.h"
38 #include "floatmatrix.h"
39 #include "floatarray.h"
40 #include "intarray.h"
41 #include "stressvector.h"
42 #include "strainvector.h"
43 #include "mathfem.h"
44 #include "contextioerr.h"
45 #include "datastream.h"
46 #include "classfactory.h"
47 
48 namespace oofem {
49 REGISTER_Material(LargeStrainMasterMaterial);
50 
51 // constructor
52 LargeStrainMasterMaterial :: LargeStrainMasterMaterial(int n, Domain *d) : StructuralMaterial(n, d)
53 {
54  slaveMat = 0;
55 }
56 
57 // destructor
58 LargeStrainMasterMaterial :: ~LargeStrainMasterMaterial()
59 { }
60 
61 // reads the model parameters from the input file
62 IRResultType
63 LargeStrainMasterMaterial :: initializeFrom(InputRecord *ir)
64 {
65  IRResultType result; // required by IR_GIVE_FIELD macro
66 
67  IR_GIVE_OPTIONAL_FIELD(ir, slaveMat, _IFT_LargeStrainMasterMaterial_slaveMat); // number of slave material
68  IR_GIVE_OPTIONAL_FIELD(ir, m, _IFT_LargeStrainMasterMaterial_m); // type of Set-Hill strain tensor
69 
70  return IRRT_OK;
71 }
72 
73 // creates a new material status corresponding to this class
74 MaterialStatus *
75 LargeStrainMasterMaterial :: CreateStatus(GaussPoint *gp) const
76 {
77  return new LargeStrainMasterMaterialStatus(1, this->giveDomain(), gp, slaveMat);
78 }
79 
80 
81 void
82 LargeStrainMasterMaterial :: giveFirstPKStressVector_3d(FloatArray &answer, GaussPoint *gp, const FloatArray &vF, TimeStep *tStep)
83 {
84  LargeStrainMasterMaterialStatus *status = static_cast< LargeStrainMasterMaterialStatus * >( this->giveStatus(gp) );
85  this->initTempStatus(gp);
86 
87  StructuralMaterial *sMat = static_cast< StructuralMaterial * >( domain->giveMaterial(slaveMat) );
88 
89  double lambda1, lambda2, lambda3, E1, E2, E3;
90  FloatArray eVals, SethHillStrainVector, stressVector, stressM;
91  FloatMatrix F, Ft, C, eVecs, SethHillStrain;
92  FloatMatrix L1, L2, T;
93  //store of deformation gradient into 3x3 matrix
94  F.beMatrixForm(vF);
95  //compute right Cauchy-Green tensor(C), its eigenvalues and eigenvectors
96  C.beTProductOf(F, F);
97  // compute eigen values and eigen vectors of C
98  C.jaco_(eVals, eVecs, 15);
99  // compute Seth - Hill's strain measure, it depends on mParameter
100  lambda1 = eVals.at(1);
101  lambda2 = eVals.at(2);
102  lambda3 = eVals.at(3);
103  if ( m == 0 ) {
104  E1 = 1. / 2. * log(lambda1);
105  E2 = 1. / 2. * log(lambda2);
106  E3 = 1. / 2. * log(lambda3);
107  } else {
108  E1 = 1. / ( 2. * m ) * ( pow(lambda1, m) - 1. );
109  E2 = 1. / ( 2. * m ) * ( pow(lambda2, m) - 1. );
110  E3 = 1. / ( 2. * m ) * ( pow(lambda3, m) - 1. );
111  }
112 
113  SethHillStrain.resize(3, 3);
114  for ( int i = 1; i < 4; i++ ) {
115  for ( int j = 1; j < 4; j++ ) {
116  SethHillStrain.at(i, j) = E1 * eVecs.at(i, 1) * eVecs.at(j, 1) + E2 *eVecs.at(i, 2) * eVecs.at(j, 2) + E3 *eVecs.at(i, 3) * eVecs.at(j, 3);
117  }
118  }
119 
120  SethHillStrainVector.beSymVectorFormOfStrain(SethHillStrain);
121  sMat->giveRealStressVector_3d(stressVector, gp, SethHillStrainVector, tStep);
122  this->constructTransformationMatrix(T, eVecs);
123 
124  stressVector.at(4) = 2 * stressVector.at(4);
125  stressVector.at(5) = 2 * stressVector.at(5);
126  stressVector.at(6) = 2 * stressVector.at(6);
127 
128 
129  stressM.beProductOf(T, stressVector);
130  stressM.at(4) = 1. / 2. * stressM.at(4);
131  stressM.at(5) = 1. / 2. * stressM.at(5);
132  stressM.at(6) = 1. / 2. * stressM.at(6);
133 
134  this->constructL1L2TransformationMatrices(L1, L2, eVals, stressM, E1, E2, E3);
135 
136  FloatMatrix junk, P, TL;
137  FloatArray secondPK;
138  junk.beProductOf(L1, T);
139  P.beTProductOf(T, junk);
140  //transformation of the stress to the 2PK stress and then to 1PK
141  stressVector.at(4) = 0.5 * stressVector.at(4);
142  stressVector.at(5) = 0.5 * stressVector.at(5);
143  stressVector.at(6) = 0.5 * stressVector.at(6);
144  secondPK.beProductOf(P, stressVector);
145  answer.beProductOf(F, secondPK); // P = F*S
146  junk.zero();
147  junk.beProductOf(L2, T);
148  TL.beTProductOf(T, junk);
149  status->setPmatrix(P);
150  status->setTLmatrix(TL);
151  status->letTempStressVectorBe(answer);
152 }
153 
154 
155 void
156 LargeStrainMasterMaterial :: constructTransformationMatrix(FloatMatrix &answer, const FloatMatrix &eVecs)
157 {
158  answer.resize(6, 6);
159  answer.at(1, 1) = eVecs.at(1, 1) * eVecs.at(1, 1);
160  answer.at(1, 2) = eVecs.at(2, 1) * eVecs.at(2, 1);
161  answer.at(1, 3) = eVecs.at(3, 1) * eVecs.at(3, 1);
162  answer.at(1, 4) = eVecs.at(2, 1) * eVecs.at(3, 1);
163  answer.at(1, 5) = eVecs.at(1, 1) * eVecs.at(3, 1);
164  answer.at(1, 6) = eVecs.at(1, 1) * eVecs.at(2, 1);
165 
166  answer.at(2, 1) = eVecs.at(1, 2) * eVecs.at(1, 2);
167  answer.at(2, 2) = eVecs.at(2, 2) * eVecs.at(2, 2);
168  answer.at(2, 3) = eVecs.at(3, 2) * eVecs.at(3, 2);
169  answer.at(2, 4) = eVecs.at(2, 2) * eVecs.at(3, 2);
170  answer.at(2, 5) = eVecs.at(1, 2) * eVecs.at(3, 2);
171  answer.at(2, 6) = eVecs.at(1, 2) * eVecs.at(2, 2);
172 
173  answer.at(3, 1) = eVecs.at(1, 3) * eVecs.at(1, 3);
174  answer.at(3, 2) = eVecs.at(2, 3) * eVecs.at(2, 3);
175  answer.at(3, 3) = eVecs.at(3, 3) * eVecs.at(3, 3);
176  answer.at(3, 4) = eVecs.at(2, 3) * eVecs.at(3, 3);
177  answer.at(3, 5) = eVecs.at(1, 3) * eVecs.at(3, 3);
178  answer.at(3, 6) = eVecs.at(1, 3) * eVecs.at(2, 3);
179 
180  answer.at(4, 1) = 2 * eVecs.at(1, 2) * eVecs.at(1, 3);
181  answer.at(4, 2) = 2 * eVecs.at(2, 2) * eVecs.at(2, 3);
182  answer.at(4, 3) = 2 * eVecs.at(3, 2) * eVecs.at(3, 3);
183  answer.at(4, 4) = eVecs.at(2, 2) * eVecs.at(3, 3) + eVecs.at(3, 2) * eVecs.at(2, 3);
184  answer.at(4, 5) = eVecs.at(1, 2) * eVecs.at(3, 3) + eVecs.at(3, 2) * eVecs.at(1, 3);
185  answer.at(4, 6) = eVecs.at(1, 2) * eVecs.at(2, 3) + eVecs.at(2, 2) * eVecs.at(1, 3);
186 
187  answer.at(5, 1) = 2 * eVecs.at(1, 1) * eVecs.at(1, 3);
188  answer.at(5, 2) = 2 * eVecs.at(2, 1) * eVecs.at(2, 3);
189  answer.at(5, 3) = 2 * eVecs.at(3, 1) * eVecs.at(3, 3);
190  answer.at(5, 4) = eVecs.at(2, 1) * eVecs.at(3, 3) + eVecs.at(3, 1) * eVecs.at(2, 3);
191  answer.at(5, 5) = eVecs.at(1, 1) * eVecs.at(3, 3) + eVecs.at(3, 1) * eVecs.at(1, 3);
192  answer.at(5, 6) = eVecs.at(1, 1) * eVecs.at(2, 3) + eVecs.at(2, 1) * eVecs.at(1, 3);
193 
194  answer.at(6, 1) = 2 * eVecs.at(1, 1) * eVecs.at(1, 2);
195  answer.at(6, 2) = 2 * eVecs.at(2, 1) * eVecs.at(2, 2);
196  answer.at(6, 3) = 2 * eVecs.at(3, 1) * eVecs.at(3, 2);
197  answer.at(6, 4) = eVecs.at(2, 1) * eVecs.at(3, 2) + eVecs.at(3, 1) * eVecs.at(2, 2);
198  answer.at(6, 5) = eVecs.at(1, 1) * eVecs.at(3, 2) + eVecs.at(3, 1) * eVecs.at(1, 2);
199  answer.at(6, 6) = eVecs.at(1, 1) * eVecs.at(2, 2) + eVecs.at(2, 1) * eVecs.at(1, 2);
200 }
201 
202 
203 void
204 LargeStrainMasterMaterial :: constructL1L2TransformationMatrices(FloatMatrix &answer1, FloatMatrix &answer2, const FloatArray &eigenValues, FloatArray &stressM, double E1, double E2, double E3)
205 {
206  double gamma12, gamma13, gamma23, gamma112, gamma221, gamma113, gamma331, gamma223, gamma332, gamma;
207  double lambda1 = eigenValues.at(1);
208  double lambda2 = eigenValues.at(2);
209  double lambda3 = eigenValues.at(3);
210  double lambda1P = pow(lambda1, m - 1);
211  double lambda2P = pow(lambda2, m - 1);
212  double lambda3P = pow(lambda3, m - 1);
213  if ( ( lambda1 == lambda2 ) && ( lambda2 == lambda3 ) ) { // three equal eigenvalues
214  gamma12 = gamma13 = gamma23 = 1. / 2. * lambda1P;
215 
216  answer2.at(1, 1) = 2 * stressM.at(1) * ( m - 1 ) * pow(lambda1, m - 2);
217  answer2.at(2, 2) = 2 * stressM.at(2) * ( m - 1 ) * pow(lambda2, m - 2);
218  answer2.at(3, 3) = 2 * stressM.at(3) * ( m - 1 ) * pow(lambda3, m - 2);
219  answer2.at(4, 4) = 1. / 2. * ( stressM.at(2) + stressM.at(3) ) * ( m - 1 ) * pow(lambda2, m - 2);
220  answer2.at(5, 5) = 1. / 2. * ( stressM.at(1) + stressM.at(3) ) * ( m - 1 ) * pow(lambda3, m - 2);
221  answer2.at(6, 6) = 1. / 2. * ( stressM.at(1) + stressM.at(2) ) * ( m - 1 ) * pow(lambda1, m - 2);
222 
223  answer2.at(1, 5) = answer2.at(5, 1) = stressM.at(5) * ( m - 1 ) * pow(lambda3, m - 2);
224  answer2.at(1, 6) = answer2.at(6, 1) = stressM.at(6) * ( m - 1 ) * pow(lambda1, m - 2);
225  answer2.at(2, 4) = answer2.at(4, 2) = stressM.at(4) * ( m - 1 ) * pow(lambda2, m - 2);
226  answer2.at(2, 6) = answer2.at(6, 2) = stressM.at(6) * ( m - 1 ) * pow(lambda1, m - 2);
227  answer2.at(3, 4) = answer2.at(4, 3) = stressM.at(4) * ( m - 1 ) * pow(lambda2, m - 2);
228  answer2.at(3, 5) = answer2.at(5, 3) = stressM.at(5) * ( m - 1 ) * pow(lambda3, m - 2);
229 
230  answer2.at(4, 5) = answer2.at(5, 4) = 1. / 2. * stressM.at(6) * ( m - 1 ) * pow(lambda1, m - 2);
231  answer2.at(4, 6) = answer2.at(6, 4) = 1. / 2. * stressM.at(5) * ( m - 1 ) * pow(lambda1, m - 2);
232  answer2.at(5, 6) = answer2.at(6, 5) = 1. / 2. * stressM.at(4) * ( m - 1 ) * pow(lambda1, m - 2);
233  } else if ( lambda1 == lambda2 ) { //two equal eigenvalues
234  gamma12 = 1. / 2. * lambda1P;
235  gamma13 = ( E1 - E3 ) / ( lambda1 - lambda3 );
236  gamma23 = ( E2 - E3 ) / ( lambda2 - lambda3 );
237  gamma113 = ( pow(lambda1, m - 1) * ( lambda1 - lambda3 ) - 2 * ( E1 - E3 ) ) / ( ( lambda1 - lambda3 ) * ( lambda1 - lambda3 ) );
238  gamma331 = ( pow(lambda3, m - 1) * ( lambda3 - lambda1 ) - 2 * ( E3 - E1 ) ) / ( ( lambda3 - lambda1 ) * ( lambda3 - lambda1 ) );
239 
240 
241  answer2.at(1, 1) = 2 * stressM.at(1) * ( m - 1 ) * pow(lambda1, m - 2);
242  answer2.at(2, 2) = 2 * stressM.at(2) * ( m - 1 ) * pow(lambda2, m - 2);
243  answer2.at(3, 3) = 2 * stressM.at(3) * ( m - 1 ) * pow(lambda3, m - 2);
244 
245  answer2.at(4, 4) = stressM.at(2) * gamma113 + stressM.at(3) * gamma331;
246  answer2.at(5, 5) = stressM.at(1) * gamma113 + stressM.at(3) * gamma331;
247  answer2.at(6, 6) = 1. / 2. * ( stressM.at(1) + stressM.at(2) ) * ( m - 1 ) * pow(lambda1, m - 2);
248 
249  answer2.at(1, 5) = answer2.at(5, 1) = 2. * stressM.at(5) * gamma113;
250  answer2.at(1, 6) = answer2.at(6, 1) = stressM.at(6) * ( m - 1 ) * pow(lambda1, m - 2);
251  answer2.at(2, 4) = answer2.at(4, 2) = 2. * stressM.at(4) * gamma113;
252  answer2.at(2, 6) = answer2.at(6, 2) = stressM.at(6) * ( m - 1 ) * pow(lambda1, m - 2);
253  answer2.at(3, 4) = answer2.at(4, 3) = 2. * stressM.at(4) * gamma331;
254  answer2.at(3, 5) = answer2.at(5, 3) = 2. * stressM.at(5) * gamma331;
255  answer2.at(4, 5) = answer2.at(5, 4) = stressM.at(6) * gamma113;
256  answer2.at(4, 6) = answer2.at(6, 4) = stressM.at(5) * gamma113;
257  answer2.at(5, 6) = answer2.at(6, 5) = stressM.at(4) * gamma113;
258  } else if ( lambda2 == lambda3 ) {
259  gamma23 = 1. / 2. * lambda2P;
260  gamma12 = ( E1 - E2 ) / ( lambda1 - lambda2 );
261  gamma13 = ( E1 - E3 ) / ( lambda1 - lambda3 );
262  gamma112 = ( pow(lambda1, m - 1) * ( lambda1 - lambda2 ) - 2 * ( E1 - E2 ) ) / ( ( lambda1 - lambda2 ) * ( lambda1 - lambda2 ) );
263  gamma221 = ( pow(lambda2, m - 1) * ( lambda2 - lambda1 ) - 2 * ( E2 - E1 ) ) / ( ( lambda2 - lambda1 ) * ( lambda2 - lambda1 ) );
264 
265  answer2.at(1, 1) = 2 * stressM.at(1) * ( m - 1 ) * pow(lambda1, m - 2);
266  answer2.at(2, 2) = 2 * stressM.at(2) * ( m - 1 ) * pow(lambda2, m - 2);
267  answer2.at(3, 3) = 2 * stressM.at(3) * ( m - 1 ) * pow(lambda3, m - 2);
268 
269  answer2.at(4, 4) = 1. / 2. * ( stressM.at(2) + stressM.at(3) ) * ( m - 1 ) * pow(lambda2, m - 2);
270  answer2.at(5, 5) = stressM.at(1) * gamma112 + stressM.at(3) * gamma221;
271  answer2.at(6, 6) = stressM.at(1) * gamma112 + stressM.at(2) * gamma221;
272 
273  answer2.at(1, 5) = answer2.at(5, 1) = 2. * stressM.at(5) * gamma112;
274  answer2.at(1, 6) = answer2.at(6, 1) = 2. * stressM.at(6) * gamma112;
275  answer2.at(2, 4) = answer2.at(4, 2) = stressM.at(4) * ( m - 1 ) * pow(lambda2, m - 2);
276  answer2.at(2, 6) = answer2.at(6, 2) = 2. * stressM.at(6) * gamma221;
277  answer2.at(3, 4) = answer2.at(4, 3) = stressM.at(4) * ( m - 1 ) * pow(lambda2, m - 2);
278  answer2.at(3, 5) = answer2.at(5, 3) = 2. * stressM.at(5) * gamma221;
279  answer2.at(4, 5) = answer2.at(5, 4) = stressM.at(6) * gamma221;
280  answer2.at(4, 6) = answer2.at(6, 4) = stressM.at(5) * gamma221;
281  answer2.at(5, 6) = answer2.at(6, 5) = stressM.at(4) * gamma221;
282  } else if ( lambda1 == lambda3 ) {
283  gamma13 = 1. / 2. * lambda1P;
284  gamma12 = ( E1 - E2 ) / ( lambda1 - lambda2 );
285  gamma23 = ( E2 - E3 ) / ( lambda2 - lambda3 );
286  gamma223 = ( pow(lambda2, m - 1) * ( lambda2 - lambda3 ) - 2 * ( E2 - E3 ) ) / ( ( lambda2 - lambda3 ) * ( lambda2 - lambda3 ) );
287  gamma332 = ( pow(lambda3, m - 1) * ( lambda3 - lambda2 ) - 2 * ( E3 - E2 ) ) / ( ( lambda3 - lambda2 ) * ( lambda3 - lambda2 ) );
288 
289  answer2.at(1, 1) = 2 * stressM.at(1) * ( m - 1 ) * pow(lambda1, m - 2);
290  answer2.at(2, 2) = 2 * stressM.at(2) * ( m - 1 ) * pow(lambda2, m - 2);
291  answer2.at(3, 3) = 2 * stressM.at(3) * ( m - 1 ) * pow(lambda3, m - 2);
292 
293  answer2.at(4, 4) = stressM.at(2) * gamma223 + stressM.at(3) * gamma332;
294  answer2.at(5, 5) = 1. / 2. * ( stressM.at(1) + stressM.at(3) ) * ( m - 1 ) * pow(lambda3, m - 2);
295  answer2.at(6, 6) = stressM.at(1) * gamma332 + stressM.at(2) * gamma223;
296 
297  answer2.at(1, 5) = answer2.at(5, 1) = stressM.at(5) * ( m - 1 ) * pow(lambda3, m - 2);
298  answer2.at(1, 6) = answer2.at(6, 1) = 2. * stressM.at(6) * gamma332;
299  answer2.at(2, 4) = answer2.at(4, 2) = 2. * stressM.at(4) * gamma223;
300  answer2.at(2, 6) = answer2.at(6, 2) = 2. * stressM.at(4) * gamma223;
301  answer2.at(3, 4) = answer2.at(4, 3) = 2. * stressM.at(4) * gamma332;
302  answer2.at(3, 5) = answer2.at(5, 3) = stressM.at(5) * ( m - 1 ) * pow(lambda3, m - 2);
303  answer2.at(4, 5) = answer2.at(5, 4) = stressM.at(6) * gamma332;
304  answer2.at(4, 6) = answer2.at(6, 4) = stressM.at(5) * gamma332;
305  answer2.at(5, 6) = answer2.at(6, 5) = stressM.at(4) * gamma332;
306  } else { //three different eigenvalues
307  gamma12 = ( E1 - E2 ) / ( lambda1 - lambda2 );
308  gamma13 = ( E1 - E3 ) / ( lambda1 - lambda3 );
309  gamma23 = ( E2 - E3 ) / ( lambda2 - lambda3 );
310 
311  gamma112 = ( pow(lambda1, m - 1) * ( lambda1 - lambda2 ) - 2 * ( E1 - E2 ) ) / ( ( lambda1 - lambda2 ) * ( lambda1 - lambda2 ) );
312  gamma221 = ( pow(lambda2, m - 1) * ( lambda2 - lambda1 ) - 2 * ( E2 - E1 ) ) / ( ( lambda2 - lambda1 ) * ( lambda2 - lambda1 ) );
313  gamma113 = ( pow(lambda1, m - 1) * ( lambda1 - lambda3 ) - 2 * ( E1 - E3 ) ) / ( ( lambda1 - lambda3 ) * ( lambda1 - lambda3 ) );
314  gamma331 = ( pow(lambda3, m - 1) * ( lambda3 - lambda1 ) - 2 * ( E3 - E1 ) ) / ( ( lambda3 - lambda1 ) * ( lambda3 - lambda1 ) );
315  gamma223 = ( pow(lambda2, m - 1) * ( lambda2 - lambda3 ) - 2 * ( E2 - E3 ) ) / ( ( lambda2 - lambda3 ) * ( lambda2 - lambda3 ) );
316  gamma332 = ( pow(lambda3, m - 1) * ( lambda3 - lambda2 ) - 2 * ( E3 - E2 ) ) / ( ( lambda3 - lambda2 ) * ( lambda3 - lambda2 ) );
317 
318  gamma = ( lambda1 * ( E2 - E3 ) + lambda2 * ( E3 - E1 ) + lambda3 * ( E1 - E2 ) ) / ( ( lambda1 - lambda2 ) * ( lambda2 - lambda3 ) * ( lambda3 - lambda1 ) );
319 
320  answer2.at(1, 1) = 2 * stressM.at(1) * ( m - 1 ) * pow(lambda1, m - 2);
321  answer2.at(2, 2) = 2 * stressM.at(2) * ( m - 1 ) * pow(lambda2, m - 2);
322  answer2.at(3, 3) = 2 * stressM.at(3) * ( m - 1 ) * pow(lambda3, m - 2);
323 
324  answer2.at(4, 4) = stressM.at(2) * gamma223 + stressM.at(3) * gamma332;
325  answer2.at(5, 5) = stressM.at(1) * gamma113 + stressM.at(3) * gamma331;
326  answer2.at(6, 6) = stressM.at(1) * gamma112 + stressM.at(2) * gamma221;
327 
328  answer2.at(1, 5) = answer2.at(5, 1) = 2. * stressM.at(5) * gamma113;
329  answer2.at(1, 6) = answer2.at(6, 1) = 2. * stressM.at(6) * gamma112;
330  answer2.at(2, 4) = answer2.at(4, 2) = 2. * stressM.at(4) * gamma223;
331  answer2.at(2, 6) = answer2.at(6, 2) = 2. * stressM.at(6) * gamma221;
332  answer2.at(3, 4) = answer2.at(4, 3) = 2. * stressM.at(4) * gamma332;
333  answer2.at(3, 5) = answer2.at(5, 3) = 2. * stressM.at(5) * gamma331;
334  answer2.at(4, 5) = answer2.at(5, 4) = 2. * stressM.at(6) * gamma;
335  answer2.at(4, 6) = answer2.at(6, 4) = 2. * stressM.at(5) * gamma;
336  answer2.at(5, 6) = answer2.at(6, 5) = 2. * stressM.at(4) * gamma;
337  }
338 
339 
340  answer1.at(1, 1) = lambda1P;
341  answer1.at(2, 2) = lambda2P;
342  answer1.at(3, 3) = lambda3P;
343  answer1.at(4, 4) = gamma23;
344  answer1.at(5, 5) = gamma13;
345  answer1.at(6, 6) = gamma12;
346 }
347 
348 void
349 LargeStrainMasterMaterial :: give3dMaterialStiffnessMatrix_dPdF(FloatMatrix &answer, MatResponseMode mode, GaussPoint *gp, TimeStep *tStep)
350 {
351  LargeStrainMasterMaterialStatus *status = static_cast< LargeStrainMasterMaterialStatus * >( this->giveStatus(gp) );
352  StructuralMaterial *sMat = static_cast< StructuralMaterial * >( domain->giveMaterial(slaveMat) );
353  FloatMatrix stiffness;
354 
355  sMat->give3dMaterialStiffnessMatrix(stiffness, mode, gp, tStep);
356  FloatMatrix junk;
358  stiffness.at(1, 4) = 2. * stiffness.at(1, 4);
359  stiffness.at(4, 1) = 2. * stiffness.at(4, 1);
360  stiffness.at(1, 5) = 2. * stiffness.at(1, 5);
361  stiffness.at(5, 1) = 2. * stiffness.at(5, 1);
362  stiffness.at(1, 6) = 2. * stiffness.at(1, 6);
363  stiffness.at(6, 1) = 2. * stiffness.at(6, 1);
364  stiffness.at(2, 4) = 2. * stiffness.at(2, 4);
365  stiffness.at(4, 2) = 2. * stiffness.at(4, 2);
366  stiffness.at(2, 5) = 2. * stiffness.at(2, 5);
367  stiffness.at(5, 2) = 2. * stiffness.at(5, 2);
368  stiffness.at(2, 6) = 2. * stiffness.at(2, 6);
369  stiffness.at(6, 2) = 2. * stiffness.at(6, 2);
370  stiffness.at(3, 4) = 2. * stiffness.at(3, 4);
371  stiffness.at(4, 3) = 2. * stiffness.at(4, 3);
372  stiffness.at(3, 5) = 2. * stiffness.at(3, 5);
373  stiffness.at(5, 3) = 2. * stiffness.at(5, 3);
374  stiffness.at(3, 6) = 2. * stiffness.at(3, 6);
375  stiffness.at(6, 3) = 2. * stiffness.at(6, 3);
376  stiffness.at(4, 4) = 4. * stiffness.at(4, 4);
377  stiffness.at(4, 5) = 4. * stiffness.at(4, 5);
378  stiffness.at(5, 4) = 4. * stiffness.at(5, 4);
379  stiffness.at(4, 6) = 4. * stiffness.at(4, 6);
380  stiffness.at(6, 4) = 4. * stiffness.at(6, 4);
381  stiffness.at(5, 5) = 4. * stiffness.at(5, 5);
382  stiffness.at(5, 6) = 4. * stiffness.at(5, 6);
383  stiffness.at(6, 5) = 4. * stiffness.at(6, 5);
384  stiffness.at(6, 6) = 4. * stiffness.at(6, 6);
386  junk.beProductOf( stiffness, status->givePmatrix() );
387  stiffness.beProductOf(status->givePmatrix(), junk);
388  stiffness.add( status->giveTLmatrix() );
389 
390  this->convert_dSdE_2_dPdF(answer, stiffness, status->giveTempStressVector(), status->giveTempFVector(), _3dMat);
391 }
392 
393 
394 int
395 LargeStrainMasterMaterial :: giveIPValue(FloatArray &answer, GaussPoint *gp, InternalStateType type, TimeStep *tStep)
396 {
397  LargeStrainMasterMaterialStatus *status = static_cast< LargeStrainMasterMaterialStatus * >( this->giveStatus(gp) );
398 
399  if ( type == IST_StressTensor ) {
400  answer = status->giveStressVector();
401  return 1;
402  } else if ( type == IST_StrainTensor ) {
403  answer = status->giveStrainVector();
404  return 1;
405  } else {
406  StructuralMaterial *sMat = static_cast< StructuralMaterial * >( domain->giveMaterial(slaveMat) );
407  return sMat->giveIPValue(answer, gp, type, tStep);
408  }
409 }
410 
411 
412 //=============================================================================
413 
414 LargeStrainMasterMaterialStatus :: LargeStrainMasterMaterialStatus(int n, Domain *d, GaussPoint *g, int s) : StructuralMaterialStatus(n, d, g),
415  Pmatrix(6, 6),
416  TLmatrix(6, 6),
417  transformationMatrix(6, 6),
418  slaveMat(s)
419 {
420  Pmatrix.beUnitMatrix();
421 }
422 
423 
424 LargeStrainMasterMaterialStatus :: ~LargeStrainMasterMaterialStatus()
425 { }
426 
427 
428 void
429 LargeStrainMasterMaterialStatus :: printOutputAt(FILE *file, TimeStep *tStep)
430 {
431  StructuralMaterial *sMat = static_cast< StructuralMaterial * >( domain->giveMaterial(slaveMat) );
432  MaterialStatus *mS = sMat->giveStatus(gp);
433 
434  mS->printOutputAt(file, tStep);
435  // StructuralMaterialStatus :: printOutputAt(file, tStep);
436 }
437 
438 
439 // initializes temporary variables based on their values at the previous equlibrium state
440 void LargeStrainMasterMaterialStatus :: initTempStatus()
441 {
442  StructuralMaterial *sMat = static_cast< StructuralMaterial * >( domain->giveMaterial(slaveMat) );
443  MaterialStatus *mS = sMat->giveStatus(gp);
444  mS->initTempStatus();
445  //StructuralMaterialStatus :: initTempStatus();
446 }
447 
448 
449 // updates internal variables when equilibrium is reached
450 void
451 LargeStrainMasterMaterialStatus :: updateYourself(TimeStep *tStep)
452 {
453  StructuralMaterial *sMat = static_cast< StructuralMaterial * >( domain->giveMaterial(slaveMat) );
454  MaterialStatus *mS = sMat->giveStatus(gp);
455  mS->updateYourself(tStep);
456  // StructuralMaterialStatus :: updateYourself(tStep);
457 }
458 
459 
460 // saves full information stored in this status
461 // temporary variables are NOT stored
462 contextIOResultType
463 LargeStrainMasterMaterialStatus :: saveContext(DataStream &stream, ContextMode mode, void *obj)
464 {
465  contextIOResultType iores;
466  StructuralMaterial *sMat = dynamic_cast< StructuralMaterial * >( domain->giveMaterial(slaveMat) );
467  MaterialStatus *mS = sMat->giveStatus(gp);
468  // save parent class status
469  if ( ( iores = mS->saveContext(stream, mode, obj) ) != CIO_OK ) {
470  THROW_CIOERR(iores);
471  }
472 
473  // write raw data
474 
475  return CIO_OK;
476 }
477 
478 
479 contextIOResultType
480 LargeStrainMasterMaterialStatus :: restoreContext(DataStream &stream, ContextMode mode, void *obj)
481 //
482 // restores full information stored in stream to this Status
483 //
484 {
485  contextIOResultType iores;
486 
487  // read parent class status
488  if ( ( iores = StructuralMaterialStatus :: restoreContext(stream, mode, obj) ) != CIO_OK ) {
489  THROW_CIOERR(iores);
490  }
491 
492  return CIO_OK; // return succes
493 }
494 } // end namespace oofem
oofem::InternalStateType
InternalStateType
Type representing the physical meaning of element or constitutive model internal variable.
Definition: internalstatetype.h:187
oofem::MaterialStatus::updateYourself
virtual void updateYourself(TimeStep *)
Update equilibrium history variables according to temp-variables.
Definition: matstatus.h:108
oofem::Material::giveStatus
virtual MaterialStatus * giveStatus(GaussPoint *gp) const
Returns material status of receiver in given integration point.
Definition: material.C:244
oofem::LargeStrainMasterMaterial::give3dMaterialStiffnessMatrix_dPdF
virtual void give3dMaterialStiffnessMatrix_dPdF(FloatMatrix &answer, MatResponseMode, GaussPoint *gp, TimeStep *tStep)
Definition: lsmastermat.C:349
oofem::IntegrationPointStatus::gp
GaussPoint * gp
Associated integration point.
Definition: integrationpointstatus.h:57
oofem::Domain
Class and object Domain.
Definition: domain.h:115
oofem::LargeStrainMasterMaterial::giveFirstPKStressVector_3d
virtual void giveFirstPKStressVector_3d(FloatArray &answer, GaussPoint *gp, const FloatArray &vF, TimeStep *tStep)
Default implementation relies on giveRealStressVector for second Piola-Kirchoff stress.
Definition: lsmastermat.C:82
oofem::FEMComponent::domain
Domain * domain
Link to domain object, useful for communicating with other FEM components.
Definition: femcmpnn.h:82
oofem::LargeStrainMasterMaterialStatus::printOutputAt
virtual void printOutputAt(FILE *file, TimeStep *tStep)
Print receiver&#39;s output to given stream.
Definition: lsmastermat.C:429
oofem::LargeStrainMasterMaterial::m
double m
Specifies the strain tensor.
Definition: lsmastermat.h:72
oofem::TL
Total Lagrange.
Definition: fmode.h:44
oofem::DataStream
The purpose of DataStream abstract class is to allow to store/restore context to different streams...
Definition: datastream.h:54
oofem::StructuralMaterial::giveIPValue
virtual int giveIPValue(FloatArray &answer, GaussPoint *gp, InternalStateType type, TimeStep *tStep)
Returns the integration point corresponding value in Reduced form.
Definition: structuralmaterial.C:2115
oofem::FloatArray::at
double & at(int i)
Coefficient access function.
Definition: floatarray.h:131
oofem::StructuralMaterialStatus
This class implements a structural material status information.
Definition: structuralms.h:65
oofem::LargeStrainMasterMaterial::initializeFrom
virtual IRResultType initializeFrom(InputRecord *ir)
Initializes receiver according to object description stored in input record.
Definition: lsmastermat.C:63
oofem::FEMComponent::saveContext
virtual contextIOResultType saveContext(DataStream &stream, ContextMode mode, void *obj=NULL)
Stores receiver state to output stream.
Definition: femcmpnn.C:51
oofem::LargeStrainMasterMaterialStatus::giveTLmatrix
const FloatMatrix & giveTLmatrix()
Definition: lsmastermat.h:122
oofem::MatResponseMode
MatResponseMode
Describes the character of characteristic material matrix.
Definition: matresponsemode.h:64
THROW_CIOERR
#define THROW_CIOERR(e)
Definition: contextioerr.h:61
oofem::FloatMatrix::beMatrixForm
void beMatrixForm(const FloatArray &aArray)
Definition: floatmatrix.C:1657
oofem::FloatMatrix::jaco_
bool jaco_(FloatArray &eval, FloatMatrix &v, int nf)
Computes eigenvalues and eigenvectors of receiver (must be symmetric) The receiver is preserved...
Definition: floatmatrix.C:1912
oofem::LargeStrainMasterMaterialStatus::setTLmatrix
void setTLmatrix(const FloatMatrix &values)
Definition: lsmastermat.h:126
oofem::MaterialStatus::initTempStatus
virtual void initTempStatus()
Initializes the temporary internal variables, describing the current state according to previously re...
Definition: matstatus.h:103
oofem::Domain::giveMaterial
Material * giveMaterial(int n)
Service for accessing particular domain material model.
Definition: domain.C:281
oofem::StructuralMaterial::give3dMaterialStiffnessMatrix
virtual void give3dMaterialStiffnessMatrix(FloatMatrix &answer, MatResponseMode mode, GaussPoint *gp, TimeStep *tStep)
Computes full 3d material stiffness matrix at given integration point, time, respecting load history ...
Definition: structuralmaterial.h:344
oofem::LargeStrainMasterMaterial::CreateStatus
virtual MaterialStatus * CreateStatus(GaussPoint *gp) const
Creates new copy of associated status and inserts it into given integration point.
Definition: lsmastermat.C:75
oofem::FloatArray::beProductOf
void beProductOf(const FloatMatrix &aMatrix, const FloatArray &anArray)
Receiver becomes the result of the product of aMatrix and anArray.
Definition: floatarray.C:676
_IFT_LargeStrainMasterMaterial_m
#define _IFT_LargeStrainMasterMaterial_m
Definition: lsmastermat.h:48
oofem::FloatMatrix::at
double at(int i, int j) const
Coefficient access function.
Definition: floatmatrix.h:176
_IFT_LargeStrainMasterMaterial_slaveMat
#define _IFT_LargeStrainMasterMaterial_slaveMat
Definition: lsmastermat.h:49
oofem::StructuralMaterialStatus::restoreContext
virtual contextIOResultType restoreContext(DataStream &stream, ContextMode mode, void *obj=NULL)
Restores the receiver state previously written in stream.
Definition: structuralms.C:157
oofem::StructuralMaterial::giveRealStressVector_3d
virtual void giveRealStressVector_3d(FloatArray &answer, GaussPoint *gp, const FloatArray &reducedE, TimeStep *tStep)
Default implementation relies on giveRealStressVector for second Piola-Kirchoff stress.
Definition: structuralmaterial.C:112
oofem::MaterialStatus
Abstract base class representing a material status information.
Definition: matstatus.h:84
oofem::LargeStrainMasterMaterialStatus::setPmatrix
void setPmatrix(const FloatMatrix &values)
Definition: lsmastermat.h:125
oofem::FloatArray
Class representing vector of real numbers.
Definition: floatarray.h:82
oofem::StructuralMaterial::convert_dSdE_2_dPdF
void convert_dSdE_2_dPdF(FloatMatrix &answer, const FloatMatrix &dSdE, const FloatArray &S, const FloatArray &F, MaterialMode matMode)
Definition: structuralmaterial.C:434
oofem::FloatMatrix
Implementation of matrix containing floating point numbers.
Definition: floatmatrix.h:94
oofem::LargeStrainMasterMaterialStatus::initTempStatus
virtual void initTempStatus()
Initializes the temporary internal variables, describing the current state according to previously re...
Definition: lsmastermat.C:440
oofem::IRResultType
IRResultType
Type defining the return values of InputRecord reading operations.
Definition: irresulttype.h:47
oofem::StructuralMaterialStatus::letTempStressVectorBe
void letTempStressVectorBe(const FloatArray &v)
Assigns tempStressVector to given vector v.
Definition: structuralms.h:135
oofem::StructuralMaterialStatus::giveStressVector
const FloatArray & giveStressVector() const
Returns the const pointer to receiver&#39;s stress vector.
Definition: structuralms.h:107
oofem::LargeStrainMasterMaterial::LargeStrainMasterMaterial
LargeStrainMasterMaterial(int n, Domain *d)
Definition: lsmastermat.C:52
oofem::FloatMatrix::resize
void resize(int rows, int cols)
Checks size of receiver towards requested bounds.
Definition: floatmatrix.C:1358
oofem::InputRecord
Class representing the general Input Record.
Definition: inputrecord.h:101
oofem::FloatMatrix::add
void add(const FloatMatrix &a)
Adds matrix to the receiver.
Definition: floatmatrix.C:1023
oofem::FloatMatrix::beTProductOf
void beTProductOf(const FloatMatrix &a, const FloatMatrix &b)
Assigns to the receiver product of .
Definition: floatmatrix.C:367
oofem::StructuralMaterialStatus::giveTempFVector
const FloatArray & giveTempFVector() const
Returns the const pointer to receiver&#39;s temporary deformation gradient vector.
Definition: structuralms.h:123
oofem::ContextMode
long ContextMode
Context mode (mask), defining the type of information written/read to/from context.
Definition: contextmode.h:43
oofem::StructuralMaterial
Abstract base class for all "structural" constitutive models.
Definition: structuralmaterial.h:113
oofem::FloatMatrix::zero
void zero()
Zeroes all coefficient of receiver.
Definition: floatmatrix.C:1326
oofem::FEMComponent::giveDomain
Domain * giveDomain() const
Definition: femcmpnn.h:100
oofem::FloatMatrix::beUnitMatrix
void beUnitMatrix()
Sets receiver to unity matrix.
Definition: floatmatrix.C:1332
oofem::REGISTER_Material
REGISTER_Material(DummyMaterial)
oofem::FloatArray::beSymVectorFormOfStrain
void beSymVectorFormOfStrain(const FloatMatrix &aMatrix)
Definition: floatarray.C:1014
oofem::LargeStrainMasterMaterialStatus::LargeStrainMasterMaterialStatus
LargeStrainMasterMaterialStatus(int n, Domain *d, GaussPoint *g, int s)
Definition: lsmastermat.C:414
IR_GIVE_OPTIONAL_FIELD
#define IR_GIVE_OPTIONAL_FIELD(__ir, __value, __id)
Macro facilitating the use of input record reading methods.
Definition: inputrecord.h:78
oofem::FloatMatrix::beProductOf
void beProductOf(const FloatMatrix &a, const FloatMatrix &b)
Assigns to the receiver product of .
Definition: floatmatrix.C:337
oofem::StructuralMaterialStatus::giveTempStressVector
const FloatArray & giveTempStressVector() const
Returns the const pointer to receiver&#39;s temporary stress vector.
Definition: structuralms.h:117
oofem
the oofem namespace is to define a context or scope in which all oofem names are defined.
oofem::LargeStrainMasterMaterial::constructTransformationMatrix
void constructTransformationMatrix(FloatMatrix &answer, const FloatMatrix &eigenVectors)
transformation matrices
Definition: lsmastermat.C:156
oofem::MaterialStatus::printOutputAt
virtual void printOutputAt(FILE *file, TimeStep *tStep)
Print receiver&#39;s output to given stream.
Definition: matstatus.h:97
oofem::LargeStrainMasterMaterial::constructL1L2TransformationMatrices
void constructL1L2TransformationMatrices(FloatMatrix &answer1, FloatMatrix &answer2, const FloatArray &eigenValues, FloatArray &stress, double E1, double E2, double E3)
Definition: lsmastermat.C:204
oofem::LargeStrainMasterMaterial::giveIPValue
virtual int giveIPValue(FloatArray &answer, GaussPoint *gp, InternalStateType type, TimeStep *tStep)
Returns the integration point corresponding value in Reduced form.
Definition: lsmastermat.C:395
oofem::Material::initTempStatus
virtual void initTempStatus(GaussPoint *gp)
Initializes temporary variables stored in integration point status at the beginning of new time step...
Definition: material.C:267
oofem::LargeStrainMasterMaterialStatus::restoreContext
virtual contextIOResultType restoreContext(DataStream &stream, ContextMode mode, void *obj=NULL)
Restores the receiver state previously written in stream.
Definition: lsmastermat.C:480
oofem::StructuralMaterialStatus::giveStrainVector
const FloatArray & giveStrainVector() const
Returns the const pointer to receiver&#39;s strain vector.
Definition: structuralms.h:105
oofem::GaussPoint
Class representing integration point in finite element program.
Definition: gausspoint.h:93
oofem::LargeStrainMasterMaterial::slaveMat
int slaveMat
&#39;slave&#39; material model number.
Definition: lsmastermat.h:70
oofem::LargeStrainMasterMaterialStatus::saveContext
virtual contextIOResultType saveContext(DataStream &stream, ContextMode mode, void *obj=NULL)
Stores receiver state to output stream.
Definition: lsmastermat.C:463
oofem::TimeStep
Class representing solution step.
Definition: timestep.h:80
oofem::LargeStrainMasterMaterialStatus::updateYourself
virtual void updateYourself(TimeStep *)
Update equilibrium history variables according to temp-variables.
Definition: lsmastermat.C:451
This page is part of the OOFEM documentation. Copyright (c) 2011 Borek Patzak
Project e-mail: info@oofem.org
Generated at Tue Jan 2 2018 20:07:29 for OOFEM by doxygen 1.8.11 written by Dimitri van Heesch, © 1997-2011