time_harmonic_fourier_decomposed_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
27// Include guards to prevent multiple inclusion of the header
28#ifndef OOMPH_FOURIER_DECOMPOSED_TIME_HARMONIC_LINEAR_ELASTICITY_ELEMENTS_HEADER
29#define OOMPH_FOURIER_DECOMPOSED_TIME_HARMONIC_LINEAR_ELASTICITY_ELEMENTS_HEADER
30
31// Config header generated by autoconfig
32#ifdef HAVE_CONFIG_H
33#include <oomph-lib-config.h>
34#endif
35
36
37#ifdef OOMPH_HAS_MPI
38#include "mpi.h"
39#endif
40
41#include <complex>
42
43
44// OOMPH-LIB headers
45#include "../generic/Qelements.h"
46#include "../generic/Telements.h"
47#include "../generic/projection.h"
48#include "../generic/error_estimator.h"
49
50
51namespace oomph
52{
53 //=======================================================================
54 /// A base class for elements that solve the Fourier decomposed (in
55 /// cylindrical polars) equations of time-harmonic linear elasticity.
56 //=======================================================================
58 : public virtual FiniteElement
59 {
60 public:
61 /// Return the index at which the i-th (i=0: r, i=1: z; i=2: theta)
62 /// real or imag unknown displacement component is stored at the nodes.
63 /// The default assignment here (u_r_real, u_z_real, ..., u_theta_imag)
64 /// is appropriate for single-physics problems.
66 const unsigned i) const
67 {
68 return std::complex<unsigned>(i, i + 3);
69 }
70
71 /// Compute vector of FE interpolated displacement u at local coordinate s
73 const Vector<double>& s, Vector<std::complex<double>>& disp) const
74 {
75 // Find number of nodes
76 unsigned n_node = nnode();
77
78 // Local shape function
79 Shape psi(n_node);
80
81 // Find values of shape function
82 shape(s, psi);
83
84 for (unsigned i = 0; i < 3; i++)
85 {
86 // Index at which the nodal value is stored
87 std::complex<unsigned> u_nodal_index =
89
90 // Initialise value of u
91 disp[i] = std::complex<double>(0.0, 0.0);
92
93 // Loop over the local nodes and sum
94 for (unsigned l = 0; l < n_node; l++)
95 {
96 const std::complex<double> u_value(
97 nodal_value(l, u_nodal_index.real()),
98 nodal_value(l, u_nodal_index.imag()));
99
100 disp[i] += u_value * psi[l];
101 }
102 }
103 }
104
105 /// Return FE interpolated displacement u[i] (i=0: r, i=1: z; i=2:
106 /// theta) at local coordinate s
108 const Vector<double>& s, const unsigned& i) const
109 {
110 // Find number of nodes
111 unsigned n_node = nnode();
112
113 // Local shape function
114 Shape psi(n_node);
115
116 // Find values of shape function
117 shape(s, psi);
118
119 // Get nodal index at which i-th velocity is stored
120 std::complex<unsigned> u_nodal_index =
122
123 // Initialise value of u
124 std::complex<double> interpolated_u(0.0, 0.0);
125
126 // Loop over the local nodes and sum
127 for (unsigned l = 0; l < n_node; l++)
128 {
129 const std::complex<double> u_value(
130 nodal_value(l, u_nodal_index.real()),
131 nodal_value(l, u_nodal_index.imag()));
132
133 interpolated_u += u_value * psi[l];
134 }
135
136 return (interpolated_u);
137 }
138
139
140 /// Function pointer to function that specifies the body force
141 /// as a function of the Cartesian coordinates and time FCT(x,b) --
142 /// x and b are Vectors!
143 typedef void (*BodyForceFctPt)(const Vector<double>& x,
145
146 /// Constructor: Set null pointers for constitutive law.
147 /// Set physical parameter values to
148 /// default values, and set body force to zero.
152 Nu_pt(0),
155 {
156 }
157
158 /// Access function for square of non-dim frequency
159 const std::complex<double>& omega_sq() const
160 {
161 return *Omega_sq_pt;
162 }
163
164 /// Access function for square of non-dim frequency
165 std::complex<double>*& omega_sq_pt()
166 {
167 return Omega_sq_pt;
168 }
169
170 /// Return the pointer to Young's modulus
171 std::complex<double>*& youngs_modulus_pt()
172 {
173 return Youngs_modulus_pt;
174 }
175
176 /// Access function to Young's modulus
177 inline std::complex<double> youngs_modulus() const
178 {
179 return (*Youngs_modulus_pt);
180 }
181
182 /// Access function for Poisson's ratio
183 std::complex<double>& nu() const
184 {
185#ifdef PARANOID
186 if (Nu_pt == 0)
187 {
188 std::ostringstream error_message;
189 error_message << "No pointer to Poisson's ratio set. Please set one!\n";
190 throw OomphLibError(error_message.str(),
191 OOMPH_CURRENT_FUNCTION,
192 OOMPH_EXCEPTION_LOCATION);
193 }
194#endif
195 return *Nu_pt;
196 }
197
198 /// Access function for pointer to Poisson's ratio
199 std::complex<double>*& nu_pt()
200 {
201 return Nu_pt;
202 }
203
204 /// Access function for Fourier wavenumber
206 {
207#ifdef PARANOID
208 if (Fourier_wavenumber_pt == 0)
209 {
210 std::ostringstream error_message;
211 error_message
212 << "No pointer to Fourier wavenumber set. Please set one!\n";
213 throw OomphLibError(error_message.str(),
214 OOMPH_CURRENT_FUNCTION,
215 OOMPH_EXCEPTION_LOCATION);
216 }
217#endif
218 return *Fourier_wavenumber_pt;
219 }
220
221 /// Access function for pointer to Fourier wavenumber
223 {
225 }
226
227 /// Access function: Pointer to body force function
229 {
230 return Body_force_fct_pt;
231 }
232
233 /// Access function: Pointer to body force function (const version)
235 {
236 return Body_force_fct_pt;
237 }
238
239 /// Evaluate body force at Eulerian coordinate x at present time
240 /// (returns zero vector if no body force function pointer has been set)
241 inline void body_force(const Vector<double>& x,
242 Vector<std::complex<double>>& b) const
243 {
244 // If no function has been set, return zero vector
245 if (Body_force_fct_pt == 0)
246 {
247 // Get spatial dimension of element
248 unsigned n = dim();
249 for (unsigned i = 0; i < n; i++)
250 {
251 b[i] = std::complex<double>(0.0, 0.0);
252 }
253 }
254 else
255 {
256 (*Body_force_fct_pt)(x, b);
257 }
258 }
259
260 /// The number of "DOF types" that degrees of freedom in this element
261 /// are sub-divided into: for now lump them all into one DOF type.
262 /// Can be adjusted later
263 unsigned ndof_types() const
264 {
265 return 1;
266 }
267
268 /// Create a list of pairs for all unknowns in this element,
269 /// so that the first entry in each pair contains the global equation
270 /// number of the unknown, while the second one contains the number
271 /// of the "DOF type" that this unknown is associated with.
272 /// (Function can obviously only be called if the equation numbering
273 /// scheme has been set up.)
275 std::list<std::pair<unsigned long, unsigned>>& dof_lookup_list) const
276 {
277 // temporary pair (used to store dof lookup prior to being added
278 // to list)
279 std::pair<unsigned long, unsigned> dof_lookup;
280
281 // number of nodes
282 const unsigned n_node = this->nnode();
283
284 // Integer storage for local unknown
285 int local_unknown = 0;
286
287 // Loop over the nodes
288 for (unsigned n = 0; n < n_node; n++)
289 {
290 // Loop over dimension (real and imag displacement components)
291 for (unsigned i = 0; i < 6; i++)
292 {
293 // If the variable is free
294 local_unknown = nodal_local_eqn(n, i);
295
296 // ignore pinned values
297 if (local_unknown >= 0)
298 {
299 // store dof lookup in temporary pair: First entry in pair
300 // is global equation number; second entry is DOF type
301 dof_lookup.first = this->eqn_number(local_unknown);
302 dof_lookup.second = 0;
303
304 // add to list
305 dof_lookup_list.push_front(dof_lookup);
306 }
307 }
308 }
309 }
310
311
312 protected:
313 /// Square of nondim frequency
314 std::complex<double>* Omega_sq_pt;
315
316 /// Pointer to the Young's modulus
317 std::complex<double>* Youngs_modulus_pt;
318
319 /// Pointer to Poisson's ratio
320 std::complex<double>* Nu_pt;
321
322 /// Pointer to Fourier wavenumber
324
325 /// Pointer to body force function
327
328 /// Static default value for squared frequency
329 static std::complex<double> Default_omega_sq_value;
330
331 /// Static default value for Young's modulus (1.0 -- for natural
332 /// scaling, i.e. all stresses have been non-dimensionalised by
333 /// the same reference Young's modulus. Setting the "non-dimensional"
334 /// Young's modulus (obtained by de-referencing Youngs_modulus_pt)
335 /// to a number larger than one means that the material is stiffer
336 /// than assumed in the non-dimensionalisation.
337 static std::complex<double> Default_youngs_modulus_value;
338 };
339
340
341 /// ////////////////////////////////////////////////////////////////////
342 /// ////////////////////////////////////////////////////////////////////
343 /// ////////////////////////////////////////////////////////////////////
344
345
346 //=======================================================================
347 /// A class for elements that solve the Fourier decomposed (in cylindrical
348 /// polars) equations of time-harmonic linear elasticity
349 //=======================================================================
352 {
353 public:
354 /// Constructor
356
357 /// Number of values required at node n.
358 unsigned required_nvalue(const unsigned& n) const
359 {
360 return 6;
361 }
362
363 /// Return the residuals for the equations (the discretised
364 /// principle of virtual displacements)
366 {
369 }
370
371
372 /// The jacobian is calculated by finite differences by default,
373 /// We need only to take finite differences w.r.t. positional variables
374 /// For this element
376 DenseMatrix<double>& jacobian)
377 {
378 // Add the contribution to the residuals
379 this
380 ->fill_in_generic_contribution_to_residuals_fourier_decomp_time_harmonic_linear_elasticity(
381 residuals, jacobian, 1);
382 }
383
384
385 /// Get strain tensor
386 void get_strain(const Vector<double>& s,
387 DenseMatrix<std::complex<double>>& strain);
388
389 /// Compute norm of solution: square of the L2 norm
390 void compute_norm(double& norm);
391
392 /// Output exact solution: r,z, u_r_real, u_z_real, ..., u_theta_imag
393 void output_fct(std::ostream& outfile,
394 const unsigned& nplot,
396
397 /// Output: r,z, u_r_real, u_z_real, ..., u_theta_imag
398 void output(std::ostream& outfile)
399 {
400 unsigned n_plot = 5;
401 output(outfile, n_plot);
402 }
403
404 /// Output: r,z, u_r_real, u_z_real, ..., u_theta_imag
405 void output(std::ostream& outfile, const unsigned& n_plot);
406
407 /// C-style output: r,z, u_r_real, u_z_real, ..., u_theta_imag
408 void output(FILE* file_pt)
409 {
410 unsigned n_plot = 5;
411 output(file_pt, n_plot);
412 }
413
414 /// Output: r,z, u_r_real, u_z_real, ..., u_theta_imag
415 void output(FILE* file_pt, const unsigned& n_plot);
416
417 /// Validate against exact solution.
418 /// Solution is provided via function pointer.
419 /// Plot at a given number of plot points and compute L2 error
420 /// and L2 norm of displacement solution over element
421 void compute_error(std::ostream& outfile,
423 double& error,
424 double& norm);
425
426
427 private:
428 /// Private helper function to compute residuals and (if requested
429 /// via flag) also the Jacobian matrix.
431 Vector<double>& residuals, DenseMatrix<double>& jacobian, unsigned flag);
432 };
433
434
435 /// /////////////////////////////////////////////////////////////////////
436 /// /////////////////////////////////////////////////////////////////////
437 /// /////////////////////////////////////////////////////////////////////
438
439
440 //===========================================================================
441 /// An Element that solves the equations of Fourier decomposed (in cylindrical
442 /// polars) time-harmonic linear elasticity, using QElements for the geometry.
443 //============================================================================
444 template<unsigned NNODE_1D>
446 : public virtual QElement<2, NNODE_1D>,
448 {
449 public:
450 /// Constructor
452 : QElement<2, NNODE_1D>(),
454 {
455 }
456
457 /// Output function
458 void output(std::ostream& outfile)
459 {
461 }
462
463 /// Output function
464 void output(std::ostream& outfile, const unsigned& n_plot)
465 {
467 n_plot);
468 }
469
470
471 /// C-style output function
472 void output(FILE* file_pt)
473 {
475 }
476
477 /// C-style output function
478 void output(FILE* file_pt, const unsigned& n_plot)
479 {
481 n_plot);
482 }
483 };
484
485
486 //============================================================================
487 /// FaceGeometry of a linear
488 /// QTimeHarmonicFourierDecomposedLinearElasticityElement element
489 //============================================================================
490 template<unsigned NNODE_1D>
493 : public virtual QElement<1, NNODE_1D>
494 {
495 public:
496 /// Constructor must call the constructor of the underlying element
497 FaceGeometry() : QElement<1, NNODE_1D>() {}
498 };
499
500
501 /// /////////////////////////////////////////////////////////////////////
502 /// /////////////////////////////////////////////////////////////////////
503 /// /////////////////////////////////////////////////////////////////////
504
505
506 //===========================================================================
507 /// An Element that solves the equations of Fourier decomposed (in cylindrical
508 /// polars) time-harmonic linear elasticity, using TElements for the geometry.
509 //============================================================================
510 template<unsigned NNODE_1D>
512 : public virtual TElement<2, NNODE_1D>,
514 public virtual ElementWithZ2ErrorEstimator
515 {
516 public:
517 /// Constructor
519 : TElement<2, NNODE_1D>(),
521 {
522 }
523
524 /// Output function
525 void output(std::ostream& outfile)
526 {
528 }
529
530 /// Output function
531 void output(std::ostream& outfile, const unsigned& n_plot)
532 {
534 n_plot);
535 }
536
537 /// C-style output function
538 void output(FILE* file_pt)
539 {
541 }
542
543 /// C-style output function
544 void output(FILE* file_pt, const unsigned& n_plot)
545 {
547 n_plot);
548 }
549
550
551 /// Number of vertex nodes in the element
552 unsigned nvertex_node() const
553 {
555 }
556
557 /// Pointer to the j-th vertex node in the element
558 Node* vertex_node_pt(const unsigned& j) const
559 {
561 }
562
563 /// Order of recovery shape functions for Z2 error estimation:
564 /// Same order as shape functions.
566 {
567 return NNODE_1D - 1;
568 }
569
570 /// Number of 'flux' terms for Z2 error estimation
572 {
573 // 3 Diagonal strain rates and 3 off diagonal terms for real and imag part
574 return 12;
575 }
576
577 /// Get 'flux' for Z2 error recovery: Upper triangular entries
578 /// in strain tensor.
580 {
581#ifdef PARANOID
582 unsigned num_entries = 12;
583 if (flux.size() != num_entries)
584 {
585 std::ostringstream error_message;
586 error_message << "The flux vector has the wrong number of entries, "
587 << flux.size() << ", whereas it should be " << num_entries
588 << std::endl;
589 throw OomphLibError(error_message.str(),
590 OOMPH_CURRENT_FUNCTION,
591 OOMPH_EXCEPTION_LOCATION);
592 }
593#endif
594
595 // Get strain matrix
597 this->get_strain(s, strain);
598
599 // Pack into flux Vector
600 unsigned icount = 0;
601
602 // Start with diagonal terms
603 for (unsigned i = 0; i < 3; i++)
604 {
605 flux[icount] = strain(i, i).real();
606 icount++;
607 flux[icount] = strain(i, i).imag();
608 icount++;
609 }
610
611 // Off diagonals row by row
612 for (unsigned i = 0; i < 3; i++)
613 {
614 for (unsigned j = i + 1; j < 3; j++)
615 {
616 flux[icount] = strain(i, j).real();
617 icount++;
618 flux[icount] = strain(i, j).imag();
619 icount++;
620 }
621 }
622 }
623 };
624
625
626 //============================================================================
627 /// FaceGeometry of a linear
628 /// TTimeHarmonicFourierDecomposedLinearElasticityElement element
629 //============================================================================
630 template<unsigned NNODE_1D>
633 : public virtual TElement<1, NNODE_1D>
634 {
635 public:
636 /// Constructor must call the constructor of the underlying element
637 FaceGeometry() : TElement<1, NNODE_1D>() {}
638 };
639
640
641 /// /////////////////////////////////////////////////////////////////
642 /// /////////////////////////////////////////////////////////////////
643 /// /////////////////////////////////////////////////////////////////
644
645
646 //==========================================================
647 /// Fourier-decomposed time-harmonic linear elasticity
648 /// upgraded to become projectable
649 //==========================================================
650 template<class TIME_HARMONIC_LINEAR_ELAST_ELEMENT>
652 : public virtual ProjectableElement<TIME_HARMONIC_LINEAR_ELAST_ELEMENT>
653 {
654 public:
655 /// Constructor [this was only required explicitly
656 /// from gcc 4.5.2 onwards...]
658
659 /// Specify the values associated with field fld.
660 /// The information is returned in a vector of pairs which comprise
661 /// the Data object and the value within it, that correspond to field fld.
662 /// In the underlying time-harmonic linear elasticity elemements the
663 /// real and complex parts of the displacements are stored
664 /// at the nodal values
666 {
667 // Create the vector
669
670 // Loop over all nodes and extract the fld-th nodal value
671 unsigned nnod = this->nnode();
672 for (unsigned j = 0; j < nnod; j++)
673 {
674 // Add the data value associated with the velocity components
675 data_values.push_back(std::make_pair(this->node_pt(j), fld));
676 }
677
678 // Return the vector
679 return data_values;
680 }
681
682 /// Number of fields to be projected: 3*dim, corresponding to
683 /// real and imag parts of the displacement components
685 {
686 return 3 * this->dim();
687 }
688
689 /// Number of history values to be stored for fld-th field.
690 /// (includes present value!)
691 unsigned nhistory_values_for_projection(const unsigned& fld)
692 {
693#ifdef PARANOID
694 if (fld > 5)
695 {
696 std::stringstream error_stream;
697 error_stream << "Elements only store six fields so fld can't be"
698 << " " << fld << std::endl;
699 throw OomphLibError(
700 error_stream.str(), OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION);
701 }
702#endif
703 return this->node_pt(0)->ntstorage();
704 }
705
706 /// Number of positional history values: Read out from
707 /// positional timestepper
708 /// (Note: count includes current value!)
710 {
711 return this->node_pt(0)->position_time_stepper_pt()->ntstorage();
712 }
713
714 /// Return Jacobian of mapping and shape functions of field fld
715 /// at local coordinate s
716 double jacobian_and_shape_of_field(const unsigned& fld,
717 const Vector<double>& s,
718 Shape& psi)
719 {
720 unsigned n_dim = this->dim();
721 unsigned n_node = this->nnode();
722 DShape dpsidx(n_node, n_dim);
723
724 // Call the derivatives of the shape functions and return
725 // the Jacobian
726 return this->dshape_eulerian(s, psi, dpsidx);
727 }
728
729
730 /// Return interpolated field fld at local coordinate s, at time
731 /// level t (t=0: present; t>0: history values)
732 double get_field(const unsigned& t,
733 const unsigned& fld,
734 const Vector<double>& s)
735 {
736 unsigned n_node = this->nnode();
737
738#ifdef PARANOID
739 unsigned n_dim = this->node_pt(0)->ndim();
740#endif
741
742 // Local shape function
743 Shape psi(n_node);
744
745 // Find values of shape function
746 this->shape(s, psi);
747
748 // Initialise value of u
749 double interpolated_u = 0.0;
750
751 // Sum over the local nodes at that time
752 for (unsigned l = 0; l < n_node; l++)
753 {
754#ifdef PARANOID
755 unsigned nvalue = this->node_pt(l)->nvalue();
756 if (nvalue != 3 * n_dim)
757 {
758 std::stringstream error_stream;
759 error_stream
760 << "Current implementation only works for non-resized nodes\n"
761 << "but nvalue= " << nvalue << "!= 3 dim = " << 3 * n_dim
762 << std::endl;
763 throw OomphLibError(error_stream.str(),
764 OOMPH_CURRENT_FUNCTION,
765 OOMPH_EXCEPTION_LOCATION);
766 }
767#endif
768 interpolated_u += this->nodal_value(t, l, fld) * psi[l];
769 }
770 return interpolated_u;
771 }
772
773
774 /// Return number of values in field fld
775 unsigned nvalue_of_field(const unsigned& fld)
776 {
777 return this->nnode();
778 }
779
780
781 /// Return local equation number of value j in field fld.
782 int local_equation(const unsigned& fld, const unsigned& j)
783 {
784#ifdef PARANOID
785 unsigned n_dim = this->node_pt(0)->ndim();
786 unsigned nvalue = this->node_pt(j)->nvalue();
787 if (nvalue != 3 * n_dim)
788 {
789 std::stringstream error_stream;
790 error_stream
791 << "Current implementation only works for non-resized nodes\n"
792 << "but nvalue= " << nvalue << "!= 3 dim = " << 3 * n_dim
793 << std::endl;
794 throw OomphLibError(
795 error_stream.str(), OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION);
796 }
797#endif
798 return this->nodal_local_eqn(j, fld);
799 }
800 };
801
802
803 //=======================================================================
804 /// Face geometry for element is the same as that for the underlying
805 /// wrapped element
806 //=======================================================================
807 template<class ELEMENT>
810 : public virtual FaceGeometry<ELEMENT>
811 {
812 public:
813 FaceGeometry() : FaceGeometry<ELEMENT>() {}
814 };
815
816
817 //=======================================================================
818 /// Face geometry of the Face Geometry for element is the same as
819 /// that for the underlying wrapped element
820 //=======================================================================
821 template<class ELEMENT>
824 : public virtual FaceGeometry<FaceGeometry<ELEMENT>>
825 {
826 public:
828 };
829
830
831} // namespace oomph
832
833
834#endif
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
unsigned nvalue() const
Return number of values stored in data object (incl pinned ones).
Definition: nodes.h:483
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
Base class for finite elements that can compute the quantities that are required for the Z2 error est...
//////////////////////////////////////////////////////////////////// ////////////////////////////////...
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
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
Nodes are derived from Data, but, in addition, have a definite (Eulerian) position in a space of a gi...
Definition: nodes.h:906
TimeStepper *& position_time_stepper_pt()
Return a pointer to the position timestepper.
Definition: nodes.h:1022
unsigned ndim() const
Return (Eulerian) spatial dimension of the node.
Definition: nodes.h:1054
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
///////////////////////////////////////////////////////////////// ///////////////////////////////////...
unsigned nhistory_values_for_projection(const unsigned &fld)
Number of history values to be stored for fld-th field. (includes present value!)
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.
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...
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 nhistory_values_for_coordinate_projection()
Number of positional history values: Read out from positional timestepper (Note: count includes curre...
int local_equation(const unsigned &fld, const unsigned &j)
Return local equation number of value j in field fld.
unsigned nfields_for_projection()
Number of fields to be projected: 3*dim, corresponding to real and imag parts of the displacement com...
ProjectableTimeHarmonicFourierDecomposedLinearElasticityElement()
Constructor [this was only required explicitly from gcc 4.5.2 onwards...].
/////////////////////////////////////////////////////////////////////// /////////////////////////////...
Definition: Qelements.h:459
///////////////////////////////////////////////////////////////////// ///////////////////////////////...
A Class for shape functions. In simple cases, the shape functions have only one index that can be tho...
Definition: shape.h:76
General TElement class.
Definition: Telements.h:1208
///////////////////////////////////////////////////////////////////// ///////////////////////////////...
unsigned nrecovery_order()
Order of recovery shape functions for Z2 error estimation: Same order as shape functions.
void get_Z2_flux(const Vector< double > &s, Vector< double > &flux)
Get 'flux' for Z2 error recovery: Upper triangular entries in strain tensor.
Node * vertex_node_pt(const unsigned &j) const
Pointer to the j-th vertex node in the element.
A base class for elements that solve the Fourier decomposed (in cylindrical polars) equations of time...
void interpolated_u_time_harmonic_fourier_decomposed_linear_elasticity(const Vector< double > &s, Vector< std::complex< double > > &disp) const
Compute vector of FE interpolated displacement u at local coordinate s.
std::complex< double > *& omega_sq_pt()
Access function for square of non-dim frequency.
const std::complex< double > & omega_sq() const
Access function for square of non-dim frequency.
static std::complex< double > Default_youngs_modulus_value
Static default value for Young's modulus (1.0 – for natural scaling, i.e. all stresses have been non-...
static std::complex< double > Default_omega_sq_value
Static default value for squared frequency.
void body_force(const Vector< double > &x, Vector< std::complex< double > > &b) const
Evaluate body force at Eulerian coordinate x at present time (returns zero vector if no body force fu...
unsigned ndof_types() const
The number of "DOF types" that degrees of freedom in this element are sub-divided into: for now lump ...
virtual std::complex< unsigned > u_index_time_harmonic_fourier_decomposed_linear_elasticity(const unsigned i) const
Return the index at which the i-th (i=0: r, i=1: z; i=2: theta) real or imag unknown displacement com...
TimeHarmonicFourierDecomposedLinearElasticityEquationsBase()
Constructor: Set null pointers for constitutive law. Set physical parameter values to default values,...
std::complex< double > *& nu_pt()
Access function for pointer to Poisson's ratio.
void(* BodyForceFctPt)(const Vector< double > &x, Vector< std::complex< double > > &b)
Function pointer to function that specifies the body force as a function of the Cartesian coordinates...
std::complex< double > interpolated_u_time_harmonic_fourier_decomposed_linear_elasticity(const Vector< double > &s, const unsigned &i) const
Return FE interpolated displacement u[i] (i=0: r, i=1: z; i=2: theta) at local coordinate s.
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...
BodyForceFctPt body_force_fct_pt() const
Access function: Pointer to body force function (const version)
//////////////////////////////////////////////////////////////////// ////////////////////////////////...
void output(std::ostream &outfile)
Output: r,z, u_r_real, u_z_real, ..., u_theta_imag.
void fill_in_contribution_to_residuals(Vector< double > &residuals)
Return the residuals for the equations (the discretised principle of virtual displacements)
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_fct(std::ostream &outfile, const unsigned &nplot, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt)
Output exact solution: r,z, u_r_real, u_z_real, ..., u_theta_imag.
unsigned required_nvalue(const unsigned &n) const
Number of values required at node n.
virtual void fill_in_generic_contribution_to_residuals_fourier_decomp_time_harmonic_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 output(FILE *file_pt)
C-style output: r,z, u_r_real, u_z_real, ..., u_theta_imag.
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 get_strain(const Vector< double > &s, DenseMatrix< std::complex< double > > &strain)
Get strain tensor.
unsigned ntstorage() const
Return the number of doubles required to represent history (one for steady)
Definition: timesteppers.h:601
//////////////////////////////////////////////////////////////////// ////////////////////////////////...