Line data Source code
1 : // MASTODON includes
2 : #include "ComputeFPIsolatorElasticity.h"
3 :
4 : // MOOSE includes
5 : #include "MooseMesh.h"
6 : #include "Assembly.h"
7 : #include "NonlinearSystem.h"
8 : #include "MooseVariable.h"
9 : #include "MathUtils.h"
10 :
11 : // libmesh includes
12 : #include "libmesh/quadrature.h"
13 : #include "libmesh/utility.h"
14 :
15 : registerMooseObject("MastodonApp", ComputeFPIsolatorElasticity);
16 :
17 : InputParameters
18 96 : ComputeFPIsolatorElasticity::validParams()
19 : {
20 96 : InputParameters params = Material::validParams();
21 96 : params.addClassDescription("Compute the forces and the stiffness matrix for a single concave "
22 : "Friction Pendulum isolator element.");
23 : // Switches
24 192 : params.addRequiredParam<bool>(
25 : "pressure_dependent",
26 : "Switch for modeling friction dependence on the instantaneous pressure.");
27 192 : params.addRequiredParam<bool>(
28 : "temperature_dependent",
29 : "Switch for modeling friction dependence on the instantaneous temperature.");
30 192 : params.addRequiredParam<bool>(
31 : "velocity_dependent",
32 : "Switch for modeling friction dependence on the instantaneous sliding velocity.");
33 : // Material properties
34 192 : params.addRequiredRangeCheckedParam<Real>(
35 : "mu_ref", "mu_ref>0.0", "Reference co-efficient of friction.");
36 192 : params.addRequiredRangeCheckedParam<Real>("p_ref", "p_ref>0.0", "Reference axial pressure.");
37 192 : params.addRequiredRangeCheckedParam<Real>(
38 : "diffusivity", "diffusivity>0.0", "Thermal diffusivity of steel.");
39 192 : params.addRequiredRangeCheckedParam<Real>(
40 : "conductivity", "conductivity>0.0", "Thermal conductivity of steel.");
41 192 : params.addRequiredRangeCheckedParam<Real>("a", "a>0.0", "Rate parameter.");
42 192 : params.addRequiredRangeCheckedParam<Real>(
43 : "r_eff", "r_eff>0.0", "Effective radius of curvature of sliding surface.");
44 192 : params.addRequiredRangeCheckedParam<Real>(
45 : "r_contact", "r_contact>0.0", "Radius of contact area at sliding surface.");
46 192 : params.addRequiredRangeCheckedParam<Real>(
47 : "uy", "uy>0.0", "Yield displacement of the bearing in shear.");
48 192 : params.addRequiredRangeCheckedParam<Real>(
49 : "gamma", "gamma>0.0", "Gamma parameter of Newmark algorithm.");
50 192 : params.addRequiredRangeCheckedParam<Real>(
51 : "beta", "beta>0.0", "Beta parameter of Newmark algorithm.");
52 192 : params.addRequiredRangeCheckedParam<unsigned int>(
53 : "unit",
54 : "8 > unit > 0",
55 : "Tag for conversion in the pressure factor computation "
56 : "when different unit systems are used. Enter "
57 : "1 for N m s C; "
58 : "2 for kN m s C; "
59 : "3 for N mm s C; "
60 : "4 for kN mm s C; "
61 : "5 for lb in s C; "
62 : "6 for kip in s C; "
63 : "7 for lb ft s C; "
64 : "8 for kip ft s C. ");
65 192 : params.addParam<Real>(
66 : "k_x",
67 192 : 10e13,
68 : "Stiffness of the bearing in translation along x axis. Defaults to 10e13 (SI units)");
69 192 : params.addParam<Real>(
70 : "k_xx",
71 192 : 10e13,
72 : "Stiffness of the bearing in rotation about x axis. Defaults to 10e13 (SI units) ");
73 192 : params.addParam<Real>(
74 : "k_yy",
75 192 : 10e13,
76 : "Stiffness of the bearing in rotation about y axis. Defaults to 10e13 (SI units) ");
77 192 : params.addParam<Real>(
78 : "k_zz",
79 192 : 10e13,
80 : "Stiffness of the bearing in rotation about z axis. Defaults to 10e13 (SI units) ");
81 192 : params.addParam<Real>(
82 192 : "tol", 1e-8, "Convergence tolerance to satisfy equilibrium of element. Defaults to 1e-8");
83 192 : params.addParam<Real>(
84 192 : "maxiter", 100, "Maximum number of iterations for equilibrium. Defaults to 100.");
85 192 : params.set<MooseEnum>("constant_on") = "ELEMENT";
86 96 : return params;
87 0 : }
88 :
89 72 : ComputeFPIsolatorElasticity::ComputeFPIsolatorElasticity(const InputParameters & parameters)
90 : : Material(parameters),
91 72 : _pressure_dependent(getParam<bool>("pressure_dependent")),
92 144 : _temperature_dependent(getParam<bool>("temperature_dependent")),
93 144 : _velocity_dependent(getParam<bool>("velocity_dependent")),
94 144 : _mu_ref(getParam<Real>("mu_ref")),
95 144 : _p_ref(getParam<Real>("p_ref")),
96 144 : _diffusivity(getParam<Real>("diffusivity")),
97 144 : _conductivity(getParam<Real>("conductivity")),
98 144 : _a(getParam<Real>("a")),
99 144 : _r_eff(getParam<Real>("r_eff")),
100 144 : _r_contact(getParam<Real>("r_contact")),
101 144 : _uy(getParam<Real>("uy")),
102 144 : _unit(getParam<unsigned int>("unit")),
103 144 : _gamma(getParam<Real>("gamma")),
104 144 : _beta(getParam<Real>("beta")),
105 144 : _k_x(getParam<Real>("k_x")),
106 144 : _k_xx(getParam<Real>("k_xx")),
107 144 : _k_yy(getParam<Real>("k_yy")),
108 144 : _k_zz(getParam<Real>("k_zz")),
109 144 : _tol(getParam<Real>("tol")),
110 144 : _maxiter(getParam<Real>("maxiter")),
111 72 : _sD(0.5),
112 144 : _local_def(getMaterialPropertyByName<ColumnMajorMatrix>("local_deformations")),
113 144 : _basic_def(getMaterialPropertyByName<ColumnMajorMatrix>("deformations")),
114 144 : _basic_def_old(getMaterialPropertyByName<ColumnMajorMatrix>("old_deformations")),
115 144 : _basic_vel(getMaterialPropertyByName<ColumnMajorMatrix>("deformation_rates")),
116 144 : _basic_vel_old(getMaterialPropertyByName<ColumnMajorMatrix>("old_deformation_rates")),
117 72 : _Fb(declareProperty<ColumnMajorMatrix>("basic_forces")),
118 72 : _Fl(declareProperty<ColumnMajorMatrix>("local_forces")),
119 72 : _Fg(declareProperty<ColumnMajorMatrix>("global_forces")),
120 72 : _Kb(declareProperty<ColumnMajorMatrix>("basic_stiffness_matrix")),
121 72 : _Kl(declareProperty<ColumnMajorMatrix>("local_stiffness_matrix")),
122 72 : _Kg(declareProperty<ColumnMajorMatrix>("global_stiffness_matrix")),
123 144 : _total_gl(getMaterialPropertyByName<ColumnMajorMatrix>("total_global_to_local_transformation")),
124 144 : _total_lb(getMaterialPropertyByName<ColumnMajorMatrix>("total_local_to_basic_transformation")),
125 144 : _length(getMaterialPropertyByName<Real>("initial_isolator_length")),
126 72 : _pi(libMesh::pi),
127 72 : _ubPlastic(declareProperty<RealVectorValue>("plastic displacements in basic system")),
128 216 : _ubPlastic_old(getMaterialPropertyOld<RealVectorValue>("plastic displacements in basic system"))
129 :
130 : {
131 : // Calculation of lateral stiffness of elastic component
132 72 : _mass_slider = _p_ref * _pi * _r_contact * _r_contact / 9.81; // Mass of the slider
133 72 : _k0 = _mass_slider * 9.81 * _mu_ref / _uy; // Initial elastic stiffness in lateral direction
134 :
135 : // Other parameters
136 72 : _kpF = 1;
137 72 : _ktF = 1;
138 72 : _kvF = 1;
139 72 : _temperature_center = 0;
140 72 : _heatflux_center = 0;
141 72 : _mu_adj = _mu_ref;
142 72 : _disp_currentstep = 0;
143 72 : _vel_currentstep = 0;
144 72 : _vec_time.push_back(0);
145 72 : _vec_heatflux.push_back(0);
146 72 : }
147 :
148 : void
149 40 : ComputeFPIsolatorElasticity::initQpStatefulProperties()
150 : {
151 40 : _ubPlastic[_qp].zero();
152 40 : }
153 :
154 : void
155 12316 : ComputeFPIsolatorElasticity::computeQpProperties()
156 : {
157 12316 : initializeFPIsolator();
158 :
159 : // Computing shear forces and stiffness terms
160 12316 : computeShear();
161 :
162 : // Compute forces in other directions
163 12316 : _Fb[_qp](0, 0) = _Kb[_qp](0, 0) * _basic_def[_qp](0, 0); // translation in basic x direction
164 12316 : _Fb[_qp](3, 0) = _Kb[_qp](3, 3) * _basic_def[_qp](3, 0); // rotation along basic x direction
165 12316 : _Fb[_qp](4, 0) = _Kb[_qp](4, 4) * _basic_def[_qp](4, 0); // rotation along basic y direction
166 12316 : _Fb[_qp](5, 0) = _Kb[_qp](5, 5) * _basic_def[_qp](5, 0); // rotation along basic z direction
167 :
168 : // Finalize forces and stiffness matrix and convert them into global co-ordinate system
169 12316 : finalize();
170 12316 : }
171 :
172 : void
173 12316 : ComputeFPIsolatorElasticity::initializeFPIsolator()
174 : {
175 : // Initialize stiffness matrices
176 12316 : _Kb[_qp].reshape(6, 6);
177 12316 : _Kb[_qp].identity();
178 12316 : _Kb[_qp](0, 0) = _k_x; // axial stiffness along x axis
179 12316 : _Kb[_qp](1, 1) = _k0; // elastic lateral stiffness along y axis
180 12316 : _Kb[_qp](2, 2) = _k0; // elastic lateral stiffness along z axis
181 12316 : _Kb[_qp](3, 3) = _k_xx; // torsional stiffness about x axis
182 12316 : _Kb[_qp](4, 4) = _k_yy; // rotational stiffness about y axis
183 12316 : _Kb[_qp](5, 5) = _k_zz; // rotational stiffness about z axis
184 :
185 : // Initialize forces
186 12316 : _Fb[_qp].reshape(6, 1); // forces in the basic system
187 12316 : _Fb[_qp].zero();
188 12316 : _Fl[_qp].reshape(12, 1); // forces in local system (6+6 = 12dofs)
189 12316 : _Fl[_qp].zero();
190 12316 : _Fg[_qp] = _Fl[_qp]; // forces in global system
191 :
192 12316 : _Fb[_qp] = _Kb[_qp] * _basic_def[_qp];
193 12316 : }
194 :
195 : void
196 12316 : ComputeFPIsolatorElasticity::computeShear()
197 : {
198 12316 : if (_t != _vec_time[_vec_time.size() - 1])
199 : {
200 : // Update last element of time vector with current analysis time
201 4300 : _vec_time.push_back(_t);
202 4300 : _vec_heatflux.push_back(0);
203 : }
204 :
205 : // Compute resultant displacement in horizontal direction
206 12316 : Real resultantU = sqrt(_basic_def[_qp](1, 0) * _basic_def[_qp](1, 0) +
207 12316 : _basic_def[_qp](2, 0) * _basic_def[_qp](2, 0));
208 :
209 : // Calculate radii in basic y- and z- direction
210 12316 : Real r_y = sqrt(pow(_r_eff, 2) - pow(_basic_def[_qp](1, 0), 2));
211 12316 : Real r_z = sqrt(pow(_r_eff, 2) - pow(_basic_def[_qp](2, 0), 2));
212 :
213 : // Calculate velocities based on current deformations according to Newmark formulation
214 12316 : Real vel1 = _basic_vel_old[_qp](1, 0) +
215 12316 : (_gamma / _beta) * (((_basic_def[_qp](1, 0) - _basic_def_old[_qp](1, 0)) / _dt) -
216 12316 : (_basic_vel_old[_qp](1, 0)));
217 12316 : Real vel2 = _basic_vel_old[_qp](2, 0) +
218 12316 : (_gamma / _beta) * (((_basic_def[_qp](2, 0) - _basic_def_old[_qp](2, 0)) / _dt) -
219 12316 : (_basic_vel_old[_qp](2, 0)));
220 :
221 : // Calculate absolute velocity
222 : Real velAbs =
223 12316 : sqrt(pow(vel1 / r_y * _basic_def[_qp](1, 0) + vel2 / r_z * _basic_def[_qp](2, 0), 2) +
224 : pow(vel1, 2) + pow(vel2, 2));
225 :
226 : // Check for uplift
227 : bool is_uplift = 0;
228 : Real kFactUplift = 1.0E-6;
229 12316 : if (_Fb[_qp](0, 0) >= 0.0)
230 : {
231 30 : _Fb[_qp](0, 0) = -0.0001 * _Fb[_qp](0, 0);
232 : is_uplift = 1;
233 : }
234 :
235 : // Unit conversion for pressure to be used in the pressure factor computation
236 : Real p_Unit_Convert; // To convert to MPa
237 12316 : if (_unit == 1)
238 : {
239 : p_Unit_Convert = 0.000001;
240 : }
241 12316 : if (_unit == 2)
242 : {
243 : p_Unit_Convert = 0.001;
244 : }
245 12316 : if (_unit == 3)
246 : {
247 : p_Unit_Convert = 1.0;
248 : }
249 12316 : if (_unit == 4)
250 : {
251 : p_Unit_Convert = 1000.0;
252 : }
253 12316 : if (_unit == 5)
254 : {
255 : p_Unit_Convert = 0.006894;
256 : }
257 12316 : if (_unit == 6)
258 : {
259 : p_Unit_Convert = 6.894;
260 : }
261 12316 : if (_unit == 7)
262 : {
263 : p_Unit_Convert = 0.00004788;
264 : }
265 12316 : if (_unit == 8)
266 : {
267 : p_Unit_Convert = 0.04788;
268 : }
269 :
270 : // Calculate normal force
271 24632 : Real N = -_Fb[_qp](0, 0);
272 :
273 : // Calculate instantaneous pressure
274 12316 : Real Inst_Pressure = fabs(N) / (_pi * _r_contact * _r_contact);
275 :
276 : // Calculate displacement increment in the current step
277 12316 : _disp_currentstep = sqrt(pow((_basic_def[_qp](1, 0) - _basic_def_old[_qp](1, 0)), 2) +
278 12316 : pow((_basic_def[_qp](2, 0) - _basic_def_old[_qp](2, 0)), 2));
279 :
280 : // Calculate sliding velocity
281 12316 : _vel_currentstep = _disp_currentstep / _dt;
282 :
283 12316 : if (velAbs > 0.0)
284 5340 : _vel_currentstep = velAbs;
285 :
286 : // Calculate pressure factor
287 12316 : if (_pressure_dependent)
288 6716 : _kpF = pow(0.7, ((Inst_Pressure - _p_ref) * p_Unit_Convert / 50.0));
289 :
290 : // Calculate temperature factor
291 12316 : if (_temperature_dependent)
292 : {
293 : // Calculate flux
294 6716 : if (resultantU < (_r_contact * sqrt(_pi / 4)))
295 6716 : _vec_heatflux[_t_step] =
296 6716 : _mu_adj * (fabs(N) / (_pi * _r_contact * _r_contact)) * _vel_currentstep;
297 :
298 : else
299 0 : _vec_heatflux[_t_step] = 0.0;
300 :
301 6716 : _heatflux_center = _vec_heatflux[_t_step];
302 :
303 : // Temperature at the sliding surface
304 : Real temperature_change = 0.0;
305 :
306 : // Evaluate the temperature integral
307 6716 : if (_t_step >= 1 && resultantU > 0.0)
308 : {
309 : Real dt_analysis = 0.000001;
310 : Real tau = 0.000001;
311 :
312 913546 : for (int jTemp = 1; jTemp <= _t_step; jTemp++)
313 : {
314 910616 : dt_analysis = _dt;
315 910616 : tau = _vec_time[jTemp];
316 :
317 910616 : temperature_change =
318 : temperature_change +
319 910616 : sqrt(_diffusivity / _pi) * _vec_heatflux[_t_step - jTemp] * dt_analysis /
320 910616 : (_conductivity * sqrt(tau)); // Solving the integral for temperature increment
321 : }
322 :
323 2930 : _temperature_center =
324 2930 : 20.0 + temperature_change; // Temperature at the center of the sliding surface
325 :
326 2930 : _ktF =
327 2930 : 0.789 * (pow(0.7, (_temperature_center / 50.0)) + 0.4); // Update the temperature factor
328 : }
329 : }
330 :
331 : // Calculate velocity factor
332 12316 : if (_velocity_dependent)
333 6716 : _kvF = 1 - 0.5 * exp(-_a * _vel_currentstep);
334 :
335 : // Calculate adjusted co-efficient of friction
336 12316 : _mu_adj = _mu_ref * _kpF * _ktF * _kvF;
337 :
338 : // Calculate shear forces and stiffness in basic y- and z-direction
339 12316 : unsigned int iter = 0;
340 : _qbOld.zero();
341 :
342 : do
343 : {
344 : // Save shear forces of last iteration
345 17342 : _qbOld(0) = _Fb[_qp](1, 0);
346 17342 : _qbOld(1) = _Fb[_qp](2, 0);
347 :
348 : // Yield force
349 17342 : Real qYield = _mu_adj * N;
350 :
351 : // Calculate stiffness of elastic components
352 17342 : Real k2y = N / r_y;
353 17342 : Real k2z = N / r_z;
354 :
355 : // Calculate initial stiffness of hysteretic components
356 17342 : Real k0y = _k0 - k2y;
357 17342 : Real k0z = _k0 - k2z;
358 :
359 : // Account for uplift
360 : kFactUplift = 1.0E-6;
361 17342 : if (is_uplift)
362 : {
363 30 : k0y = kFactUplift * k0y;
364 30 : k0z = kFactUplift * k0z;
365 : }
366 :
367 : // Trial shear forces of hysteretic component
368 : _qTrial.zero();
369 17342 : _qTrial(0) = k0y * (_basic_def[_qp](1, 0) - _ubPlastic_old[_qp](0));
370 17342 : _qTrial(1) = k0z * (_basic_def[_qp](2, 0) - _ubPlastic_old[_qp](1));
371 :
372 : // Compute yield criterion of hysteretic component
373 17342 : Real qTrialNorm = _qTrial.norm();
374 17342 : Real Y = qTrialNorm - qYield;
375 :
376 : // Elastic step -> no updates required
377 17342 : if (Y <= 0.0)
378 : {
379 : // Set shear forces
380 9558 : _Fb[_qp](1, 0) = _qTrial(0) - N * _local_def[_qp](5, 0) + k2y * _basic_def[_qp](1, 0);
381 9558 : _Fb[_qp](2, 0) = _qTrial(1) + N * _local_def[_qp](4, 0) + k2z * _basic_def[_qp](2, 0);
382 :
383 : // Set tangent stiffnesses
384 9558 : _Kb[_qp](1, 1) = _Kb[_qp](2, 2) = _k0;
385 9558 : _Kb[_qp](1, 2) = _Kb[_qp](2, 1) = 0.0;
386 :
387 : // Account for uplift
388 9558 : if (is_uplift)
389 : {
390 30 : _Kb[_qp](1, 1) = kFactUplift * _k0;
391 30 : _Kb[_qp](2, 2) = kFactUplift * _k0;
392 : }
393 :
394 9558 : _ubPlastic[_qp](0) = _ubPlastic_old[_qp](0);
395 9558 : _ubPlastic[_qp](1) = _ubPlastic_old[_qp](1);
396 : }
397 :
398 : // Plastic step -> return mapping
399 : else
400 : {
401 : // Compute consistency parameters
402 7784 : Real dGammaY = Y / k0y;
403 7784 : Real dGammaZ = Y / k0z;
404 :
405 : // Update plastic displacements
406 7784 : _ubPlastic[_qp](0) = _ubPlastic_old[_qp](0) + dGammaY * _qTrial(0) / qTrialNorm;
407 7784 : _ubPlastic[_qp](1) = _ubPlastic_old[_qp](1) + dGammaZ * _qTrial(1) / qTrialNorm;
408 :
409 : // Set shear forces
410 7784 : _Fb[_qp](1, 0) = qYield * _qTrial(0) / qTrialNorm - N * _local_def[_qp](5, 0) +
411 7784 : k2y * _basic_def[_qp](1, 0);
412 7784 : _Fb[_qp](2, 0) = qYield * _qTrial(1) / qTrialNorm + N * _local_def[_qp](4, 0) +
413 7784 : k2z * _basic_def[_qp](2, 0);
414 :
415 : // Account for uplift
416 7784 : Real k0Plastic = _k0;
417 :
418 7784 : if (is_uplift)
419 0 : k0Plastic = kFactUplift * _k0;
420 :
421 : // Set tangent stiffness
422 : Real D = pow(qTrialNorm, 3);
423 :
424 7784 : _Kb[_qp](1, 1) = k0Plastic * (1.0 - (qYield * _qTrial(1) * _qTrial(1) / D)) + k2y;
425 7784 : _Kb[_qp](1, 2) = -qYield * k0Plastic * _qTrial(0) * _qTrial(1) / D;
426 7784 : _Kb[_qp](2, 1) = _Kb[_qp](1, 2);
427 7784 : _Kb[_qp](2, 2) = k0Plastic * (1.0 - (qYield * _qTrial(0) * _qTrial(0) / D)) + k2z;
428 : }
429 :
430 17342 : iter++;
431 :
432 34684 : } while (sqrt(pow(_Fb[_qp](1, 0) - _qbOld(0), 2) + pow(_Fb[_qp](2, 0) - _qbOld(1), 2) >= _tol) &&
433 5026 : (iter <= _maxiter));
434 :
435 : // Issue warning if iteration did not converge
436 12316 : if (iter >= _maxiter)
437 : {
438 0 : mooseError("Error: ComputeFPIsolatorElasticity::computeshearforce() - did not find the shear "
439 : "force after ",
440 : iter,
441 : " iterations and norm: ",
442 0 : sqrt(pow(_Fb[_qp](1, 0) - _qbOld(0), 2) + pow(_Fb[_qp](2, 0) - _qbOld(1), 2)),
443 : ".\n");
444 : }
445 12316 : }
446 :
447 : void
448 12316 : ComputeFPIsolatorElasticity::addPDeltaEffects()
449 : {
450 : // Add P-Delta moment stiffness terms
451 12316 : Real Ls = (1 - _sD) * _length[_qp];
452 :
453 12316 : _Kl[_qp](5, 1) -= _Fb[_qp](0, 0);
454 12316 : _Kl[_qp](5, 7) += _Fb[_qp](0, 0);
455 12316 : _Kl[_qp](5, 11) -= _Fb[_qp](0, 0) * Ls;
456 12316 : _Kl[_qp](11, 11) -= _Fb[_qp](0, 0) * Ls;
457 12316 : _Kl[_qp](4, 2) += _Fb[_qp](0, 0);
458 12316 : _Kl[_qp](4, 8) -= _Fb[_qp](0, 0);
459 12316 : _Kl[_qp](4, 10) -= _Fb[_qp](0, 0) * Ls;
460 12316 : _Kl[_qp](10, 10) -= _Fb[_qp](0, 0) * Ls;
461 :
462 : // Add V-Delta torsion stiffness terms
463 12316 : _Kl[_qp](3, 1) += _Fb[_qp](2, 0);
464 12316 : _Kl[_qp](3, 2) -= _Fb[_qp](1, 0);
465 12316 : _Kl[_qp](3, 7) -= _Fb[_qp](2, 0);
466 12316 : _Kl[_qp](3, 8) += _Fb[_qp](1, 0);
467 12316 : _Kl[_qp](3, 10) += _Fb[_qp](1, 0) * Ls;
468 12316 : _Kl[_qp](3, 11) += _Fb[_qp](2, 0) * Ls;
469 12316 : _Kl[_qp](9, 10) -= _Fb[_qp](1, 0) * Ls;
470 12316 : _Kl[_qp](9, 11) -= _Fb[_qp](2, 0) * Ls;
471 12316 : }
472 :
473 : void
474 12316 : ComputeFPIsolatorElasticity::finalize()
475 : {
476 : // Convert forces from basic to local to global coordinate system
477 12316 : _Fl[_qp] = _total_lb[_qp].transpose() * _Fb[_qp]; // local forces
478 12316 : _Fg[_qp] = _total_gl[_qp].transpose() * _Fl[_qp]; // global forces
479 :
480 : // Convert stiffness matrix from basic to local coordinate system
481 12316 : _Kl[_qp] = _total_lb[_qp].transpose() * _Kb[_qp] * _total_lb[_qp];
482 :
483 : // add P-∆ and V-∆ effects to local stiffness matrix
484 12316 : addPDeltaEffects();
485 :
486 : // Converting stiffness matrix from local to global coordinate system
487 12316 : _Kg[_qp] = _total_gl[_qp].transpose() * _Kl[_qp] * _total_gl[_qp];
488 12316 : }
|
