MAST
Multidisciplinary-design Adaptation and Sensitivity Toolkit (MAST)
solid_1d_section_element_property_card.cpp
Go to the documentation of this file.
1 /*
2  * MAST: Multidisciplinary-design Adaptation and Sensitivity Toolkit
3  * Copyright (C) 2013-2020 Manav Bhatia and MAST authors
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
24 #include "base/elem_base.h"
25 
26 
27 
28 namespace MAST {
29  namespace Solid1DSectionProperty {
30 
31  class Area: public MAST::FieldFunction<Real> {
32  public:
34  const MAST::FieldFunction<Real>& hz):
35  MAST::FieldFunction<Real>("Area"),
36  _hy(hy),
37  _hz(hz) {
38  _functions.insert(&hy);
39  _functions.insert(&hz);
40  }
41 
42  virtual ~Area() { }
43 
44  virtual void operator() (const libMesh::Point& p,
45  const Real t,
46  Real& m) const {
47  Real hy, hz;
48  _hy(p, t, hy);
49  _hz(p, t, hz);
50 
51  m = hy*hz;
52  }
53 
54  virtual void derivative ( const MAST::FunctionBase& f,
55  const libMesh::Point& p,
56  const Real t,
57  Real& m) const {
58  Real hy, hz, dhy, dhz;
59  _hy(p, t, hy); _hy.derivative( f, p, t, dhy);
60  _hz(p, t, hz); _hz.derivative( f, p, t, dhz);
61 
62  m = dhy*hz + hy*dhz;
63  }
64 
65  protected:
66 
68  };
69 
70 
71 
72  class TorsionalConstant: public MAST::FieldFunction<Real> {
73  public:
75  const MAST::FieldFunction<Real>& hz):
76  MAST::FieldFunction<Real>("TorsionalConstant"),
77  _hy(hy),
78  _hz(hz) {
79  _functions.insert(&hy);
80  _functions.insert(&hz);
81  }
82 
83 
84  virtual ~TorsionalConstant() { }
85 
86  virtual void operator() (const libMesh::Point& p,
87  const Real t,
88  Real& m) const {
89  Real hy, hz, a, b;
90  _hy(p, t, hy);
91  _hz(p, t, hz);
92 
93 
94  // shorter side is b, and longer side is a
95  if (hy > hz) {
96  a = hy;
97  b = hz;
98  }
99  else {
100  a = hz;
101  b = hy;
102  }
103 
104  m = a*pow(b,3)*(1./3.-.21*b/a*(1.-pow(b/a,4)/12.));
105  }
106 
107 
108  virtual void derivative ( const MAST::FunctionBase& f,
109  const libMesh::Point& p,
110  const Real t,
111  Real& m) const {
112  Real hy, hz, dhy, dhz, a, b, da, db;
113  _hy(p, t, hy); _hy.derivative( f, p, t, dhy);
114  _hz(p, t, hz); _hz.derivative( f, p, t, dhz);
115 
116  // shorter side is b, and longer side is a
117  if (hy > hz) {
118  a = hy; da = dhy;
119  b = hz; db = dhz;
120  }
121  else {
122  a = hz; da = dhz;
123  b = hy; db = dhy;
124  }
125 
126  m =
127  da*pow(b,3)*(1./3.-.21*b/a*(1.-pow(b/a,4)/12.)) +
128  a*3.*pow(b,2)*db*(1./3.-.21*b/a*(1.-pow(b/a,4)/12.)) +
129  a*pow(b,3)*(-.21*db/a*(1.-pow(b/a,4)/12.) +
130  (.21*b/pow(a,2)*da*(1.-pow(b/a,4)/12.)) +
131  (-.21*b/a*(-4.*pow(b,3)*db/pow(a,4)/12.+
132  4.*pow(b,4)/pow(a,5)*da/12.)));
133  }
134 
135  protected:
136 
138  };
139 
140 
141 
142  class PolarInertia: public MAST::FieldFunction<Real> {
143  public:
145  const MAST::FieldFunction<Real>& hz,
146  const MAST::FieldFunction<Real>& hy_offset,
147  const MAST::FieldFunction<Real>& hz_offset):
148  MAST::FieldFunction<Real>("PolarInertia"),
149  _hy(hy),
150  _hz(hz),
151  _hy_offset(hy_offset),
152  _hz_offset(hz_offset) {
153  _functions.insert(&hy);
154  _functions.insert(&hz);
155  _functions.insert(&hy_offset);
156  _functions.insert(&hz_offset);
157  }
158 
159  virtual ~PolarInertia() { }
160 
161  virtual void operator() (const libMesh::Point& p,
162  const Real t,
163  Real& m) const {
164  Real hy, hz, offy, offz;
165  _hy(p, t, hy);
166  _hz(p, t, hz);
167  _hy_offset(p, t, offy);
168  _hz_offset(p, t, offz);
169 
170  m = hy*hz*((pow(hy,2) + pow(hz,2))/12. +
171  pow(offy,2) + pow(offz,2));
172  }
173 
174  virtual void derivative ( const MAST::FunctionBase& f,
175  const libMesh::Point& p,
176  const Real t,
177  Real& m) const {
178  Real hy, hz, dhy, dhz, offy, offz, doffy, doffz;
179  _hy (p, t, hy); _hy.derivative( f, p, t, dhy);
180  _hz (p, t, hz); _hz.derivative( f, p, t, dhz);
181  _hy_offset (p, t, offy); _hy_offset.derivative( f, p, t, doffy);
182  _hz_offset (p, t, offz); _hz_offset.derivative( f, p, t, doffz);
183 
184 
185  m =
186  (dhy*hz + hy*dhz) * ((pow(hy,2) + pow(hz,2))/12. +
187  pow(offy,2) + pow(offz,2)) +
188  2.*hy*hz*((hy*dhy + hz*dhz)/12. +
189  offy*doffy + offz*doffz);
190  }
191 
192  protected:
193 
194  const MAST::FieldFunction<Real>& _hy, &_hz, &_hy_offset, &_hz_offset;
195  };
196 
197 
198 
203  class AreaYMoment: public MAST::FieldFunction<Real> {
204  public:
206  const MAST::FieldFunction<Real>& hz,
207  const MAST::FieldFunction<Real>& hz_offset):
208  MAST::FieldFunction<Real>("AreaYMoment"),
209  _hy(hy),
210  _hz(hz),
211  _hz_offset(hz_offset) {
212  _functions.insert(&hy);
213  _functions.insert(&hz);
214  _functions.insert(&hz_offset);
215  }
216 
217 
218  virtual ~AreaYMoment() { }
219 
220  virtual void operator() (const libMesh::Point& p,
221  const Real t,
222  Real& m) const {
223  Real hy, hz, off;
224  _hy(p, t, hy);
225  _hz(p, t, hz);
226  _hz_offset(p, t, off);
227 
228  m = hy*hz*off;
229  }
230 
231 
232  virtual void derivative ( const MAST::FunctionBase& f,
233  const libMesh::Point& p,
234  const Real t,
235  Real& m) const {
236  Real hy, hz, off, dhy, dhz, doff;
237  _hy (p, t, hy); _hy.derivative( f, p, t, dhy);
238  _hz (p, t, hz); _hz.derivative( f, p, t, dhz);
239  _hz_offset (p, t, off); _hz_offset.derivative( f, p, t, doff);
240 
241  m = dhy*hz*off + hy*dhz*off + hy*hz*doff;
242  }
243 
244  protected:
245 
247  };
248 
249 
250 
255  class AreaZMoment: public MAST::FieldFunction<Real> {
256  public:
258  const MAST::FieldFunction<Real>& hz,
259  const MAST::FieldFunction<Real>& hy_offset):
260  MAST::FieldFunction<Real>("AreaZMoment"),
261  _hy(hy),
262  _hz(hz),
263  _hy_offset(hy_offset) {
264  _functions.insert(&hy);
265  _functions.insert(&hz);
266  _functions.insert(&hy_offset);
267  }
268 
269  virtual ~AreaZMoment() { }
270 
271  virtual void operator() (const libMesh::Point& p,
272  const Real t,
273  Real& m) const {
274  Real hy, hz, off;
275  _hy(p, t, hy);
276  _hz(p, t, hz);
277  _hy_offset(p, t, off);
278 
279  m = hy*hz*off;
280  }
281 
282 
283  virtual void derivative ( const MAST::FunctionBase& f,
284  const libMesh::Point& p,
285  const Real t,
286  Real& m) const {
287  Real hy, hz, off, dhy, dhz, doff;
288  _hy(p, t, hy); _hy.derivative( f, p, t, dhy);
289  _hz(p, t, hz); _hz.derivative( f, p, t, dhz);
290  _hy_offset(p, t, off); _hy_offset.derivative( f, p, t, doff);
291 
292  m = dhy*hz*off + hy*dhz*off + hy*hz*doff;
293  }
294 
295  protected:
296 
297  const MAST::FieldFunction<Real>& _hy, &_hz, &_hy_offset;
298  };
299 
300 
301 
311  class AreaInertiaMatrix: public MAST::FieldFunction<RealMatrixX> {
312  public:
314  const MAST::FieldFunction<Real>& hz,
315  const MAST::FieldFunction<Real>& hy_offset,
316  const MAST::FieldFunction<Real>& hz_offset):
317  MAST::FieldFunction<RealMatrixX>("AreaInertiaMatrix"),
318  _hy(hy),
319  _hz(hz),
320  _hy_offset(hy_offset),
321  _hz_offset(hz_offset) {
322  _functions.insert(&hy);
323  _functions.insert(&hz);
324  _functions.insert(&hy_offset);
325  _functions.insert(&hz_offset);
326  }
327 
328  virtual ~AreaInertiaMatrix() { }
329 
330  virtual void operator() (const libMesh::Point& p,
331  const Real t,
332  RealMatrixX& m) const {
333  Real hy, hz, offy, offz;
334  m = RealMatrixX::Zero(2,2);
335  _hy(p, t, hy);
336  _hz(p, t, hz);
337  _hy_offset(p, t, offy);
338  _hz_offset(p, t, offz);
339 
340  m(0,0) = hz*pow(hy,3)/12. + hy*hz*pow(offy,2) ; // Izz for v-bending
341  m(0,1) = hy*hz*offy*offz;
342  m(1,0) = m(0,1);
343  m(1,1) = hy*pow(hz,3)/12. + hy*hz*pow(offz,2) ; // Iyy for w-bending
344  }
345 
346 
347  virtual void derivative ( const MAST::FunctionBase& f,
348  const libMesh::Point& p,
349  const Real t,
350  RealMatrixX& m) const {
351  Real hy, hz, offy, offz, dhy, dhz, doffy, doffz;
352  m = RealMatrixX::Zero(2,2);
353  _hy(p, t, hy); _hy.derivative( f, p, t, dhy);
354  _hz(p, t, hz); _hz.derivative( f, p, t, dhz);
355  _hy_offset(p, t, offy); _hy_offset.derivative( f, p, t, doffy);
356  _hz_offset(p, t, offz); _hz_offset.derivative( f, p, t, doffz);
357 
358 
359  m(0,0) = dhz*pow(hy,3)/12. + hz*pow(hy,2)/4.*dhy +
360  dhy*hz*pow(offy,2) + hy*dhz*pow(offy,2) + 2.*hy*hz*offy*doffy ;
361  m(0,1) = dhy*hz*offy*offz + hy*dhz*offy*offz +
362  hy*hz*doffy*offz + hy*hz*offy*doffz;
363  m(1,0) = m(0,1);
364  m(1,1) = dhy*pow(hz,3)/12. + hy*pow(hz,2)/4.*dhz +
365  dhy*hz*pow(offz,2) + hy*dhz*pow(offz,2) + 2.*hy*hz*offz*doffz ;
366  }
367 
368  protected:
369 
370  const MAST::FieldFunction<Real>& _hy, &_hz, &_hy_offset, &_hz_offset;
371  };
372 
373 
375  {
376  public:
378  MAST::FieldFunction<Real>("WarpingConstant")
379  {
380  }
381 
383  {
384  }
385 
386  virtual void operator() (const libMesh::Point& p,
387  const Real t,
388  Real& m) const
389  {
390  m = 0.0;
391  }
392 
393  virtual void derivative (const MAST::FunctionBase& f,
394  const libMesh::Point& p,
395  const Real t,
396  Real& m) const
397  {
398  m = 0.0;
399  }
400  };
401 
402 
403  class ShearCoefficientMatrix: public MAST::FieldFunction<RealMatrixX> {
404  public:
406  const MAST::FieldFunction<Real>& Kyy):
407  MAST::FieldFunction<RealMatrixX>("ShearCoefficientMatrix"),
408  _Kzz(Kzz),
409  _Kyy(Kyy)
410  {
411  _functions.insert(&Kzz);
412  _functions.insert(&Kyy);
413  }
414 
416 
417  virtual void operator() (const libMesh::Point& p,
418  const Real t,
419  RealMatrixX& m) const {
420  Real Kzz;
421  Real Kyy;
422  _Kzz(p, t, Kzz);
423  _Kyy(p, t, Kyy);
424  m = RealMatrixX::Zero(2,2);
425  m(0,0) = Kzz;
426  m(1,1) = Kyy;
427  }
428 
429  virtual void derivative (const MAST::FunctionBase& f,
430  const libMesh::Point& p,
431  const Real t,
432  RealMatrixX& m) const {
433  Real dKzz, dKyy;
434  _Kyy.derivative(f, p, t, dKyy);
435  _Kzz.derivative(f, p, t, dKzz);
436  m = RealMatrixX::Zero(2,2);
437  m(0,0) = dKzz;
438  m(1,1) = dKyy;
439  }
440  protected:
443  };
444 
445 
446  class ExtensionStiffnessMatrix: public MAST::FieldFunction<RealMatrixX> {
447  public:
449  const MAST::FieldFunction<Real>& A,
450  const MAST::FieldFunction<Real>& J);
451 
453 
454  virtual void operator() (const libMesh::Point& p,
455  const Real t,
456  RealMatrixX& m) const;
457 
458 
459 
460  virtual void derivative ( const MAST::FunctionBase& f,
461  const libMesh::Point& p,
462  const Real t,
463  RealMatrixX& m) const;
464 
465  protected:
466 
469  };
470 
471 
472 
474  public:
476  const MAST::FieldFunction<Real>& A_y_moment,
477  const MAST::FieldFunction<Real>& A_z_moment);
478 
480 
481  virtual void operator() (const libMesh::Point& p,
482  const Real t,
483  RealMatrixX& m) const;
484 
485 
486 
487  virtual void derivative ( const MAST::FunctionBase& f,
488  const libMesh::Point& p,
489  const Real t,
490  RealMatrixX& m) const;
491 
492  protected:
493 
496  };
497 
498 
499  class BendingStiffnessMatrix: public MAST::FieldFunction<RealMatrixX> {
500  public:
503 
505 
506  virtual void operator() (const libMesh::Point& p,
507  const Real t,
508  RealMatrixX& m) const;
509 
510 
511 
512  virtual void derivative ( const MAST::FunctionBase& f,
513  const libMesh::Point& p,
514  const Real t,
515  RealMatrixX& m) const;
516 
517  protected:
518 
521  };
522 
523 
524  class TransverseStiffnessMatrix: public MAST::FieldFunction<RealMatrixX> {
525  public:
527  const MAST::FieldFunction<Real>& A,
528  const MAST::FieldFunction<RealMatrixX>& Kappa):
529  MAST::FieldFunction<RealMatrixX>("TransverseStiffnessMatrix1D"),
530  _material_stiffness(mat),
531  _A(A),
532  _Kappa(Kappa)
533  {
534  _functions.insert(&mat);
535  _functions.insert(&A);
536  _functions.insert(&Kappa);
537  }
538 
540 
541  virtual void operator() (const libMesh::Point& p,
542  const Real t,
543  RealMatrixX& m) const {
544  Real A;
545  RealMatrixX Kappa;
546  _A(p, t, A);
547  _Kappa(p, t, Kappa);
548  _material_stiffness(p, t, m);
549 
550  m *= A;
551 
552  m(0,0) *= Kappa(0,0);
553  m(1,1) *= Kappa(1,1);
554  m(0,1) *= Kappa(0,1);
555  m(1,0) *= Kappa(1,0);
556 
561  }
562 
563  virtual void derivative ( const MAST::FunctionBase& f,
564  const libMesh::Point& p,
565  const Real t,
566  RealMatrixX& m) const {
567  RealMatrixX G, dG, dm, Kappa, dKappa;
568  Real A, dA;
569  _A (p, t, A); _A.derivative( f, p, t, dA);
570  _material_stiffness (p, t, G); _material_stiffness.derivative( f, p, t, dG);
571  _Kappa (p, t, Kappa);
572  _Kappa.derivative(f, p, t, dKappa);
573 
574  m(0,0) = (dG(0,0) * Kappa(0,0) * A) + (G(0,0) * dKappa(0,0) * A) + (G(0,0) * Kappa(0,0) * dA);
575  m(1,1) = (dG(1,1) * Kappa(1,1) * A) + (G(1,1) * dKappa(1,1) * A) + (G(1,1) * Kappa(1,1) * dA);
576  m(0,1) = (dG(0,1) * Kappa(0,1) * A) + (G(0,1) * dKappa(0,1) * A) + (G(0,1) * Kappa(0,1) * dA);
577  m(1,0) = (dG(1,0) * Kappa(1,0) * A) + (G(1,0) * dKappa(1,0) * A) + (G(1,0) * Kappa(1,0) * dA);
578  }
579 
580  protected:
581 
585  };
586 
587 
588  class InertiaMatrix: public MAST::FieldFunction<RealMatrixX> {
589  public:
591  const MAST::FieldFunction<Real>& A,
592  const MAST::FieldFunction<Real>& A_y_moment,
593  const MAST::FieldFunction<Real>& A_z_moment,
594  const MAST::FieldFunction<Real>& Ip,
596 
597  virtual ~InertiaMatrix() { }
598 
599  virtual void operator() (const libMesh::Point& p,
600  const Real t,
601  RealMatrixX& m) const;
602 
603 
604 
605  virtual void derivative ( const MAST::FunctionBase& f,
606  const libMesh::Point& p,
607  const Real t,
608  RealMatrixX& m) const;
609 
610  protected:
611 
612  const MAST::FieldFunction<Real>& _rho, &_A, &_A_y_moment, &_A_z_moment, &_Ip;
614  };
615 
616 
617 
618  class ThermalExpansionAMatrix: public MAST::FieldFunction<RealMatrixX> {
619  public:
621  const MAST::FieldFunction<RealMatrixX>& mat_expansion,
622  const MAST::FieldFunction<Real>& A);
623 
625 
626  virtual void operator() (const libMesh::Point& p,
627  const Real t,
628  RealMatrixX& m) const;
629 
630 
631 
632  virtual void derivative ( const MAST::FunctionBase& f,
633  const libMesh::Point& p,
634  const Real t,
635  RealMatrixX& m) const;
636 
637  protected:
638 
642  };
643 
644 
645 
646  class ThermalExpansionBMatrix: public MAST::FieldFunction<RealMatrixX> {
647  public:
649  const MAST::FieldFunction<RealMatrixX>& mat_expansion,
650  const MAST::FieldFunction<Real>& A_y_moment,
651  const MAST::FieldFunction<Real>& A_z_moment);
652 
654 
655  virtual void operator() (const libMesh::Point& p,
656  const Real t,
657  RealMatrixX& m) const;
658 
659 
660 
661  virtual void derivative ( const MAST::FunctionBase& f,
662  const libMesh::Point& p,
663  const Real t,
664  RealMatrixX& m) const;
665 
666  protected:
667 
671  };
672 
673 
674 
675 
677  public MAST::FieldFunction<RealMatrixX> {
678  public:
681  const MAST::FieldFunction<Real>& A);
682 
683  virtual ~PrestressAMatrix() { }
684 
685  virtual void operator() (const libMesh::Point& p,
686  const Real t,
687  RealMatrixX& m) const;
688 
689 
690 
691  virtual void derivative ( const MAST::FunctionBase& f,
692  const libMesh::Point& p,
693  const Real t,
694  RealMatrixX& m) const;
695 
696  //virtual void convert_to_vector(const RealMatrixX& m, DenseRealVector& v) const;
697 
698  protected:
699 
702  };
703 
704 
705 
707  public MAST::FieldFunction<RealMatrixX> {
708  public:
711  const MAST::FieldFunction<Real>& A_y_moment,
712  const MAST::FieldFunction<Real>& A_z_moment);
713 
714  virtual ~PrestressBMatrix() { }
715 
716  virtual void operator() (const libMesh::Point& p,
717  const Real t,
718  RealMatrixX& m) const;
719 
720 
721 
722  virtual void derivative ( const MAST::FunctionBase& f,
723  const libMesh::Point& p,
724  const Real t,
725  RealMatrixX& m) const;
726 
727  //virtual void convert_to_vector(const RealMatrixX& m, DenseRealVector& v) const;
728 
729  protected:
730 
733  };
734 
735 
736 
737 
739  public MAST::FieldFunction<RealMatrixX> {
740 
741  public:
742 
744  const MAST::FieldFunction<Real>& A);
745 
746  virtual ~ThermalConductanceMatrix();
747 
748  virtual void operator() (const libMesh::Point& p,
749  const Real t,
750  RealMatrixX& m) const;
751 
752 
753 
754  virtual void derivative ( const MAST::FunctionBase& f,
755  const libMesh::Point& p,
756  const Real t,
757  RealMatrixX& m) const;
758 
759  protected:
760 
762 
764  };
765 
766 
767 
768 
770  public MAST::FieldFunction<RealMatrixX> {
771 
772  public:
773 
775  const MAST::FieldFunction<Real>& h);
776 
777  virtual ~ThermalCapacitanceMatrix();
778 
779  virtual void operator() (const libMesh::Point& p,
780  const Real t,
781  RealMatrixX& m) const;
782 
783 
784 
785  virtual void derivative ( const MAST::FunctionBase& f,
786  const libMesh::Point& p,
787  const Real t,
788  RealMatrixX& m) const;
789 
790  protected:
791 
793 
795  };
796 
797 
798 
799  }
800 
801 
802 }
803 
804 
808  const MAST::FieldFunction<Real>& A,
809  const MAST::FieldFunction<Real>& J):
810 MAST::FieldFunction<RealMatrixX> ("ExtensionStiffnessMatrix1D"),
811 _material_stiffness(mat),
812 _A(A),
813 _J(J) {
814  _functions.insert(&mat);
815  _functions.insert(&A);
816  _functions.insert(&J);
817 }
818 
819 
820 
821 
822 
823 void
826  const Real t,
827  RealMatrixX& m) const {
828  // [C]*h
829  Real A, J;
830  _A(p, t, A);
831  _J(p, t, J);
832  _material_stiffness(p, t, m);
833  m.row(0) *= A;
834  m.row(1) *= J;
835 }
836 
837 
838 
839 
840 void
843  const libMesh::Point& p,
844  const Real t,
845  RealMatrixX& m) const {
846  RealMatrixX dm;
847  m = RealMatrixX::Zero(2,2);
848  Real A, J, dA, dJ;
849  _A(p, t, A); _A.derivative( f, p, t, dA);
850  _J(p, t, J); _J.derivative( f, p, t, dJ);
851  _material_stiffness(p, t, m); _material_stiffness.derivative( f, p, t, dm);
852 
853  // [C]*dh
854  m.row(0) *= dA;
855  m.row(1) *= dJ;
856 
857  // += [dC]*h
858  dm.row(0) *= A;
859  dm.row(1) *= J;
860  m += dm;
861 }
862 
863 
864 
865 
866 
867 
870  const MAST::FieldFunction<Real>& A_y_moment,
871  const MAST::FieldFunction<Real>& A_z_moment):
872 MAST::FieldFunction<RealMatrixX> ("ExtensionBendingStiffnessMatrix1D"),
874 _A_y_moment(A_y_moment),
875 _A_z_moment(A_z_moment) {
876  _functions.insert(&mat);
877  _functions.insert(&A_y_moment);
878  _functions.insert(&A_z_moment);
879 }
880 
881 
882 
883 void
886  const Real t,
887  RealMatrixX& m) const {
888  Real Ay, Az;
889  _A_y_moment(p, t, Ay);
890  _A_z_moment(p, t, Az);
891  _material_stiffness(p, t, m);
892 
893  m(0,1) = m(0,0)*Ay; // coupling of u and w bending (== theta_y)
894  m(0,0) *= Az; // coupling of u and v bending (== theta_z)
895 
896  m.row(1) *= 0; // no coupling for torsion for symmetic sections
897 }
898 
899 
900 
901 
902 
903 void
906  const libMesh::Point& p,
907  const Real t,
908  RealMatrixX& m) const {
909  RealMatrixX dm;
910  Real Ay, Az, dAy, dAz;
911  _A_y_moment(p, t, Ay); _A_y_moment.derivative( f, p, t, dAy);
912  _A_z_moment(p, t, Az); _A_z_moment.derivative( f, p, t, dAz);
913  _material_stiffness(p, t, m); _material_stiffness.derivative( f, p, t, dm);
914 
915  m(0,1) = m(0,0)*dAy; // coupling of u and w bending (== theta_y)
916  m(0,0) *= dAz; // coupling of u and v bending (== theta_z)
917  m.row(1) *= 0; // no coupling for torsion for symmetic sections
918 
919  dm(0,1) = dm(0,0)*Ay;
920  dm(0,0) *= Az;
921  dm.row(1)*= 0.;
922  m += dm;
923 }
924 
925 
926 
927 
931 MAST::FieldFunction<RealMatrixX> ("BendingStiffnessMatrix1D"),
933 _I(I) {
934  _functions.insert(&mat);
935  _functions.insert(&I);
936 }
937 
938 
939 
940 void
942 BendingStiffnessMatrix::operator() (const libMesh::Point& p,
943  const Real t,
944  RealMatrixX& m) const {
945  RealMatrixX mat;
946  _I(p, t, m);
947  _material_stiffness(p, t, mat);
948 
949  // E*I
950  m *= mat(0,0); // scale the inertia matrix with modulus of elasticity
951 }
952 
953 
954 
955 
956 
957 
958 
959 
960 void
963  const libMesh::Point& p,
964  const Real t,
965  RealMatrixX& m) const {
966  RealMatrixX mat, dmat, dm;
967  _I(p, t, m); _I.derivative( f, p, t, dm);
968  _material_stiffness(p, t, mat); _material_stiffness.derivative( f, p, t, dmat);
969 
970  // dE*I
971  m *= dmat(0,0); // scale the inertia matrix with modulus of elasticity
972 
973  // E*dI
974  m += mat(0,0)*dm; // scale the inertia matrix with modulus of elasticity
975 }
976 
977 
978 
979 
980 
983  const MAST::FieldFunction<Real>& A,
984  const MAST::FieldFunction<Real>& A_y_moment,
985  const MAST::FieldFunction<Real>& A_z_moment,
986  const MAST::FieldFunction<Real>& Ip,
988 MAST::FieldFunction<RealMatrixX>("InertiaMatrix1D"),
989 _rho(rho),
990 _A(A),
991 _A_y_moment(A_y_moment),
992 _A_z_moment(A_z_moment),
993 _Ip(Ip),
994 _I(I) {
995  _functions.insert(&_rho);
996  _functions.insert(&_A);
997  _functions.insert(&_A_y_moment);
998  _functions.insert(&_A_z_moment);
999  _functions.insert(&_Ip);
1000  _functions.insert(&_I);
1001 }
1002 
1003 
1004 
1005 
1006 void
1008 InertiaMatrix::operator() (const libMesh::Point& p,
1009  const Real t,
1010  RealMatrixX& m) const {
1011  m = RealMatrixX::Zero(6, 6);
1012  RealMatrixX I;
1013  Real rho, A, Ay, Az, Ip;
1014  _rho(p, t, rho);
1015  _A(p, t, A);
1016  _A_y_moment(p, t, Ay);
1017  _A_z_moment(p, t, Az);
1018  _Ip(p, t, Ip);
1019  _I(p, t, I);
1020 
1021  // translation velocities
1022  m(0,0) = A; m(1,1) = A; m(2,2) = A;
1023 
1024  // torsion
1025  m(3,3) = Ip;
1026 
1027  // rotational velocities
1028  // theta-y rotation
1029  m(0,4) = Ay; m(4,0) = Ay;
1030  m(4,4) = I(1,1); // I11 is defined about y-y axis for theta-y
1031 
1032  // theta-z rotation
1033  m(0,5) = Az; m(5,0) = Az;
1034  m(5,5) = I(0,0); // I00 is defined about z-z axis for theta-z
1035 
1036  m(4,5) = m(5,4) = I(0,1);
1037 
1038  m *= rho;
1039 }
1040 
1041 
1042 
1043 
1044 
1045 
1046 
1047 
1048 void
1051  const libMesh::Point& p,
1052  const Real t,
1053  RealMatrixX& m) const {
1054  RealMatrixX dm;
1055  m = RealMatrixX::Zero(6, 6); dm = RealMatrixX::Zero(6, 6);
1056  RealMatrixX I, dI;
1057  Real rho, A, Ay, Az, Ip, drho, dA, dAy, dAz, dIp;
1058  _rho(p, t, rho); _rho.derivative( f, p, t, drho);
1059  _A(p, t, A); _A.derivative( f, p, t, dA);
1060  _A_y_moment(p, t, Ay); _A_y_moment.derivative( f, p, t, dAy);
1061  _A_z_moment(p, t, Az); _A_z_moment.derivative( f, p, t, dAz);
1062  _Ip(p, t, Ip); _Ip.derivative( f, p, t, dIp);
1063  _I(p, t, I); _I.derivative( f, p, t, dI);
1064 
1065  // translation velocities
1066  m(0,0) = A; m(1,1) = A; m(2,2) = A;
1067  dm(0,0) = dA; dm(1,1) = dA; dm(2,2) = dA;
1068 
1069  // torsion
1070  m(3,3) = Ip;
1071  dm(3,3) = dIp;
1072 
1073  // rotational velocities
1074  // theta-y rotation
1075  m(0,4) = Ay; m(4,0) = Ay;
1076  m(4,4) = I(1,1);
1077  dm(0,4) = dAy; dm(4,0) = dAy;
1078  dm(4,4) = dI(1,1);
1079 
1080  // theta-z rotation
1081  m(0,5) = Az; m(5,0) = Az;
1082  m(5,5) = I(0,0);
1083  dm(0,5) = dAz; dm(5,0) = dAz; // v-displacement
1084  dm(5,5) = dI(0,0);
1085 
1086  m(4,5) = m(5,4) = I(0,1);
1087  dm(4,5) = dm(5,4) = dI(0,1);
1088 
1089  m *= drho;
1090  m += rho*dm;
1091 }
1092 
1093 
1094 
1095 
1098  const MAST::FieldFunction<RealMatrixX>& mat_expansion,
1099  const MAST::FieldFunction<Real>& A):
1100 MAST::FieldFunction<RealMatrixX>("ThermalExpansionAMatrix1D"),
1101 _material_stiffness(mat_stiff),
1102 _material_expansion(mat_expansion),
1103 _A(A) {
1104  _functions.insert(&mat_stiff);
1105  _functions.insert(&mat_expansion);
1106  _functions.insert(&_A);
1107 }
1108 
1109 
1110 
1111 
1112 void
1115  const Real t,
1116  RealMatrixX& m) const {
1117  Real A;
1118  RealMatrixX at;
1119  _A(p, t, A);
1120  _material_stiffness(p, t, m);
1121  _material_expansion(p, t, at);
1122 
1123  m *= at;
1124  m *= A;
1125 }
1126 
1127 
1128 
1129 
1130 
1131 void
1134  const libMesh::Point& p,
1135  const Real t,
1136  RealMatrixX& m) const {
1137  Real A, dA;
1138  RealMatrixX m1, at, dat, dm;
1139  _A(p, t, A); _A.derivative( f, p, t, dA);
1140  _material_stiffness(p, t, m1); _material_stiffness.derivative( f, p, t, dm);
1141  _material_expansion(p, t, at); _material_expansion.derivative( f, p, t, dat);
1142 
1143  m=m1;
1144 
1145  m *= at;
1146  m *= dA;
1147 
1148  m1 *= dat;
1149  dm *= at;
1150  m1 += dm;
1151 
1152  m += A*m1;
1153 }
1154 
1155 
1156 
1157 
1160  const MAST::FieldFunction<RealMatrixX>& mat_expansion,
1161  const MAST::FieldFunction<Real>& A_y_moment,
1162  const MAST::FieldFunction<Real>& A_z_moment):
1163 MAST::FieldFunction<RealMatrixX>("ThermalExpansionBMatrix1D"),
1164 _material_stiffness(mat_stiff),
1165 _material_expansion(mat_expansion),
1166 _A_y_moment(A_y_moment),
1167 _A_z_moment(A_z_moment) {
1168  _functions.insert(&mat_stiff);
1169  _functions.insert(&mat_expansion);
1170  _functions.insert(&_A_y_moment);
1171  _functions.insert(&_A_z_moment);
1172 }
1173 
1174 
1175 
1176 
1177 void
1180  const Real t,
1181  RealMatrixX& m) const {
1182  Real Ay, Az;
1183  RealMatrixX at;
1184  _A_y_moment(p, t, Ay);
1185  _A_z_moment(p, t, Az);
1186  _material_stiffness(p, t, m);
1187  _material_expansion(p, t, at);
1188 
1189  m *= at;
1190  m(1,0) = Ay * m(0,0); // for w-displacement, area moment about Y-axis
1191  m(0,0) *= Az; // for v-displacement, area moment about Z-axis
1192 }
1193 
1194 
1195 
1196 
1197 
1198 
1199 void
1202  const libMesh::Point& p,
1203  const Real t,
1204  RealMatrixX& m) const {
1205  Real Ay, Az, dAy, dAz;
1206  RealMatrixX at, dat, m1, dm;
1207  _A_y_moment(p, t, Ay); _A_y_moment.derivative( f, p, t, dAy);
1208  _A_z_moment(p, t, Az); _A_z_moment.derivative( f, p, t, dAz);
1209  _material_stiffness(p, t, m1); _material_stiffness.derivative( f, p, t, dm);
1210  _material_expansion(p, t, at); _material_expansion.derivative( f, p, t, dat);
1211 
1212  m = m1;
1213  m *= at;
1214  m(1,0) = dAy * m(0,0);
1215  m(0,0) *= dAz;
1216 
1217  m1 *= dat;
1218  dm *= at;
1219  m1 += dm;
1220  m1(1,0) = Ay * m1(0,0);
1221  m1(0,0) *= Az;
1222 
1223  m += m1;
1224 }
1225 
1226 
1227 
1228 
1232  const MAST::FieldFunction<Real>& A):
1233 MAST::FieldFunction<RealMatrixX>("PrestressAMatrix1D"),
1234 _prestress(prestress),
1235 _T(T),
1236 _A(A) {
1237  _functions.insert(&prestress);
1238  _functions.insert(&T);
1239  _functions.insert(&A);
1240 }
1241 
1242 
1243 
1244 
1245 void
1247 PrestressAMatrix::operator() (const libMesh::Point& p,
1248  const Real t,
1249  RealMatrixX& m) const {
1250  RealMatrixX s, T;
1251  m = RealMatrixX::Zero(2, 2);
1252  Real A;
1253  _A(p, t, A);
1254  _prestress(p, t, s);
1255  _T(p, t, T);
1256 
1257  // convert the stress to the local coordinate
1258  s *= T;
1259  s = T.transpose() * s;
1260 
1261  m(0,0) = s(0,0)*A; // only sigma_xx is applied, and torsion is neglected
1262 }
1263 
1264 
1265 
1266 
1267 
1268 void
1271  const libMesh::Point& p,
1272  const Real t,
1273  RealMatrixX& m) const {
1274  RealMatrixX s, ds, T, dT;
1275  m = RealMatrixX::Zero(2, 2);
1276  Real A, dA;
1277  _A(p, t, A); _A.derivative( f, p, t, dA);
1278  _prestress(p, t, s); _prestress.derivative( f, p, t, ds);
1279  _T(p, t, T); _T.derivative( f, p, t, dT);
1280 
1281  // convert the stress to the local coordinate
1282  s *= T;
1283  s = T.transpose() * s;
1284 
1285  // ds = dT^T s T + T^T s dT + T^T ds T
1286  RealMatrixX tmp;
1287  ds * T;
1288  ds = T.transpose()*ds;
1289 
1290  tmp = T.transpose() * s * dT + dT.transpose() * s * T;
1291  ds += tmp;
1292 
1293  m(0,0) = s(0,0)*dA + ds(0,0)*A; // only sigma_xx is applied, and torsion is neglected
1294 }
1295 
1296 
1297 
1298 
1302  const MAST::FieldFunction<Real>& A_y_moment,
1303  const MAST::FieldFunction<Real>& A_z_moment):
1304 MAST::FieldFunction<RealMatrixX>("PrestressBMatrix1D"),
1305 _prestress(prestress),
1306 _T(T),
1307 _A_y_moment(A_y_moment),
1308 _A_z_moment(A_z_moment) {
1309  _functions.insert(&prestress);
1310  _functions.insert(&T);
1311  _functions.insert(&A_y_moment);
1312  _functions.insert(&A_z_moment);
1313 }
1314 
1315 
1316 
1317 
1318 void
1320 PrestressBMatrix::operator() (const libMesh::Point& p,
1321  const Real t,
1322  RealMatrixX& m) const {
1323  RealMatrixX s, T;
1324  m = RealMatrixX::Zero(2, 2);
1325  Real Ay, Az;
1326  _A_y_moment(p, t, Ay);
1327  _A_z_moment(p, t, Az);
1328  _prestress(p, t, s);
1329  _T(p, t, T);
1330 
1331  // convert the stress to the local coordinate
1332  s = T.transpose() * s * T;
1333 
1334  // only sigma_xx is applied, and torsion is neglected
1335  m(0,0) = s(0,0)*Az;
1336  m(0,1) = s(0,0)*Ay;
1337 }
1338 
1339 
1340 
1341 
1342 
1343 
1344 void
1347  const libMesh::Point& p,
1348  const Real t,
1349  RealMatrixX& m) const {
1350  RealMatrixX s, ds, T, dT;
1351  m = RealMatrixX::Zero(2, 2);
1352  Real Ay, Az, dAy, dAz;
1353  _A_y_moment(p, t, Ay); _A_y_moment.derivative( f, p, t, dAy);
1354  _A_z_moment(p, t, Az); _A_z_moment.derivative( f, p, t, dAz);
1355  _prestress(p, t, s); _prestress.derivative( f, p, t, ds);
1356  _T(p, t, T); _T.derivative( f, p, t, dT);
1357 
1358  // convert the stress to the local coordinate
1359  s = T.transpose() * s * T;
1360 
1361  // ds = dT^T s T + T^T s dT + T^T ds T
1362  RealMatrixX tmp;
1363  ds = (T.transpose() * ds * T +
1364  T.transpose() * s * dT +
1365  dT.transpose() * s * T);
1366 
1367  // only sigma_xx is applied, and torsion is neglected
1368  m(0,0) = (s(0,0)*dAz + ds(0,0)*Az);
1369  m(0,1) = s(0,0)*dAy + ds(0,0)*Ay;
1370 }
1371 
1372 
1373 
1374 
1375 
1376 
1379  const MAST::FieldFunction<Real>& A):
1380 MAST::FieldFunction<RealMatrixX>("ThermalConductanceMatrix"),
1381 _mat_cond(mat_cond),
1382 _A(A) {
1383  _functions.insert(&mat_cond);
1384  _functions.insert(&A);
1385 }
1386 
1387 
1390 
1391 
1392 void
1394 operator() (const libMesh::Point& p,
1395  const Real t,
1396  RealMatrixX& m) const {
1397 
1398  m = RealMatrixX::Zero(1, 1);
1399  Real A;
1400  _mat_cond(p, t, m);
1401  _A(p, t, A);
1402 
1403  m *= A;
1404 }
1405 
1406 
1407 
1408 void
1410  const libMesh::Point& p,
1411  const Real t,
1412  RealMatrixX& m) const {
1413  m = RealMatrixX::Zero(1, 1);
1414  RealMatrixX dm;
1415  Real A, dA;
1416  _mat_cond(p, t, m);
1417  _mat_cond.derivative( f, p, t, dm);
1418  _A(p, t, A);
1419  _A.derivative( f, p, t, dA);
1420 
1421  m *= dA;
1422  m += dm*A;
1423 }
1424 
1425 
1426 
1427 
1430  const MAST::FieldFunction<Real>& h):
1431 MAST::FieldFunction<RealMatrixX>("ThermalCapacitanceMatrix"),
1432 _mat_cap(mat_cap),
1433 _h(h) {
1434  _functions.insert(&mat_cap);
1435  _functions.insert(&h);
1436 }
1437 
1438 
1439 
1440 
1443 
1444 
1445 void
1447 operator() (const libMesh::Point& p,
1448  const Real t,
1449  RealMatrixX& m) const {
1450 
1451  m = RealMatrixX::Zero(1, 1);
1452  Real h;
1453  _mat_cap(p, t, m);
1454  _h(p, t, h);
1455 
1456  m *= h;
1457 }
1458 
1459 
1460 
1461 void
1463  const libMesh::Point& p,
1464  const Real t,
1465  RealMatrixX& m) const {
1466  m = RealMatrixX::Zero(1, 1);
1467  RealMatrixX dm;
1468  Real h, dh;
1469  _mat_cap(p, t, m);
1470  _mat_cap.derivative( f, p, t, dm);
1471  _h(p, t, h);
1472  _h.derivative( f, p, t, dh);
1473 
1474  m *= dh;
1475  m += dm*h;
1476 }
1477 
1478 
1479 
1480 
1481 bool
1484  return _material->depends_on(f) || // check if the material property depends on the function
1485  MAST::ElementPropertyCardBase::depends_on(f); // check with this property card
1486 }
1487 
1488 
1489 
1492 
1493  libmesh_assert(_initialized);
1494  return *_A;
1495 }
1496 
1497 
1498 
1501 
1502  libmesh_assert(_initialized);
1503  return *_J;
1504 }
1505 
1506 
1509 
1510  libmesh_assert(_initialized);
1511  return *_Ip;
1512 }
1513 
1514 
1517 
1518  libmesh_assert(_initialized);
1519  return *_Ay;
1520 }
1521 
1522 
1525 
1526  libmesh_assert(_initialized);
1527  return *_Az;
1528 }
1529 
1530 
1533 
1534  libmesh_assert(_initialized);
1535  return *_AI;
1536 }
1537 
1538 
1541 
1542  libmesh_assert(_initialized);
1543  return *_Kappa;
1544 }
1545 
1546 
1549 
1550  libmesh_assert(_initialized);
1551  return *_Gamma;
1552 }
1553 
1554 
1557 
1558  libmesh_assert(_initialized);
1559  return *_A;
1560 }
1561 
1562 
1563 
1566 
1567  libmesh_assert(_initialized);
1568  return *_J;
1569 }
1570 
1571 
1574 
1575  libmesh_assert(_initialized);
1576  return *_Ip;
1577 }
1578 
1579 
1582 
1583  libmesh_assert(_initialized);
1584  return *_Ay;
1585 }
1586 
1587 
1590 
1591  libmesh_assert(_initialized);
1592  return *_Az;
1593 }
1594 
1595 
1598 
1599  libmesh_assert(_initialized);
1600  return *_AI;
1601 }
1602 
1603 
1606 
1607  libmesh_assert(_initialized);
1608  return *_Kappa;
1609 }
1610 
1611 
1614 
1615  libmesh_assert(_initialized);
1616  return *_Gamma;
1617 }
1618 
1619 
1620 void
1622 
1623  _A.reset();
1624  _Ay.reset();
1625  _Az.reset();
1626  _J.reset();
1627  _Ip.reset();
1628  _AI.reset();
1629  _Gamma.reset();
1630  _Kappa.reset();
1631 
1632  _initialized = false;
1633 }
1634 
1635 
1636 
1637 void
1639 
1640  libmesh_assert(!_initialized);
1641 
1643  &hy = this->get<MAST::FieldFunction<Real> >("hy"),
1644  &hz = this->get<MAST::FieldFunction<Real> >("hz"),
1645  &Kappazz = this->get<MAST::FieldFunction<Real> >("Kappazz"),
1646  &Kappayy = this->get<MAST::FieldFunction<Real> >("Kappayy"),
1647  &hy_off = this->get<MAST::FieldFunction<Real> >("hy_off"),
1648  &hz_off = this->get<MAST::FieldFunction<Real> >("hz_off");
1649 
1650  _A.reset(new MAST::Solid1DSectionProperty::Area(hy,
1651  hz));
1653  hz,
1654  hz_off));
1656  hz,
1657  hy_off));
1659  hz));
1661  hz,
1662  hy_off,
1663  hz_off));
1665  hz,
1666  hy_off,
1667  hz_off));
1668 
1669  _Kappa.reset(new MAST::Solid1DSectionProperty::ShearCoefficientMatrix(Kappazz, Kappayy));
1670 
1672 
1673  _initialized = true;
1674 }
1675 
1676 
1677 std::unique_ptr<MAST::FieldFunction<RealMatrixX> >
1680 
1681  // make sure that the init method has been called on the card
1682  libmesh_assert(_initialized);
1683 
1686  (_material->stiffness_matrix(1), *_A, *_J);
1687 
1688  return std::unique_ptr<MAST::FieldFunction<RealMatrixX> > (rval);
1689 }
1690 
1691 
1692 std::unique_ptr<MAST::FieldFunction<RealMatrixX> >
1695 
1696  // make sure that the init method has been called on the card
1697  libmesh_assert(_initialized);
1698 
1701  (_material->stiffness_matrix(1), *_A, *_J);
1702 
1703  return std::unique_ptr<MAST::FieldFunction<RealMatrixX> > (rval);
1704 }
1705 
1706 
1707 std::unique_ptr<MAST::FieldFunction<RealMatrixX> >
1710 
1711  // make sure that the init method has been called on the card
1712  libmesh_assert(_initialized);
1713 
1716  (_material->stiffness_matrix(1), *_Ay, *_Az);
1717 
1718  return std::unique_ptr<MAST::FieldFunction<RealMatrixX> > (rval);
1719 }
1720 
1721 
1722 std::unique_ptr<MAST::FieldFunction<RealMatrixX> >
1725 
1726  // make sure that the init method has been called on the card
1727  libmesh_assert(_initialized);
1728 
1731  (_material->stiffness_matrix(1), *_Ay, *_Az);
1732 
1733  return std::unique_ptr<MAST::FieldFunction<RealMatrixX> > (rval);
1734 }
1735 
1736 
1737 std::unique_ptr<MAST::FieldFunction<RealMatrixX> >
1740 
1741  // make sure that the init method has been called on the card
1742  libmesh_assert(_initialized);
1743 
1746  (_material->stiffness_matrix(1), *_AI);
1747 
1748  return std::unique_ptr<MAST::FieldFunction<RealMatrixX> > (rval);
1749 }
1750 
1751 
1752 std::unique_ptr<MAST::FieldFunction<RealMatrixX> >
1755 
1756  // make sure that the init method has been called on the card
1757  libmesh_assert(_initialized);
1758 
1761  (_material->stiffness_matrix(1), *_AI);
1762 
1763  return std::unique_ptr<MAST::FieldFunction<RealMatrixX> > (rval);
1764 }
1765 
1766 
1767 std::unique_ptr<MAST::FieldFunction<RealMatrixX> >
1770 
1771  libmesh_error();
1772 
1773  return std::unique_ptr<MAST::FieldFunction<RealMatrixX> > (nullptr);
1774 }
1775 
1776 
1777 std::unique_ptr<MAST::FieldFunction<RealMatrixX> >
1780 
1781  libmesh_error();
1782 
1783  return std::unique_ptr<MAST::FieldFunction<RealMatrixX> > (nullptr);
1784 }
1785 
1786 
1787 std::unique_ptr<MAST::FieldFunction<RealMatrixX> >
1790 
1791  // make sure that the init method has been called on the card
1792  libmesh_assert(_initialized);
1793 
1796  (_material->get<FieldFunction<Real> >("rho"),
1797  *_A,
1798  *_Ay,
1799  *_Az,
1800  *_Ip,
1801  *_AI);
1802 
1803  return std::unique_ptr<MAST::FieldFunction<RealMatrixX> > (rval);
1804 }
1805 
1806 
1807 std::unique_ptr<MAST::FieldFunction<RealMatrixX> >
1810 
1811  // make sure that the init method has been called on the card
1812  libmesh_assert(_initialized);
1813 
1816  (_material->get<FieldFunction<Real> >("rho"),
1817  *_A,
1818  *_Ay,
1819  *_Az,
1820  *_Ip,
1821  *_AI);
1822 
1823  return std::unique_ptr<MAST::FieldFunction<RealMatrixX> > (rval);
1824 }
1825 
1826 
1827 std::unique_ptr<MAST::FieldFunction<RealMatrixX> >
1830 
1831  // make sure that the init method has been called on the card
1832  libmesh_assert(_initialized);
1833 
1836  (_material->stiffness_matrix(1),
1837  _material->thermal_expansion_matrix(1),
1838  *_A);
1839 
1840  return std::unique_ptr<MAST::FieldFunction<RealMatrixX> > (rval);
1841 }
1842 
1843 
1844 std::unique_ptr<MAST::FieldFunction<RealMatrixX> >
1847 
1848  // make sure that the init method has been called on the card
1849  libmesh_assert(_initialized);
1850 
1853  (_material->stiffness_matrix(1),
1854  _material->thermal_expansion_matrix(1),
1855  *_A);
1856 
1857  return std::unique_ptr<MAST::FieldFunction<RealMatrixX> > (rval);
1858 }
1859 
1860 
1861 std::unique_ptr<MAST::FieldFunction<RealMatrixX> >
1864 
1865 
1866  // make sure that the init method has been called on the card
1867  libmesh_assert(_initialized);
1868 
1871  (_material->stiffness_matrix(1),
1872  _material->thermal_expansion_matrix(1),
1873  *_Ay,
1874  *_Az);
1875 
1876  return std::unique_ptr<MAST::FieldFunction<RealMatrixX> > (rval);
1877 }
1878 
1879 
1880 std::unique_ptr<MAST::FieldFunction<RealMatrixX> >
1883 
1884 
1885  // make sure that the init method has been called on the card
1886  libmesh_assert(_initialized);
1887 
1890  (_material->stiffness_matrix(1),
1891  _material->thermal_expansion_matrix(1),
1892  *_Ay,
1893  *_Az);
1894 
1895  return std::unique_ptr<MAST::FieldFunction<RealMatrixX> > (rval);
1896 }
1897 
1898 
1899 std::unique_ptr<MAST::FieldFunction<RealMatrixX> >
1902 
1903 
1904  // make sure that the init method has been called on the card
1905  libmesh_assert(_initialized);
1906 
1909  (_material->transverse_shear_stiffness_matrix(),
1910  *_A, *_Kappa);
1911 
1912  return std::unique_ptr<MAST::FieldFunction<RealMatrixX> > (rval);
1913 }
1914 
1915 
1916 std::unique_ptr<MAST::FieldFunction<RealMatrixX> >
1919 
1920 
1921  // make sure that the init method has been called on the card
1922  libmesh_assert(_initialized);
1923 
1926  (_material->transverse_shear_stiffness_matrix(),
1927  *_A, *_Kappa);
1928 
1929  return std::unique_ptr<MAST::FieldFunction<RealMatrixX> > (rval);
1930 }
1931 
1932 
1933 std::unique_ptr<MAST::FieldFunction<RealMatrixX> >
1936 
1937  // make sure that the init method has been called on the card
1938  libmesh_assert(_initialized);
1939 
1941  // TODO: figure out the interface for prestress and T matrix
1942  libmesh_assert(false);
1943  // = new MAST::Solid1DSectionProperty::PrestressAMatrix
1944  //(this->get<MAST::FieldFunction<RealMatrixX> >("prestress"),
1945  // e.local_elem().T_matrix(),
1946  // *_A);
1947 
1948  return std::unique_ptr<MAST::FieldFunction<RealMatrixX> > (rval);
1949 }
1950 
1951 
1952 std::unique_ptr<MAST::FieldFunction<RealMatrixX> >
1955 
1956  // make sure that the init method has been called on the card
1957  libmesh_assert(_initialized);
1958 
1960  // TODO: figure out the interface for prestress and T matrix
1961  libmesh_assert(false);
1962  // = new MAST::Solid1DSectionProperty::PrestressAMatrix
1963  //(this->get<MAST::FieldFunction<RealMatrixX> >("prestress"),
1964  // e.local_elem().T_matrix(),
1965  // *_A);
1966 
1967  return std::unique_ptr<MAST::FieldFunction<RealMatrixX> > (rval);
1968 }
1969 
1970 
1971 std::unique_ptr<MAST::FieldFunction<RealMatrixX> >
1974 
1975  // make sure that the init method has been called on the card
1976  libmesh_assert(_initialized);
1977 
1979  // TODO: figure out the interface for prestress and T matrix
1980  libmesh_assert(false);
1981  // = new MAST::Solid1DSectionProperty::PrestressBMatrix
1982  //(this->get<MAST::FieldFunction<RealMatrixX> >("prestress"),
1983  // e.local_elem().T_matrix(),
1984  // *_Ay,
1985  // *_Az);
1986 
1987  return std::unique_ptr<MAST::FieldFunction<RealMatrixX> > (rval);
1988 }
1989 
1990 
1991 std::unique_ptr<MAST::FieldFunction<RealMatrixX> >
1994 
1995  // make sure that the init method has been called on the card
1996  libmesh_assert(_initialized);
1997 
1999  // TODO: figure out the interface for prestress and T matrix
2000  libmesh_assert(false);
2001  // = new MAST::Solid1DSectionProperty::PrestressBMatrix
2002  //(this->get<MAST::FieldFunction<RealMatrixX> >("prestress"),
2003  // e.local_elem().T_matrix(),
2004  // *_Ay,
2005  // *_Az);
2006 
2007  return std::unique_ptr<MAST::FieldFunction<RealMatrixX> > (rval);
2008 }
2009 
2010 
2011 std::unique_ptr<MAST::FieldFunction<RealMatrixX> >
2014 
2015  // make sure that the init method has been called on the card
2016  libmesh_assert(_initialized);
2017 
2020  (_material->conductance_matrix(1),
2021  *_A);
2022 
2023  return std::unique_ptr<MAST::FieldFunction<RealMatrixX> > (rval);
2024 }
2025 
2026 
2027 std::unique_ptr<MAST::FieldFunction<RealMatrixX> >
2030 
2031  // make sure that the init method has been called on the card
2032  libmesh_assert(_initialized);
2033 
2036  (_material->conductance_matrix(1),
2037  *_A);
2038 
2039  return std::unique_ptr<MAST::FieldFunction<RealMatrixX> > (rval);
2040 }
2041 
2042 
2043 std::unique_ptr<MAST::FieldFunction<RealMatrixX> >
2046 
2047  // make sure that the init method has been called on the card
2048  libmesh_assert(_initialized);
2049 
2052  (_material->capacitance_matrix(1),
2053  *_A);
2054 
2055  return std::unique_ptr<MAST::FieldFunction<RealMatrixX> > (rval);
2056 }
2057 
2058 
2059 std::unique_ptr<MAST::FieldFunction<RealMatrixX> >
2062 
2063  // make sure that the init method has been called on the card
2064  libmesh_assert(_initialized);
2065 
2068  (_material->capacitance_matrix(1),
2069  *_A);
2070 
2071  return std::unique_ptr<MAST::FieldFunction<RealMatrixX> > (rval);
2072 }
2073 
2074 
2077 section(const MAST::ElementBase& e) const {
2078 
2079  return _A.get();
2080 }
2081 
2082 
2085 section() const {
2086 
2087  return _A.get();
2088 }
virtual void operator()(const libMesh::Point &p, const Real t, Real &m) const
calculates the value of the function at the specified point, p, and time, t, and returns it in v...
calculates the area moment about the Y-axis due to an offset along the Z-axis
BendingStiffnessMatrix(const MAST::FieldFunction< RealMatrixX > &mat, const MAST::FieldFunction< RealMatrixX > &I)
ThermalCapacitanceMatrix(const MAST::FieldFunction< RealMatrixX > &mat_cond, const MAST::FieldFunction< Real > &h)
virtual std::unique_ptr< MAST::FieldFunction< RealMatrixX > > transverse_shear_stiffness_matrix() const
Area(const MAST::FieldFunction< Real > &hy, const MAST::FieldFunction< Real > &hz)
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 std::unique_ptr< MAST::FieldFunction< RealMatrixX > > stiffness_D_matrix() const
virtual const MAST::FieldFunction< Real > & Ip() const
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 const MAST::FieldFunction< Real > * section() const
ExtensionBendingStiffnessMatrix(const MAST::FieldFunction< RealMatrixX > &mat, const MAST::FieldFunction< Real > &A_y_moment, const MAST::FieldFunction< Real > &A_z_moment)
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...
ExtensionStiffnessMatrix(const MAST::FieldFunction< RealMatrixX > &mat, const MAST::FieldFunction< Real > &A, const MAST::FieldFunction< Real > &J)
AreaZMoment(const MAST::FieldFunction< Real > &hy, const MAST::FieldFunction< Real > &hz, const MAST::FieldFunction< Real > &hy_offset)
virtual void derivative(const MAST::FunctionBase &f, const libMesh::Point &p, const Real t, Real &m) const
calculates the value of the derivative of function with respect to the function f at the specified po...
virtual std::unique_ptr< MAST::FieldFunction< RealMatrixX > > prestress_A_matrix() const
PrestressBMatrix(const MAST::FieldFunction< RealMatrixX > &prestress, const MAST::FieldFunction< RealMatrixX > &T, const MAST::FieldFunction< Real > &A_y_moment, const MAST::FieldFunction< Real > &A_z_moment)
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...
std::set< const MAST::FunctionBase * > _functions
set of functions that this function depends on
virtual std::unique_ptr< MAST::FieldFunction< RealMatrixX > > thermal_capacitance_matrix() const
InertiaMatrix(const MAST::FieldFunction< Real > &rho, const MAST::FieldFunction< Real > &A, const MAST::FieldFunction< Real > &A_y_moment, const MAST::FieldFunction< Real > &A_z_moment, const MAST::FieldFunction< Real > &Ip, const MAST::FieldFunction< RealMatrixX > &I)
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, ValType &v) const
calculates the value of the function derivative and returns it in v.
libMesh::Real Real
virtual std::unique_ptr< MAST::FieldFunction< RealMatrixX > > inertia_matrix() const
PolarInertia(const MAST::FieldFunction< Real > &hy, const MAST::FieldFunction< Real > &hz, const MAST::FieldFunction< Real > &hy_offset, const MAST::FieldFunction< Real > &hz_offset)
virtual void derivative(const MAST::FunctionBase &f, const libMesh::Point &p, const Real t, Real &m) const
calculates the value of the derivative of function with respect to the function f at the specified po...
virtual const MAST::FieldFunction< RealMatrixX > & I() const
calculates the area moment about the Z-axis due to an offset along the Y-axis
virtual bool depends_on(const MAST::FunctionBase &f) const
returns true if the property card depends on the function f
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 > > stiffness_B_matrix() 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 > > stiffness_A_matrix() const
ThermalExpansionAMatrix(const MAST::FieldFunction< RealMatrixX > &mat_stiff, const MAST::FieldFunction< RealMatrixX > &mat_expansion, const MAST::FieldFunction< Real > &A)
TorsionalConstant(const MAST::FieldFunction< Real > &hy, const MAST::FieldFunction< Real > &hz)
Matrix< Real, Dynamic, Dynamic > RealMatrixX
virtual std::unique_ptr< MAST::FieldFunction< RealMatrixX > > prestress_B_matrix() const
virtual const MAST::FieldFunction< Real > & Az() const
virtual const MAST::FieldFunction< Real > & Gam() const
AreaInertiaMatrix(const MAST::FieldFunction< Real > &hy, const MAST::FieldFunction< Real > &hz, const MAST::FieldFunction< Real > &hy_offset, const MAST::FieldFunction< Real > &hz_offset)
virtual std::unique_ptr< MAST::FieldFunction< RealMatrixX > > thermal_expansion_A_matrix() const
ThermalExpansionBMatrix(const MAST::FieldFunction< RealMatrixX > &mat_stiff, const MAST::FieldFunction< RealMatrixX > &mat_expansion, const MAST::FieldFunction< Real > &A_y_moment, const MAST::FieldFunction< Real > &A_z_moment)
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 std::unique_ptr< MAST::FieldFunction< RealMatrixX > > thermal_conductance_matrix() const
calculates the 2x2 matrix of area inertia for the section with individual entries as ...
virtual const MAST::FieldFunction< Real > & A() const
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 std::unique_ptr< MAST::FieldFunction< RealMatrixX > > thermal_expansion_B_matrix() const
ThermalConductanceMatrix(const MAST::FieldFunction< RealMatrixX > &mat_cond, const MAST::FieldFunction< Real > &A)
virtual void derivative(const MAST::FunctionBase &f, const libMesh::Point &p, const Real t, Real &m) const
calculates the value of the derivative of function with respect to the function f at the specified po...
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 const MAST::FieldFunction< Real > & Ay() 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 derivative(const MAST::FunctionBase &f, const libMesh::Point &p, const Real t, Real &m) const
calculates the value of the derivative of function with respect to the function f at the specified po...
TransverseStiffnessMatrix(const MAST::FieldFunction< RealMatrixX > &mat, const MAST::FieldFunction< Real > &A, const MAST::FieldFunction< RealMatrixX > &Kappa)
virtual void derivative(const MAST::FunctionBase &f, const libMesh::Point &p, const Real t, Real &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 void derivative(const MAST::FunctionBase &f, const libMesh::Point &p, const Real t, Real &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 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...
PrestressAMatrix(const MAST::FieldFunction< RealMatrixX > &prestress, const MAST::FieldFunction< RealMatrixX > &T, const MAST::FieldFunction< Real > &A)
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...
AreaYMoment(const MAST::FieldFunction< Real > &hy, const MAST::FieldFunction< Real > &hz, const MAST::FieldFunction< Real > &hz_offset)
virtual std::unique_ptr< MAST::FieldFunction< RealMatrixX > > damping_matrix() const
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 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 const MAST::FieldFunction< RealMatrixX > & Kap() const
ShearCoefficientMatrix(const MAST::FieldFunction< Real > &Kzz, const MAST::FieldFunction< Real > &Kyy)
virtual const MAST::FieldFunction< Real > & J() const
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 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 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 bool depends_on(const MAST::FunctionBase &f) const
returns true if the property card depends on the function f
This is the base class for elements that implement calculation of finite element quantities over the ...
Definition: elem_base.h:72