MAST
orthotropic_element_property_card_3D.cpp
Go to the documentation of this file.
1 /*
2  * MAST: Multidisciplinary-design Adaptation and Sensitivity Toolkit
3  * Copyright (C) 2013-2019 Manav Bhatia
4  *
5  * This library is free software; you can redistribute it and/or
6  * modify it under the terms of the GNU Lesser General Public
7  * License as published by the Free Software Foundation; either
8  * version 2.1 of the License, or (at your option) any later version.
9  *
10  * This library is distributed in the hope that it will be useful,
11  * but WITHOUT ANY WARRANTY; without even the implied warranty of
12  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
13  * Lesser General Public License for more details.
14  *
15  * You should have received a copy of the GNU Lesser General Public
16  * License along with this library; if not, write to the Free Software
17  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
18  */
19 
20 // MAST includes
25 
26 
27 namespace MAST {
28 
29  namespace OrthotropicProperty3D {
30 
31 
33  public MAST::FieldFunction<RealMatrixX> {
34  public:
36  const MAST::CoordinateBase& orient);
37 
38  virtual ~StiffnessMatrix() { }
39 
40  virtual void operator() (const libMesh::Point& p,
41  const Real t,
42  RealMatrixX& m) const;
43 
44 
45  virtual void derivative ( const MAST::FunctionBase& f,
46  const libMesh::Point& p,
47  const Real t,
48  RealMatrixX& m) const;
49 
50  protected:
51 
53 
55  };
56 
57 
58 
59 
60  class InertiaMatrix: public MAST::FieldFunction<RealMatrixX> {
61  public:
63 
64  virtual ~InertiaMatrix() { }
65 
66  virtual void operator() (const libMesh::Point& p,
67  const Real t,
68  RealMatrixX& m) const;
69 
70  virtual void derivative ( const MAST::FunctionBase& f,
71  const libMesh::Point& p,
72  const Real t,
73  RealMatrixX& m) const;
74 
75  protected:
76 
78  };
79 
80 
81 
82 
83  class ThermalExpansionMatrix: public MAST::FieldFunction<RealMatrixX> {
84  public:
86  const MAST::FieldFunction<RealMatrixX>& mat_expansion,
87  const MAST::CoordinateBase& orient);
88 
90 
91  virtual void operator() (const libMesh::Point& p,
92  const Real t,
93  RealMatrixX& m) const;
94 
95  virtual void derivative ( const MAST::FunctionBase& f,
96  const libMesh::Point& p,
97  const Real t,
98  RealMatrixX& m) const;
99 
100  protected:
101 
105  };
106 
107 
108 
109 
111  public MAST::FieldFunction<RealMatrixX> {
112  public:
114  const MAST::CoordinateBase& orient);
115 
116  virtual ~PrestressAMatrix() { }
117 
118  virtual void operator() (const libMesh::Point& p,
119  const Real t,
120  RealMatrixX& m) const;
121 
122 
123  virtual void derivative ( const MAST::FunctionBase& f,
124  const libMesh::Point& p,
125  const Real t,
126  RealMatrixX& m) const;
127 
128  //virtual void convert_to_vector(const RealMatrixX& m, DenseRealVector& v) const;
129 
130  protected:
131 
134  };
135 
136 
137 
139  public MAST::FieldFunction<RealMatrixX> {
140 
141  public:
142 
144  const MAST::CoordinateBase& orient);
145 
146  virtual ~ThermalConductanceMatrix();
147 
148  virtual void operator() (const libMesh::Point& p,
149  const Real t,
150  RealMatrixX& m) const;
151 
152 
153  virtual void derivative ( const MAST::FunctionBase& f,
154  const libMesh::Point& p,
155  const Real t,
156  RealMatrixX& m) const;
157 
158  protected:
159 
161 
163  };
164 
165 
166 
167 
169  public MAST::FieldFunction<RealMatrixX> {
170 
171  public:
172 
174 
175  virtual ~ThermalCapacitanceMatrix();
176 
177  virtual void operator() (const libMesh::Point& p,
178  const Real t,
179  RealMatrixX& m) const;
180 
181 
182  virtual void derivative ( const MAST::FunctionBase& f,
183  const libMesh::Point& p,
184  const Real t,
185  RealMatrixX& m) const;
186 
187  protected:
188 
190  };
191 
192  }
193 
194 }
195 
196 
197 bool
199  return _material->depends_on(f) || // check if the material property depends on the function
200  MAST::ElementPropertyCardBase::depends_on(f); // check with this property card
201 }
202 
203 
204 
207  const MAST::CoordinateBase& orient):
208 MAST::FieldFunction<RealMatrixX> ("StiffnessMatrix3D"),
210 _orient(orient) {
211 
212  _functions.insert(&mat);
213  _functions.insert(&orient);
214 }
215 
216 
217 void
219 StiffnessMatrix::operator() (const libMesh::Point& p,
220  const Real t,
221  RealMatrixX& m) const {
223  A = RealMatrixX::Zero(3, 3),
224  Tinv = RealMatrixX::Zero(6, 6);
225 
226  _material_stiffness (p, t, m);
227  _orient (p, t, A);
228  _orient.stress_strain_transformation_matrix(A.transpose(), Tinv);
229 
230  // vk' = vj ej.ek' = vj Ajk = A^T v
231  // v = A vk'
232  // sij' = skl ek.ei' el.ej' = skl Aki Alj = A^T s A
233  // s' = T s
234  // s = Tinv s'
235  // s' = C' e'
236  // T s = C' Rinv T R s
237  // s = Tinv C Rinv T R e
238  // C = Tinv C Rinv T R
239  // T R scales last three columns by 1/2
240  // Rinv T scales last three rows by 2
241  // therefore, Rinv T R scales top right 3x3 block by 1/2,
242  // and bottom left 3x3 block by 2.
243  // Also, Rinv T R = T^{-T}
244  m = Tinv * m * Tinv.transpose();
245 }
246 
247 
248 
249 void
252  const libMesh::Point& p,
253  const Real t,
254  RealMatrixX& m) const {
255 
257  dm = RealMatrixX::Zero(6, 6),
258  A = RealMatrixX::Zero(3, 3),
259  dA = RealMatrixX::Zero(3, 3),
260  Tinv = RealMatrixX::Zero(6, 6),
261  dTinv= RealMatrixX::Zero(6, 6);
262 
263  _material_stiffness (p, t, m);
264  _material_stiffness.derivative( f, p, t, dm);
265  _orient (p, t, A);
266  _orient.stress_strain_transformation_matrix(A.transpose(), Tinv);
267  _orient.derivative (f, p, t, dA);
269  dA.transpose(),
270  dTinv);
271 
272  m =
273  dTinv * m * Tinv.transpose() +
274  Tinv * dm * Tinv.transpose() +
275  Tinv * m * dTinv.transpose();
276 }
277 
278 
279 
280 
283 MAST::FieldFunction<RealMatrixX>("InertiaMatrix3D"),
284 _material_inertia(mat) {
285  _functions.insert(&mat);
286 }
287 
288 
289 
290 void
292 InertiaMatrix::operator() (const libMesh::Point& p,
293  const Real t,
294  RealMatrixX& m) const {
295  // this only returns the material inertia
296  _material_inertia(p, t, m);
297 }
298 
299 
300 
301 
302 void
305  const libMesh::Point& p,
306  const Real t,
307  RealMatrixX& m) const {
308 
309  _material_inertia.derivative( f, p, t, m);
310 }
311 
312 
313 
314 
317  const MAST::FieldFunction<RealMatrixX>& mat_expansion,
318  const MAST::CoordinateBase& orient):
319 MAST::FieldFunction<RealMatrixX>("ThermalExpansionMatrix3D"),
320 _material_stiffness(mat_stiff),
321 _material_expansion(mat_expansion),
322 _orient(orient) {
323 
324  _functions.insert(&mat_stiff);
325  _functions.insert(&mat_expansion);
326  _functions.insert(&orient);
327 }
328 
329 
330 
331 
332 void
334 ThermalExpansionMatrix::operator() (const libMesh::Point& p,
335  const Real t,
336  RealMatrixX& m) const {
338  mat = RealMatrixX::Zero(6,1),
339  A = RealMatrixX::Zero(3, 3),
340  Tinv = RealMatrixX::Zero(6, 6);
341 
342  _material_stiffness (p, t, m);
343  _material_expansion (p, t, mat);
344  _orient (p, t, A);
345  _orient.stress_strain_transformation_matrix(A.transpose(), Tinv);
346 
347  // epsilon' = T^{-T} epsilon
348  // epsilon = T^T epsilon'
349  // C = Tinv C' T^{-T}
350  // epsilon_delta_T = C epsilon
351  // = Tinv C' T^{-T} T^T epsilon'
352  // = Tinv C' alpha' delta_temp
353  // hence,
354  // mat = Tinv C' alpha'
355  m = Tinv * m * mat;
356 }
357 
358 
359 
360 
361 
362 
363 void
366  const libMesh::Point& p,
367  const Real t,
368  RealMatrixX& m) const {
369 
371  dm = RealMatrixX::Zero(6, 6),
372  mat = RealMatrixX::Zero(6, 1),
373  dmat = RealMatrixX::Zero(6, 1),
374  A = RealMatrixX::Zero(3, 3),
375  dA = RealMatrixX::Zero(3, 3),
376  Tinv = RealMatrixX::Zero(6, 6),
377  dTinv= RealMatrixX::Zero(6, 6);
378 
379  _material_stiffness (p, t, m);
380  _material_stiffness.derivative( f, p, t, dm);
381  _material_expansion (p, t, mat);
382  _material_expansion.derivative( f, p, t, dmat);
383  _orient (p, t, A);
384  _orient.stress_strain_transformation_matrix(A.transpose(), Tinv);
385  _orient.derivative (f, p, t, dA);
387  dA.transpose(),
388  dTinv);
389 
390  m =
391  dTinv * m * mat +
392  Tinv * dm * mat +
393  Tinv * m * dmat;
394 }
395 
396 
397 
398 
401  const MAST::CoordinateBase& orient):
402 MAST::FieldFunction<RealMatrixX>("PrestressAMatrix3D"),
403 _prestress(prestress),
404 _orient(orient) {
405  _functions.insert(&prestress);
406  _functions.insert(&orient);
407 }
408 
409 
410 
411 
412 void
414 PrestressAMatrix::operator() (const libMesh::Point& p,
415  const Real t,
416  RealMatrixX& m) const {
417  _prestress (p, t, m);
418  libmesh_error();
419 }
420 
421 
422 
423 
424 
425 
426 void
429  const libMesh::Point& p,
430  const Real t,
431  RealMatrixX& m) const {
432  _prestress.derivative( f, p, t, m);
433 }
434 
435 
436 
439  const MAST::CoordinateBase& orient):
440 MAST::FieldFunction<RealMatrixX>("ThermalConductanceMatrix"),
441 _mat_cond(mat_cond),
442 _orient(orient) {
443 
444  _functions.insert(&mat_cond);
445  _functions.insert(&orient);
446 }
447 
448 
449 
452 
453 
454 void
456 operator() (const libMesh::Point& p,
457  const Real t,
458  RealMatrixX& m) const {
459 
461  A = RealMatrixX::Zero(3, 3);
462 
463  _mat_cond (p, t, m);
464  _orient (p, t, A);
465 
466  m = A.transpose() * m * A;
467 }
468 
469 
470 
471 
472 void
474  const libMesh::Point& p,
475  const Real t,
476  RealMatrixX& m) const {
478  dm = RealMatrixX::Zero(3, 3),
479  A = RealMatrixX::Zero(3, 3),
480  dA = RealMatrixX::Zero(3, 3);
481 
482  _mat_cond (p, t, m);
483  _orient (p, t, A);
484  _mat_cond.derivative( f, p, t, dm);
485  _orient.derivative ( f, p, t, dA);
486 
487  m =
488  dA.transpose() * m * A +
489  A.transpose() * dm * A +
490  A.transpose() * m * dA;
491 }
492 
493 
494 
495 
496 
497 
500 MAST::FieldFunction<RealMatrixX>("ThermalCapacitanceMatrix"),
501 _mat_cap(mat_cap) {
502 
503  _functions.insert(&mat_cap);
504 }
505 
506 
507 
510 
511 
512 void
514 operator() (const libMesh::Point& p,
515  const Real t,
516  RealMatrixX& m) const {
517 
518  _mat_cap(p, t, m);
519 }
520 
521 
522 
523 void
525  const libMesh::Point& p,
526  const Real t,
527  RealMatrixX& m) const {
528  _mat_cap.derivative( f, p, t, m);
529 }
530 
531 
532 
535 
536 
537 
538 void
541 
542  _orient = &orient;
543 }
544 
545 
546 std::unique_ptr<MAST::FieldFunction<RealMatrixX> >
549 
552  (_material->stiffness_matrix(3), *_orient);
553 
554  return std::unique_ptr<MAST::FieldFunction<RealMatrixX> >(rval);
555 }
556 
557 
558 
559 std::unique_ptr<MAST::FieldFunction<RealMatrixX> >
562 
563  libmesh_assert(false);
564 
565  return std::unique_ptr<MAST::FieldFunction<RealMatrixX> >(nullptr);
566 }
567 
568 
569 
570 
571 std::unique_ptr<MAST::FieldFunction<RealMatrixX> >
574 
575  libmesh_assert(false);
576 
577  return std::unique_ptr<MAST::FieldFunction<RealMatrixX> >(nullptr);
578 }
579 
580 
581 
582 std::unique_ptr<MAST::FieldFunction<RealMatrixX> >
585 
586  libmesh_assert(false);
587 
588  return std::unique_ptr<MAST::FieldFunction<RealMatrixX> >(nullptr);
589 }
590 
591 
592 std::unique_ptr<MAST::FieldFunction<RealMatrixX> >
595 
598  (_material->inertia_matrix(3));
599 
600  return std::unique_ptr<MAST::FieldFunction<RealMatrixX> >(rval);
601 }
602 
603 
604 
605 std::unique_ptr<MAST::FieldFunction<RealMatrixX> >
608 
611  (_material->stiffness_matrix(3),
612  _material->thermal_expansion_matrix(3),
613  *_orient);
614 
615  return std::unique_ptr<MAST::FieldFunction<RealMatrixX> >(rval);
616 }
617 
618 
619 std::unique_ptr<MAST::FieldFunction<RealMatrixX> >
622 
623  // for 3D elements, there is no difference between the A and B matrices
624  return this->thermal_expansion_A_matrix(e);
625 }
626 
627 
628 
629 std::unique_ptr<MAST::FieldFunction<RealMatrixX> >
632 
633  libmesh_assert(false);
634 
635  return std::unique_ptr<MAST::FieldFunction<RealMatrixX> >(nullptr);
636 }
637 
638 
639 std::unique_ptr<MAST::FieldFunction<RealMatrixX> >
642 
643  libmesh_assert(false);
644 
645  return std::unique_ptr<MAST::FieldFunction<RealMatrixX> >(nullptr);
646 }
647 
648 
649 std::unique_ptr<MAST::FieldFunction<RealMatrixX> >
652 
653  libmesh_assert(false);
654 
655  return std::unique_ptr<MAST::FieldFunction<RealMatrixX> >(nullptr);
656 }
657 
658 
659 std::unique_ptr<MAST::FieldFunction<RealMatrixX> >
662 
665  (_material->conductance_matrix(3), *_orient);
666 
667  return std::unique_ptr<MAST::FieldFunction<RealMatrixX> >(rval);
668 }
669 
670 
671 std::unique_ptr<MAST::FieldFunction<RealMatrixX> >
674 
677  (_material->capacitance_matrix(3));
678 
679  return std::unique_ptr<MAST::FieldFunction<RealMatrixX> >(rval);
680 }
681 
682 
683 
const MAST::FieldFunction< RealMatrixX > & _material_stiffness
virtual void operator()(const libMesh::Point &p, const Real t, RealMatrixX &m) const
calculates the value of the function at the specified point, p, and time, t, and returns it in v...
virtual std::unique_ptr< MAST::FieldFunction< RealMatrixX > > inertia_matrix(const MAST::ElementBase &e) const
virtual std::unique_ptr< MAST::FieldFunction< RealMatrixX > > thermal_conductance_matrix(const MAST::ElementBase &e) const
ThermalConductanceMatrix(const MAST::FieldFunction< RealMatrixX > &mat_cond, const MAST::CoordinateBase &orient)
virtual void derivative(const MAST::FunctionBase &f, const libMesh::Point &p, const Real t, RealMatrixX &m) const
calculates the value of the derivative of function with respect to the function f at the specified po...
virtual void derivative(const MAST::FunctionBase &f, const libMesh::Point &p, const Real t, RealMatrixX &m) const
calculates the value of the derivative of function with respect to the function f at the specified po...
virtual bool depends_on(const MAST::FunctionBase &f) const
returns true if the property card depends on the function f
virtual void derivative(const MAST::FunctionBase &f, const libMesh::Point &p, const Real t, RealMatrixX &m) const
calculates the value of the derivative of function with respect to the function f at the specified po...
PrestressAMatrix(const MAST::FieldFunction< RealMatrixX > &prestress, const MAST::CoordinateBase &orient)
ThermalExpansionMatrix(const MAST::FieldFunction< RealMatrixX > &mat_stiff, const MAST::FieldFunction< RealMatrixX > &mat_expansion, const MAST::CoordinateBase &orient)
virtual std::unique_ptr< MAST::FieldFunction< RealMatrixX > > thermal_expansion_B_matrix(const MAST::ElementBase &e) const
std::set< const MAST::FunctionBase * > _functions
set of functions that this function depends on
libMesh::Real Real
virtual void operator()(const libMesh::Point &p, const Real t, RealMatrixX &m) const
calculates the value of the function at the specified point, p, and time, t, and returns it in v...
const MAST::FieldFunction< RealMatrixX > & _material_inertia
virtual std::unique_ptr< MAST::FieldFunction< RealMatrixX > > thermal_capacitance_matrix(const MAST::ElementBase &e) const
virtual std::unique_ptr< MAST::FieldFunction< RealMatrixX > > stiffness_A_matrix(const MAST::ElementBase &e) const
virtual bool depends_on(const MAST::FunctionBase &f) const
returns true if the property card depends on the function f
virtual void derivative(const MAST::FunctionBase &f, const libMesh::Point &p, const Real t, RealMatrixX &m) const
calculates the value of the derivative of function with respect to the function f at the specified po...
virtual void derivative(const MAST::FunctionBase &f, ValType &v) const
calculates the value of the function derivative and returns it in v.
virtual void operator()(const libMesh::Point &p, const Real t, RealMatrixX &m) const
calculates the value of the function at the specified point, p, and time, t, and returns it in v...
virtual void derivative(const MAST::FunctionBase &f, const libMesh::Point &p, const Real t, RealMatrixX &m) const
calculates the value of the derivative of function with respect to the function f at the specified po...
virtual void derivative(const MAST::FunctionBase &f, const libMesh::Point &p, const Real t, RealMatrixX &m) const
calculates the value of the derivative of function with respect to the function f at the specified po...
Matrix< Real, Dynamic, Dynamic > RealMatrixX
StiffnessMatrix(const MAST::FieldFunction< RealMatrixX > &mat, const MAST::CoordinateBase &orient)
This creates the base class for functions that have a saptial and temporal dependence, and provide sensitivity operations with respect to the functions and parameters.
virtual void operator()(const libMesh::Point &p, const Real t, RealMatrixX &m) const
calculates the value of the function at the specified point, p, and time, t, and returns it in v...
ThermalCapacitanceMatrix(const MAST::FieldFunction< RealMatrixX > &mat_cond)
void stress_strain_transformation_matrix(const RealMatrixX &T, RealMatrixX &mat) const
prepares the matrix mat that transforms stress and strain tensors represented in a 6x1 vector from th...
InertiaMatrix(const MAST::FieldFunction< RealMatrixX > &mat)
virtual std::unique_ptr< MAST::FieldFunction< RealMatrixX > > damping_matrix(const MAST::ElementBase &e) const
void stress_strain_transformation_matrix_sens(const RealMatrixX &T, const RealMatrixX &dT, RealMatrixX &mat) const
virtual std::unique_ptr< MAST::FieldFunction< RealMatrixX > > prestress_A_matrix(MAST::ElementBase &e) const
virtual void operator()(const libMesh::Point &p, const Real t, RealMatrixX &m) const
calculates the value of the function at the specified point, p, and time, t, and returns it in v...
virtual std::unique_ptr< MAST::FieldFunction< RealMatrixX > > prestress_B_matrix(MAST::ElementBase &e) const
virtual std::unique_ptr< MAST::FieldFunction< RealMatrixX > > stiffness_D_matrix(const MAST::ElementBase &e) const
virtual void operator()(const libMesh::Point &p, const Real t, RealMatrixX &m) const
calculates the value of the function at the specified point, p, and time, t, and returns it in v...
virtual void set_orientation(const MAST::CoordinateBase &orient)
sets the orientation coordinate system for this section.
virtual std::unique_ptr< MAST::FieldFunction< RealMatrixX > > stiffness_B_matrix(const MAST::ElementBase &e) const
Provides the transformation matrix T to transform vector from the orientation provided in this matrix...
virtual std::unique_ptr< MAST::FieldFunction< RealMatrixX > > transverse_shear_stiffness_matrix(const MAST::ElementBase &e) const
virtual std::unique_ptr< MAST::FieldFunction< RealMatrixX > > thermal_expansion_A_matrix(const MAST::ElementBase &e) const
This is the base class for elements that implement calculation of finite element quantities over the ...
Definition: elem_base.h:72