linear_elasticity_elements.h
Go to the documentation of this file.
1// LIC// ====================================================================
2// LIC// This file forms part of oomph-lib, the object-oriented,
3// LIC// multi-physics finite-element library, available
4// LIC// at http://www.oomph-lib.org.
5// LIC//
6// LIC// Copyright (C) 2006-2022 Matthias Heil and Andrew Hazel
7// LIC//
8// LIC// This library is free software; you can redistribute it and/or
9// LIC// modify it under the terms of the GNU Lesser General Public
10// LIC// License as published by the Free Software Foundation; either
11// LIC// version 2.1 of the License, or (at your option) any later version.
12// LIC//
13// LIC// This library is distributed in the hope that it will be useful,
14// LIC// but WITHOUT ANY WARRANTY; without even the implied warranty of
15// LIC// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
16// LIC// Lesser General Public License for more details.
17// LIC//
18// LIC// You should have received a copy of the GNU Lesser General Public
19// LIC// License along with this library; if not, write to the Free Software
20// LIC// Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
21// LIC// 02110-1301 USA.
22// LIC//
23// LIC// The authors may be contacted at oomph-lib@maths.man.ac.uk.
24// LIC//
25// LIC//====================================================================
26// Header file for general linear elasticity elements
27
28// Include guards to prevent multiple inclusion of the header
29#ifndef OOMPH_LINEAR_ELASTICITY_ELEMENTS_HEADER
30#define OOMPH_LINEAR_ELASTICITY_ELEMENTS_HEADER
31
32// Config header generated by autoconfig
33#ifdef HAVE_CONFIG_H
34#include <oomph-lib-config.h>
35#endif
36
37
38// OOMPH-LIB headers
39#include "../generic/Qelements.h"
40#include "../generic/mesh.h"
41#include "../generic/hermite_elements.h"
42#include "./elasticity_tensor.h"
43#include "../generic/projection.h"
44
45namespace oomph
46{
47 //=======================================================================
48 /// A base class for elements that solve the equations of linear
49 /// elasticity in Cartesian coordinates.
50 /// Combines a few generic functions that are shared by
51 /// LinearElasticityEquations
52 /// and LinearElasticityEquationsWithPressure (hierher: The latter
53 /// don't exist yet but will be written as soon as somebody needs them...)
54 //=======================================================================
55 template<unsigned DIM>
57 {
58 public:
59 /// Return the index at which the i-th unknown displacement
60 /// component is stored. The default value, i, is appropriate for
61 /// single-physics problems.
62 virtual inline unsigned u_index_linear_elasticity(const unsigned i) const
63 {
64 return i;
65 }
66
67 /// d^2u/dt^2 at local node n
68 double d2u_dt2_linear_elasticity(const unsigned& n, const unsigned& i) const
69 {
70 // Get the timestepper
72
73 // Storage for the derivative - initialise to 0
74 double d2u_dt2 = 0.0;
75
76 // If we are doing an unsteady solve then calculate the derivative
78 {
79 // Get the nodal index
80 const unsigned u_nodal_index = u_index_linear_elasticity(i);
81
82 // Get the number of values required to represent history
83 const unsigned n_time = time_stepper_pt->ntstorage();
84
85 // Loop over history values
86 for (unsigned t = 0; t < n_time; t++)
87 {
88 // Add the contribution to the derivative
89 d2u_dt2 +=
90 time_stepper_pt->weight(2, t) * nodal_value(t, n, u_nodal_index);
91 }
92 }
93
94 return d2u_dt2;
95 }
96
97 /// Compute vector of FE interpolated displacement u at local coordinate s
99 Vector<double>& disp) const
100 {
101 // Find number of nodes
102 unsigned n_node = nnode();
103
104 // Local shape function
105 Shape psi(n_node);
106
107 // Find values of shape function
108 shape(s, psi);
109
110 for (unsigned i = 0; i < DIM; i++)
111 {
112 // Index at which the nodal value is stored
113 unsigned u_nodal_index = u_index_linear_elasticity(i);
114
115 // Initialise value of u
116 disp[i] = 0.0;
117
118 // Loop over the local nodes and sum
119 for (unsigned l = 0; l < n_node; l++)
120 {
121 disp[i] += nodal_value(l, u_nodal_index) * psi[l];
122 }
123 }
124 }
125
126 /// Return FE interpolated displacement u[i] at local coordinate s
128 const unsigned& i) const
129 {
130 // Find number of nodes
131 unsigned n_node = nnode();
132
133 // Local shape function
134 Shape psi(n_node);
135
136 // Find values of shape function
137 shape(s, psi);
138
139 // Get nodal index at which i-th velocity is stored
140 unsigned u_nodal_index = u_index_linear_elasticity(i);
141
142 // Initialise value of u
143 double interpolated_u = 0.0;
144
145 // Loop over the local nodes and sum
146 for (unsigned l = 0; l < n_node; l++)
147 {
148 interpolated_u += nodal_value(l, u_nodal_index) * psi[l];
149 }
150
151 return (interpolated_u);
152 }
153
154
155 /// Function pointer to function that specifies the body force
156 /// as a function of the Cartesian coordinates and time FCT(t,x,b) --
157 /// x and b are Vectors!
158 typedef void (*BodyForceFctPt)(const double& t,
159 const Vector<double>& x,
160 Vector<double>& b);
161
162 /// Constructor: Set null pointers for constitutive law and for
163 /// isotropic growth function. Set physical parameter values to
164 /// default values, switch on inertia and set body force to zero.
168 Unsteady(true),
170 {
171 }
172
173 /// Return the pointer to the elasticity_tensor
175 {
177 }
178
179 /// Access function to the entries in the elasticity tensor
180 inline double E(const unsigned& i,
181 const unsigned& j,
182 const unsigned& k,
183 const unsigned& l) const
184 {
185 return (*Elasticity_tensor_pt)(i, j, k, l);
186 }
187
188 /// Access function for timescale ratio (nondim density)
189 const double& lambda_sq() const
190 {
191 return *Lambda_sq_pt;
192 }
193
194 /// Access function for pointer to timescale ratio (nondim density)
195 double*& lambda_sq_pt()
196 {
197 return Lambda_sq_pt;
198 }
199
200 /// Access function: Pointer to body force function
202 {
203 return Body_force_fct_pt;
204 }
205
206 /// Access function: Pointer to body force function (const version)
208 {
209 return Body_force_fct_pt;
210 }
211
212
213 /// Switch on solid inertia
215 {
216 Unsteady = true;
217 }
218
219 /// Switch off solid inertia
221 {
222 Unsteady = false;
223 }
224
225 /// Access function to flag that switches inertia on/off (const version)
227 {
228 return Unsteady;
229 }
230
231 /// Pin the element's redundant solid pressures (needed for refinement)
233
234 /// Loop over all elements in Vector (which typically contains
235 /// all the elements in a refineable solid mesh) and pin the nodal solid
236 /// pressure degrees of freedom that are not being used. Function uses
237 /// the member function
238 /// - \c LinearElasticityEquationsBase<DIM>::
239 /// pin_elemental_redundant_nodal_pressure_dofs()
240 /// .
241 /// which is empty by default and should be implemented for
242 /// elements with nodal solid pressure degrees of freedom
243 /// (e.g. linear elasticity elements with continuous pressure
244 /// interpolation.)
246 const Vector<GeneralisedElement*>& element_pt)
247 {
248 // Loop over all elements and call the function that pins their
249 // unused nodal solid pressure data
250 unsigned n_element = element_pt.size();
251 for (unsigned e = 0; e < n_element; e++)
252 {
253 dynamic_cast<LinearElasticityEquationsBase<DIM>*>(element_pt[e])
255 }
256 }
257
258 /// Return the Cauchy stress tensor, as calculated
259 /// from the elasticity tensor at specified local coordinate
260 /// Virtual so separaete versions can (and must!) be provided
261 /// for displacement and pressure-displacement formulations.
262 virtual void get_stress(const Vector<double>& s,
263 DenseMatrix<double>& sigma) const = 0;
264
265 /// Return the strain tensor
266 void get_strain(const Vector<double>& s, DenseMatrix<double>& strain) const;
267
268 /// Evaluate body force at Eulerian coordinate x at present time
269 /// (returns zero vector if no body force function pointer has been set)
270 inline void body_force(const Vector<double>& x, Vector<double>& b) const
271 {
272 // If no function has been set, return zero vector
273 if (Body_force_fct_pt == 0)
274 {
275 // Get spatial dimension of element
276 unsigned n = dim();
277 for (unsigned i = 0; i < n; i++)
278 {
279 b[i] = 0.0;
280 }
281 }
282 else
283 {
284 // Get time from timestepper of first node (note that this must
285 // work -- body force only makes sense for elements that can be
286 // deformed and given that the deformation of solid finite elements
287 // is controlled by their nodes, nodes must exist!)
288 double time = node_pt(0)->time_stepper_pt()->time_pt()->time();
289
290 // Now evaluate the body force
291 (*Body_force_fct_pt)(time, x, b);
292 }
293 }
294
295
296 /// The number of "DOF types" that degrees of freedom in this element
297 /// are sub-divided into: for now lump them all into one DOF type.
298 /// Can be adjusted later
299 unsigned ndof_types() const
300 {
301 return 2;
302 // return 1;
303 }
304
305 /// Create a list of pairs for all unknowns in this element,
306 /// so that the first entry in each pair contains the global equation
307 /// number of the unknown, while the second one contains the number
308 /// of the "DOF types" that this unknown is associated with.
309 /// (Function can obviously only be called if the equation numbering
310 /// scheme has been set up.)
312 std::list<std::pair<unsigned long, unsigned>>& dof_lookup_list) const
313 {
314 // temporary pair (used to store dof lookup prior to being added
315 // to list)
316 std::pair<unsigned long, unsigned> dof_lookup;
317
318 // number of nodes
319 const unsigned n_node = this->nnode();
320
321 // Integer storage for local unknown
322 int local_unknown = 0;
323
324 // Loop over the nodes
325 for (unsigned n = 0; n < n_node; n++)
326 {
327 // Loop over dimension
328 for (unsigned i = 0; i < DIM; i++)
329 {
330 // If the variable is free
331 local_unknown = nodal_local_eqn(n, i);
332
333 // ignore pinned values
334 if (local_unknown >= 0)
335 {
336 // store dof lookup in temporary pair: First entry in pair
337 // is global equation number; second entry is dof type
338 dof_lookup.first = this->eqn_number(local_unknown);
339 dof_lookup.second = i;
340 // dof_lookup.second = DIM;
341
342 // add to list
343 dof_lookup_list.push_front(dof_lookup);
344 }
345 }
346 }
347 }
348
349
350 protected:
351 /// Pointer to the elasticity tensor
353
354 /// Timescale ratio (non-dim. density)
356
357 /// Flag that switches inertia on/off
359
360 /// Pointer to body force function
362
363 /// Static default value for timescale ratio (1.0 -- for natural scaling)
365 };
366
367
368 /// ////////////////////////////////////////////////////////////////////
369 /// ////////////////////////////////////////////////////////////////////
370 /// ////////////////////////////////////////////////////////////////////
371
372
373 //=======================================================================
374 /// A class for elements that solve the equations of linear elasticity
375 /// in cartesian coordinates.
376 //=======================================================================
377 template<unsigned DIM>
379 : public virtual LinearElasticityEquationsBase<DIM>
380 {
381 public:
382 /// Constructor
384
385 /// Number of values required at node n.
386 unsigned required_nvalue(const unsigned& n) const
387 {
388 return DIM;
389 }
390
391 /// Return the residuals for the solid equations (the discretised
392 /// principle of virtual displacements)
394 {
397 }
398
399 /// The jacobian is calculated by finite differences by default,
400 /// We need only to take finite differences w.r.t. positional variables
401 /// For this element
403 DenseMatrix<double>& jacobian)
404 {
405 // Add the contribution to the residuals
407 residuals, jacobian, 1);
408 }
409
410 /// Return the Cauchy stress tensor, as calculated
411 /// from the elasticity tensor at specified local coordinate
412 void get_stress(const Vector<double>& s, DenseMatrix<double>& sigma) const;
413
414
415 /// Output exact solution x,y,[z],u,v,[w]
416 void output_fct(std::ostream& outfile,
417 const unsigned& nplot,
419
420 /// Output exact solution x,y,[z],u,v,[w] (unsteady version)
421 void output_fct(std::ostream& outfile,
422 const unsigned& nplot,
423 const double& time,
425
426 /// Output: x,y,[z],u,v,[w]
427 void output(std::ostream& outfile)
428 {
429 unsigned n_plot = 5;
430 output(outfile, n_plot);
431 }
432
433 /// Output: x,y,[z],u,v,[w]
434 void output(std::ostream& outfile, const unsigned& n_plot);
435
436
437 /// C-style output: x,y,[z],u,v,[w]
438 void output(FILE* file_pt)
439 {
440 unsigned n_plot = 5;
441 output(file_pt, n_plot);
442 }
443
444 /// Output: x,y,[z],u,v,[w]
445 void output(FILE* file_pt, const unsigned& n_plot);
446
447 /// Validate against exact solution.
448 /// Solution is provided via function pointer.
449 /// Plot at a given number of plot points and compute L2 error
450 /// and L2 norm of displacement solution over element
451 void compute_error(std::ostream& outfile,
453 double& error,
454 double& norm);
455
456 /// Validate against exact solution.
457 /// Solution is provided via function pointer.
458 /// Plot at a given number of plot points and compute L2 error
459 /// and L2 norm of displacement solution over element
460 /// (unsteady version)
461 void compute_error(std::ostream& outfile,
463 const double& time,
464 double& error,
465 double& norm);
466
467 private:
468 /// Private helper function to compute residuals and (if requested
469 /// via flag) also the Jacobian matrix.
471 Vector<double>& residuals, DenseMatrix<double>& jacobian, unsigned flag);
472 };
473
474
475 /// /////////////////////////////////////////////////////////////////////
476 /// /////////////////////////////////////////////////////////////////////
477 /// /////////////////////////////////////////////////////////////////////
478
479
480 //===========================================================================
481 /// An Element that solves the equations of linear elasticity
482 /// in Cartesian coordinates, using QElements for the geometry
483 //============================================================================
484 template<unsigned DIM, unsigned NNODE_1D>
485 class QLinearElasticityElement : public virtual QElement<DIM, NNODE_1D>,
486 public virtual LinearElasticityEquations<DIM>
487 {
488 public:
489 /// Constructor
491 : QElement<DIM, NNODE_1D>(), LinearElasticityEquations<DIM>()
492 {
493 }
494
495 /// Output exact solution x,y,[z],u,v,[w]
496 void output_fct(std::ostream& outfile,
497 const unsigned& nplot,
499 {
500 LinearElasticityEquations<DIM>::output_fct(outfile, nplot, exact_soln_pt);
501 }
502
503 /// Output function
504 void output(std::ostream& outfile)
505 {
507 }
508
509 /// Output function
510 void output(std::ostream& outfile, const unsigned& n_plot)
511 {
513 }
514
515
516 /// C-style output function
517 void output(FILE* file_pt)
518 {
520 }
521
522 /// C-style output function
523 void output(FILE* file_pt, const unsigned& n_plot)
524 {
526 }
527 };
528
529
530 //============================================================================
531 /// FaceGeometry of a linear 2D QLinearElasticityElement element
532 //============================================================================
533 template<>
535 : public virtual QElement<1, 2>
536 {
537 public:
538 /// Constructor must call the constructor of the underlying solid element
539 FaceGeometry() : QElement<1, 2>() {}
540 };
541
542
543 //============================================================================
544 /// FaceGeometry of a quadratic 2D QLinearElasticityElement element
545 //============================================================================
546 template<>
548 : public virtual QElement<1, 3>
549 {
550 public:
551 /// Constructor must call the constructor of the underlying element
552 FaceGeometry() : QElement<1, 3>() {}
553 };
554
555
556 //============================================================================
557 /// FaceGeometry of a cubic 2D QLinearElasticityElement element
558 //============================================================================
559 template<>
561 : public virtual QElement<1, 4>
562 {
563 public:
564 /// Constructor must call the constructor of the underlying element
565 FaceGeometry() : QElement<1, 4>() {}
566 };
567
568
569 //============================================================================
570 /// FaceGeometry of a linear 3D QLinearElasticityElement element
571 //============================================================================
572 template<>
574 : public virtual QElement<2, 2>
575 {
576 public:
577 /// Constructor must call the constructor of the underlying element
578 FaceGeometry() : QElement<2, 2>() {}
579 };
580
581 //============================================================================
582 /// FaceGeometry of a quadratic 3D QLinearElasticityElement element
583 //============================================================================
584 template<>
586 : public virtual QElement<2, 3>
587 {
588 public:
589 /// Constructor must call the constructor of the underlying element
590 FaceGeometry() : QElement<2, 3>() {}
591 };
592
593
594 //============================================================================
595 /// FaceGeometry of a cubic 3D QLinearElasticityElement element
596 //============================================================================
597 template<>
599 : public virtual QElement<2, 4>
600 {
601 public:
602 /// Constructor must call the constructor of the underlying element
603 FaceGeometry() : QElement<2, 4>() {}
604 };
605
606 /// /////////////////////////////////////////////////////////////////
607 /// /////////////////////////////////////////////////////////////////
608 /// /////////////////////////////////////////////////////////////////
609
610
611 //==========================================================
612 /// Linear elasticity upgraded to become projectable
613 //==========================================================
614 template<class LINEAR_ELAST_ELEMENT>
616 : public virtual ProjectableElement<LINEAR_ELAST_ELEMENT>
617 {
618 public:
619 /// Constructor [this was only required explicitly
620 /// from gcc 4.5.2 onwards...]
622
623
624 /// Specify the values associated with field fld.
625 /// The information is returned in a vector of pairs which comprise
626 /// the Data object and the value within it, that correspond to field fld.
627 /// In the underlying linear elasticity elements the
628 /// the displacements are stored at the nodal values
630 {
631 // Create the vector
633
634 // Loop over all nodes and extract the fld-th nodal value
635 unsigned nnod = this->nnode();
636 for (unsigned j = 0; j < nnod; j++)
637 {
638 // Add the data value associated with the displacement components
639 data_values.push_back(std::make_pair(this->node_pt(j), fld));
640 }
641
642 // Return the vector
643 return data_values;
644 }
645
646 /// Number of fields to be projected: dim, corresponding to
647 /// the displacement components
649 {
650 return this->dim();
651 }
652
653 /// Number of history values to be stored for fld-th field.
654 /// (includes present value!)
655 unsigned nhistory_values_for_projection(const unsigned& fld)
656 {
657#ifdef PARANOID
658 if (fld > 1)
659 {
660 std::stringstream error_stream;
661 error_stream << "Elements only store two fields so fld can't be"
662 << " " << fld << std::endl;
663 throw OomphLibError(
664 error_stream.str(), OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION);
665 }
666#endif
667 return this->node_pt(0)->ntstorage();
668 }
669
670 /// Number of positional history values: Read out from
671 /// positional timestepper (Note: count includes current value!)
673 {
674 return this->node_pt(0)->position_time_stepper_pt()->ntstorage();
675 }
676
677 /// Return Jacobian of mapping and shape functions of field fld
678 /// at local coordinate s
679 double jacobian_and_shape_of_field(const unsigned& fld,
680 const Vector<double>& s,
681 Shape& psi)
682 {
683 unsigned n_dim = this->dim();
684 unsigned n_node = this->nnode();
685 DShape dpsidx(n_node, n_dim);
686
687 // Call the derivatives of the shape functions and return
688 // the Jacobian
689 return this->dshape_eulerian(s, psi, dpsidx);
690 }
691
692
693 /// Return interpolated field fld at local coordinate s, at time
694 /// level t (t=0: present; t>0: history values)
695 double get_field(const unsigned& t,
696 const unsigned& fld,
697 const Vector<double>& s)
698 {
699 unsigned n_node = this->nnode();
700
701 /* #ifdef PARANOID */
702 /* unsigned n_dim=this->node_pt(0)->ndim(); */
703 /* #endif */
704
705 // Local shape function
706 Shape psi(n_node);
707
708 // Find values of shape function
709 this->shape(s, psi);
710
711 // Initialise value of u
712 double interpolated_u = 0.0;
713
714 // Sum over the local nodes at that time
715 for (unsigned l = 0; l < n_node; l++)
716 {
717 // over-zealous I think. This will quietly do the right thing
718 // even if there are additional degrees of freedom floating around.
719 /* #ifdef PARANOID */
720 /* unsigned nvalue=this->node_pt(l)->nvalue(); */
721 /* if (nvalue!=n_dim) */
722 /* { */
723 /* std::stringstream error_stream; */
724 /* error_stream */
725 /* << "Current implementation only works for non-resized
726 * nodes\n" */
727 /* << "but nvalue= " << nvalue << "!= dim = " << n_dim <<
728 * std::endl; */
729 /* throw OomphLibError( */
730 /* error_stream.str(), */
731 /* OOMPH_CURRENT_FUNCTION, */
732 /* OOMPH_EXCEPTION_LOCATION); */
733 /* } */
734 /* #endif */
735 interpolated_u += this->nodal_value(t, l, fld) * psi[l];
736 }
737 return interpolated_u;
738 }
739
740
741 /// Return number of values in field fld
742 unsigned nvalue_of_field(const unsigned& fld)
743 {
744 return this->nnode();
745 }
746
747
748 /// Return local equation number of value j in field fld.
749 int local_equation(const unsigned& fld, const unsigned& j)
750 {
751 // over-zealous I think. This will quietly do the right thing
752 // even if there are additional degrees of freedom floating around.
753 /* #ifdef PARANOID */
754 /* unsigned n_dim=this->node_pt(0)->ndim(); */
755 /* unsigned nvalue=this->node_pt(j)->nvalue(); */
756 /* if (nvalue!=n_dim) */
757 /* { */
758 /* std::stringstream error_stream; */
759 /* error_stream */
760 /* << "Current implementation only works for non-resized nodes\n"
761 */
762 /* << "but nvalue= " << nvalue << "!= dim = " << n_dim <<
763 * std::endl; */
764 /* throw OomphLibError( */
765 /* error_stream.str(), */
766 /* OOMPH_CURRENT_FUNCTION, */
767 /* OOMPH_EXCEPTION_LOCATION); */
768 /* } */
769 /* #endif */
770 return this->nodal_local_eqn(j, fld);
771 }
772 };
773
774
775 //=======================================================================
776 /// Face geometry for element is the same as that for the underlying
777 /// wrapped element
778 //=======================================================================
779 template<class ELEMENT>
781 : public virtual FaceGeometry<ELEMENT>
782 {
783 public:
784 FaceGeometry() : FaceGeometry<ELEMENT>() {}
785 };
786
787
788 //=======================================================================
789 /// Face geometry of the Face Geometry for element is the same as
790 /// that for the underlying wrapped element
791 //=======================================================================
792 template<class ELEMENT>
794 : public virtual FaceGeometry<FaceGeometry<ELEMENT>>
795 {
796 public:
798 };
799
800
801} // namespace oomph
802
803#endif
e
Definition: cfortran.h:571
static char t char * s
Definition: cfortran.h:568
cstr elem_len * i
Definition: cfortran.h:603
char t
Definition: cfortran.h:568
A Class for the derivatives of shape functions The class design is essentially the same as Shape,...
Definition: shape.h:278
TimeStepper *& time_stepper_pt()
Return the pointer to the timestepper.
Definition: nodes.h:238
unsigned ntstorage() const
Return total number of doubles stored per value to record time history of each value (one for steady ...
Definition: nodes.cc:879
A base class that represents the fourth-rank elasticity tensor defined such that.
FaceGeometry()
Constructor must call the constructor of the underlying solid element.
FaceGeometry()
Constructor must call the constructor of the underlying element.
FaceGeometry()
Constructor must call the constructor of the underlying element.
FaceGeometry()
Constructor must call the constructor of the underlying element.
FaceGeometry()
Constructor must call the constructor of the underlying element.
FaceGeometry()
Constructor must call the constructor of the underlying element.
//////////////////////////////////////////////////////////////////// ////////////////////////////////...
Definition: elements.h:4998
A general Finite Element class.
Definition: elements.h:1313
double nodal_value(const unsigned &n, const unsigned &i) const
Return the i-th value stored at local node n. Produces suitably interpolated values for hanging nodes...
Definition: elements.h:2593
virtual void shape(const Vector< double > &s, Shape &psi) const =0
Calculate the geometric shape functions at local coordinate s. This function must be overloaded for e...
int nodal_local_eqn(const unsigned &n, const unsigned &i) const
Return the local equation number corresponding to the i-th value at the n-th local node.
Definition: elements.h:1432
unsigned dim() const
Return the spatial dimension of the element, i.e. the number of local coordinates required to paramet...
Definition: elements.h:2611
unsigned nnode() const
Return the number of nodes.
Definition: elements.h:2210
void(* SteadyExactSolutionFctPt)(const Vector< double > &, Vector< double > &)
Function pointer for function that computes vector-valued steady "exact solution" as .
Definition: elements.h:1759
double dshape_eulerian(const Vector< double > &s, Shape &psi, DShape &dpsidx) const
Compute the geometric shape functions and also first derivatives w.r.t. global coordinates at local c...
Definition: elements.cc:3298
Node *& node_pt(const unsigned &n)
Return a pointer to the local node n.
Definition: elements.h:2175
void(* UnsteadyExactSolutionFctPt)(const double &, const Vector< double > &, Vector< double > &)
Function pointer for function that computes Vector-valued time-dependent function as .
Definition: elements.h:1765
unsigned long eqn_number(const unsigned &ieqn_local) const
Return the global equation number corresponding to the ieqn_local-th local equation number.
Definition: elements.h:704
static DenseMatrix< double > Dummy_matrix
Empty dense matrix used as a dummy argument to combined residual and jacobian functions in the case w...
Definition: elements.h:227
TimeStepper *& time_stepper_pt()
Access function for pointer to time stepper: Null if object is not time-dependent.
Definition: geom_objects.h:192
A base class for elements that solve the equations of linear elasticity in Cartesian coordinates....
ElasticityTensor * Elasticity_tensor_pt
Pointer to the elasticity tensor.
BodyForceFctPt body_force_fct_pt() const
Access function: Pointer to body force function (const version)
void disable_inertia()
Switch off solid inertia.
static double Default_lambda_sq_value
Static default value for timescale ratio (1.0 – for natural scaling)
double * Lambda_sq_pt
Timescale ratio (non-dim. density)
virtual void get_stress(const Vector< double > &s, DenseMatrix< double > &sigma) const =0
Return the Cauchy stress tensor, as calculated from the elasticity tensor at specified local coordina...
ElasticityTensor *& elasticity_tensor_pt()
Return the pointer to the elasticity_tensor.
static void pin_redundant_nodal_solid_pressures(const Vector< GeneralisedElement * > &element_pt)
Loop over all elements in Vector (which typically contains all the elements in a refineable solid mes...
double interpolated_u_linear_elasticity(const Vector< double > &s, const unsigned &i) const
Return FE interpolated displacement u[i] at local coordinate s.
BodyForceFctPt & body_force_fct_pt()
Access function: Pointer to body force function.
double *& lambda_sq_pt()
Access function for pointer to timescale ratio (nondim density)
void get_strain(const Vector< double > &s, DenseMatrix< double > &strain) const
Return the strain tensor.
void body_force(const Vector< double > &x, Vector< double > &b) const
Evaluate body force at Eulerian coordinate x at present time (returns zero vector if no body force fu...
BodyForceFctPt Body_force_fct_pt
Pointer to body force function.
bool is_inertia_enabled() const
Access function to flag that switches inertia on/off (const version)
void interpolated_u_linear_elasticity(const Vector< double > &s, Vector< double > &disp) const
Compute vector of FE interpolated displacement u at local coordinate s.
bool Unsteady
Flag that switches inertia on/off.
const double & lambda_sq() const
Access function for timescale ratio (nondim density)
void enable_inertia()
Switch on solid inertia.
double E(const unsigned &i, const unsigned &j, const unsigned &k, const unsigned &l) const
Access function to the entries in the elasticity tensor.
void(* BodyForceFctPt)(const double &t, const Vector< double > &x, Vector< double > &b)
Function pointer to function that specifies the body force as a function of the Cartesian coordinates...
void get_dof_numbers_for_unknowns(std::list< std::pair< unsigned long, unsigned > > &dof_lookup_list) const
Create a list of pairs for all unknowns in this element, so that the first entry in each pair contain...
unsigned ndof_types() const
The number of "DOF types" that degrees of freedom in this element are sub-divided into: for now lump ...
virtual unsigned u_index_linear_elasticity(const unsigned i) const
Return the index at which the i-th unknown displacement component is stored. The default value,...
virtual void pin_elemental_redundant_nodal_solid_pressures()
Pin the element's redundant solid pressures (needed for refinement)
double d2u_dt2_linear_elasticity(const unsigned &n, const unsigned &i) const
d^2u/dt^2 at local node n
LinearElasticityEquationsBase()
Constructor: Set null pointers for constitutive law and for isotropic growth function....
//////////////////////////////////////////////////////////////////// ////////////////////////////////...
void output_fct(std::ostream &outfile, const unsigned &nplot, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt)
Output exact solution x,y,[z],u,v,[w].
unsigned required_nvalue(const unsigned &n) const
Number of values required at node n.
void compute_error(std::ostream &outfile, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt, double &error, double &norm)
Validate against exact solution. Solution is provided via function pointer. Plot at a given number of...
void output(std::ostream &outfile)
Output: x,y,[z],u,v,[w].
void output(FILE *file_pt)
C-style output: x,y,[z],u,v,[w].
void get_stress(const Vector< double > &s, DenseMatrix< double > &sigma) const
Return the Cauchy stress tensor, as calculated from the elasticity tensor at specified local coordina...
virtual void fill_in_generic_contribution_to_residuals_linear_elasticity(Vector< double > &residuals, DenseMatrix< double > &jacobian, unsigned flag)
Private helper function to compute residuals and (if requested via flag) also the Jacobian matrix.
void fill_in_contribution_to_jacobian(Vector< double > &residuals, DenseMatrix< double > &jacobian)
The jacobian is calculated by finite differences by default, We need only to take finite differences ...
void fill_in_contribution_to_residuals(Vector< double > &residuals)
Return the residuals for the solid equations (the discretised principle of virtual displacements)
TimeStepper *& position_time_stepper_pt()
Return a pointer to the position timestepper.
Definition: nodes.h:1022
An OomphLibError object which should be thrown when an run-time error is encountered....
Wrapper class for projectable elements. Adds "projectability" to the underlying ELEMENT.
Definition: projection.h:183
///////////////////////////////////////////////////////////////// ///////////////////////////////////...
double get_field(const unsigned &t, const unsigned &fld, const Vector< double > &s)
Return interpolated field fld at local coordinate s, at time level t (t=0: present; t>0: history valu...
int local_equation(const unsigned &fld, const unsigned &j)
Return local equation number of value j in field fld.
double jacobian_and_shape_of_field(const unsigned &fld, const Vector< double > &s, Shape &psi)
Return Jacobian of mapping and shape functions of field fld at local coordinate s.
unsigned nhistory_values_for_projection(const unsigned &fld)
Number of history values to be stored for fld-th field. (includes present value!)
unsigned nvalue_of_field(const unsigned &fld)
Return number of values in field fld.
ProjectableLinearElasticityElement()
Constructor [this was only required explicitly from gcc 4.5.2 onwards...].
Vector< std::pair< Data *, unsigned > > data_values_of_field(const unsigned &fld)
Specify the values associated with field fld. The information is returned in a vector of pairs which ...
unsigned nfields_for_projection()
Number of fields to be projected: dim, corresponding to the displacement components.
unsigned nhistory_values_for_coordinate_projection()
Number of positional history values: Read out from positional timestepper (Note: count includes curre...
/////////////////////////////////////////////////////////////////////// /////////////////////////////...
Definition: Qelements.h:459
///////////////////////////////////////////////////////////////////// ///////////////////////////////...
void output(FILE *file_pt, const unsigned &n_plot)
C-style output function.
void output(FILE *file_pt)
C-style output function.
void output(std::ostream &outfile, const unsigned &n_plot)
Output function.
void output_fct(std::ostream &outfile, const unsigned &nplot, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt)
Output exact solution x,y,[z],u,v,[w].
void output(std::ostream &outfile)
Output function.
A Class for shape functions. In simple cases, the shape functions have only one index that can be tho...
Definition: shape.h:76
////////////////////////////////////////////////////////////////////// //////////////////////////////...
Definition: timesteppers.h:231
unsigned ntstorage() const
Return the number of doubles required to represent history (one for steady)
Definition: timesteppers.h:601
virtual double weight(const unsigned &i, const unsigned &j) const
Access function for j-th weight for the i-th derivative.
Definition: timesteppers.h:594
Time *const & time_pt() const
Access function for the pointer to time (const version)
Definition: timesteppers.h:572
bool is_steady() const
Flag to indicate if a timestepper has been made steady (possibly temporarily to switch off time-depen...
Definition: timesteppers.h:389
double & time()
Return the current value of the continuous time.
Definition: timesteppers.h:123
//////////////////////////////////////////////////////////////////// ////////////////////////////////...