biharmonic_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-2023 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 #ifndef OOMPH_BIHARMONIC_ELEMENTS_HEADER
27 #define OOMPH_BIHARMONIC_ELEMENTS_HEADER
28 
29 // Config header generated by autoconfig
30 #ifdef HAVE_CONFIG_H
31 #include <oomph-lib-config.h>
32 #endif
33 
34 #ifdef OOMPH_HAS_MPI
35 // mpi headers
36 #include "mpi.h"
37 #endif
38 
39 // Generic C++ headers
40 #include <iostream>
41 #include <math.h>
42 
43 // oomph-lib headers
44 #include "../generic/matrices.h"
45 #include "../generic/elements.h"
46 #include "../generic/hermite_elements.h"
47 
48 
49 namespace oomph
50 {
51  //=============================================================================
52  /// Biharmonic Equation Class - contains the equations
53  //=============================================================================
54  template<unsigned DIM>
55  class BiharmonicEquations : public virtual FiniteElement
56  {
57  public:
58  /// source function type definition
59  typedef void (*SourceFctPt)(const Vector<double>& x, double& u);
60 
61 
62  /// Constructor (must initialise the Source_fct_pt to null)
64 
65 
67 
68  /// Access function: Nodal function value at local node n
69  /// Uses suitably interpolated value for hanging nodes.
70  virtual double u(const unsigned& n, const unsigned& k) const = 0;
71 
72 
73  /// gets source strength
74  virtual void get_source(const unsigned& ipt,
75  const Vector<double>& x,
76  double& source) const
77  {
78  // if source function is not provided, i.e. zero, then return zero
79  if (Source_fct_pt == 0)
80  {
81  source = 0.0;
82  }
83 
84  // else get source strength
85  else
86  {
87  (*Source_fct_pt)(x, source);
88  }
89  }
90 
91 
92  /// wrapper function, adds contribution to residual
94  {
95  // create a dummy matrix
96  DenseDoubleMatrix dummy(1);
97 
98  // call the generic residuals functions with flag set to zero
100  }
101 
102 
103  /// wrapper function, adds contribution to residual and generic
105  DenseMatrix<double>& jacobian)
106  {
107  // call generic routine with flag set to 1
109  }
110 
111 
112  /// output with nplot points
113  void output(std::ostream& outfile, const unsigned& nplot)
114  {
115  // Vector of local coordinates
116  Vector<double> s(DIM);
117 
118  // Tecplot header info
119  outfile << this->tecplot_zone_string(nplot);
120 
121  // Loop over plot points
122  unsigned num_plot_points = this->nplot_points(nplot);
123  for (unsigned iplot = 0; iplot < num_plot_points; iplot++)
124  {
125  // Get local coordinates of plot point
126  this->get_s_plot(iplot, nplot, s);
127 
128  for (unsigned i = 0; i < DIM; i++)
129  {
130  outfile << this->interpolated_x(s, i) << " ";
131  }
132 
133  outfile << interpolated_u_biharmonic(s) << std::endl;
134  }
135 
136  // Write tecplot footer (e.g. FE connectivity lists)
137  this->write_tecplot_zone_footer(outfile, nplot);
138  }
139 
140 
141  /// Output at default number of plot points
142  void output(std::ostream& outfile)
143  {
144  FiniteElement::output(outfile);
145  }
146 
147  /// C-style output
148  void output(FILE* file_pt)
149  {
150  FiniteElement::output(file_pt);
151  }
152 
153  /// C_style output at n_plot points
154  void output(FILE* file_pt, const unsigned& n_plot)
155  {
156  FiniteElement::output(file_pt, n_plot);
157  }
158 
159 
160  /// output fluid velocity field
161  void output_fluid_velocity(std::ostream& outfile, const unsigned& nplot)
162  {
163  // Vector of local coordinates
164  Vector<double> s(DIM);
165 
166  // Tecplot header info
167  outfile << this->tecplot_zone_string(nplot);
168 
169  // Loop over plot points
170  unsigned num_plot_points = this->nplot_points(nplot);
171  for (unsigned iplot = 0; iplot < num_plot_points; iplot++)
172  {
173  // Get local coordinates of plot point
174  this->get_s_plot(iplot, nplot, s);
175 
176  for (unsigned i = 0; i < DIM; i++)
177  {
178  outfile << this->interpolated_x(s, i) << " ";
179  }
180 
181  Vector<double> dudx(2, 0.0);
182  interpolated_dudx(s, dudx);
183 
184  outfile << dudx[1] << " " << -dudx[0] << std::endl;
185  }
186 
187  // Write tecplot footer (e.g. FE connectivity lists)
188  this->write_tecplot_zone_footer(outfile, nplot);
189  }
190 
191 
192  /// output with nplot points
194  {
195  // Find out how many nodes there are
196  unsigned n_node = this->nnode();
197 
198  // Find out how many values there
199  unsigned n_value = this->node_pt(0)->nvalue();
200 
201  // set up memory for shape functions
202  Shape psi(n_node, n_value);
203  DShape dpsidx(n_node, n_value, DIM);
204 
205  // evaluate dpsidx at local coordinate s
206  dshape_eulerian(s, psi, dpsidx);
207 
208  // initialise storage for d2u_interpolated
209  dudx[0] = 0.0;
210  dudx[1] = 0.0;
211 
212  // loop over nodes, degrees of freedom, and dimension to calculate
213  // interpolated_d2u
214  for (unsigned n = 0; n < n_node; n++)
215  {
216  for (unsigned k = 0; k < n_value; k++)
217  {
218  for (unsigned d = 0; d < DIM; d++)
219  {
220  dudx[d] += this->node_pt(n)->value(k) * dpsidx(n, k, d);
221  }
222  }
223  }
224  }
225 
226 
227  /// output analytic solution
228  void output_fct(std::ostream& outfile,
229  const unsigned& nplot,
231  {
232  // Vector of local coordinates
233  Vector<double> s(DIM);
234 
235  // Vector for coordintes
236  Vector<double> x(DIM);
237 
238  // Tecplot header info
239  outfile << this->tecplot_zone_string(nplot);
240 
241  // Exact solution Vector (here a scalar)
242  Vector<double> exact_soln(1);
243 
244  // Loop over plot points
245  unsigned num_plot_points = this->nplot_points(nplot);
246  for (unsigned iplot = 0; iplot < num_plot_points; iplot++)
247  {
248  // Get local coordinates of plot point
249  this->get_s_plot(iplot, nplot, s);
250 
251  // Get x position as Vector
252  this->interpolated_x(s, x);
253 
254  // Get exact solution at this point
255  (*exact_soln_pt)(x, exact_soln);
256 
257  // Output x,y,...,u_exact
258  for (unsigned i = 0; i < DIM; i++)
259  {
260  outfile << x[i] << " ";
261  }
262  outfile << exact_soln[0] << std::endl;
263  }
264 
265  // Write tecplot footer (e.g. FE connectivity lists)
266  this->write_tecplot_zone_footer(outfile, nplot);
267  }
268 
269  /// Output exact solution specified via function pointer
270  /// at a given time and at a given number of plot points.
271  /// Function prints as many components as are returned in solution Vector.
272  /// Implement broken FiniteElement base class version
273  void output_fct(std::ostream& outfile,
274  const unsigned& nplot,
275  const double& time,
277  {
278  FiniteElement::output_fct(outfile, nplot, time, exact_soln_pt);
279  }
280 
281 
282  /// computes the error
283  void compute_error(std::ostream& outfile,
285  double& error,
286  double& norm)
287  {
288  // Initialise
289  error = 0.0;
290  norm = 0.0;
291 
292  // Vector of local coordinates
293  Vector<double> s(DIM);
294 
295  // Vector for coordintes
296  Vector<double> x(DIM);
297 
298  // Find out how many nodes there are in the element
299  unsigned n_node = this->nnode();
300 
301  Shape psi(n_node);
302 
303  // Set the value of n_intpt
304  unsigned n_intpt = this->integral_pt()->nweight();
305 
306  // Tecplot header info
307  outfile << this->tecplot_zone_string(3);
308 
309  // Tecplot
310  // outfile << "ZONE" << std::endl;
311 
312  // Exact solution Vector (here a scalar)
313  Vector<double> exact_soln(1);
314 
315  // Loop over the integration points
316  for (unsigned ipt = 0; ipt < n_intpt; ipt++)
317  {
318  // Assign values of s
319  for (unsigned i = 0; i < DIM; i++)
320  {
321  s[i] = this->integral_pt()->knot(ipt, i);
322  }
323 
324  // Get the integral weight
325  double w = this->integral_pt()->weight(ipt);
326 
327  // Get jacobian of mapping
328  double J = this->J_eulerian(s);
329 
330  // Premultiply the weights and the Jacobian
331  double W = w * J;
332 
333  // Get x position as Vector
334  this->interpolated_x(s, x);
335 
336  // Get FE function value
337  double u_fe = interpolated_u_biharmonic(s);
338 
339  // Get exact solution at this point
340  (*exact_soln_pt)(x, exact_soln);
341 
342  // Output x,y,...,error
343  for (unsigned i = 0; i < DIM; i++)
344  {
345  outfile << x[i] << " ";
346  }
347  outfile << exact_soln[0] << " " << fabs(exact_soln[0] - u_fe)
348  << std::endl;
349 
350  // Add to error and norm
351  norm += exact_soln[0] * exact_soln[0] * W;
352  error += (exact_soln[0] - u_fe) * (exact_soln[0] - u_fe) * W;
353  }
354 
355  this->write_tecplot_zone_footer(outfile, 3);
356  }
357 
358 
359  /// Plot the error when compared against a given time-dependent exact
360  /// solution \f$ {\bf f}(t,{\bf x}) \f$. Also calculates the norm of the
361  /// error and that of the exact solution. Call broken base-class version.
362  void compute_error(std::ostream& outfile,
364  const double& time,
365  double& error,
366  double& norm)
367  {
368  FiniteElement::compute_error(outfile, exact_soln_pt, time, error, norm);
369  }
370 
371  /// calculates interpolated u at s
373  {
374  // initialise storage for u_interpolated
375  double uu = 0;
376 
377  // number of nodes
378  unsigned n_node = this->nnode();
379 
380  // number of degrees of freedom per node
381  unsigned n_value = this->node_pt(0)->nvalue();
382 
383  // set up memory for shape functions
384  Shape psi(n_node, n_value);
385 
386  // find shape fn at position s
387  this->shape(s, psi);
388 
389  // calculate interpolated u
390  for (unsigned n = 0; n < n_node; n++)
391  {
392  for (unsigned k = 0; k < n_value; k++)
393  {
394  uu += u(n, k) * psi(n, k);
395  }
396  }
397 
398  // return interpolated_u
399  return uu;
400  }
401 
402 
403  /// self test wrapper
404  unsigned self_test()
405  {
406  bool passed = true;
407 
408  // Check lower-level stuff
409  if (FiniteElement::self_test() != 0)
410  {
411  passed = false;
412  }
413 
414  // Return verdict
415  if (passed)
416  {
417  return 0;
418  }
419  else
420  {
421  return 1;
422  }
423  }
424 
425 
426  /// return number of second derivate degrees of freedom
427  unsigned get_d2_dof()
428  {
429  return d2_dof[DIM - 1];
430  }
431 
432 
433  /// The number of "DOF types" that degrees of freedom in this element
434  /// are sub-divided into (for block preconditioning)
435  unsigned ndof_types() const
436  {
437  return this->required_nvalue(1);
438  }
439 
440 
441  /// Create a list of pairs for all unknowns in this element,
442  /// so that the first entry in each pair contains the global equation
443  /// number of the unknown, while the second one contains the number
444  /// of the "DOF types" that this unknown is associated with.
445  /// (Function can obviously only be called if the equation numbering
446  /// scheme has been set up.) (for block preconditioning)
448  std::list<std::pair<unsigned long, unsigned>>& dof_lookup_list) const
449  {
450  // number of nodes
451  int n_node = this->nnode();
452 
453  // number of degrees of freedom at each node
454  int n_value = this->node_pt(0)->nvalue();
455 
456  // temporary pair (used to store dof lookup prior to being added to list
457  std::pair<unsigned long, unsigned> dof_lookup;
458 
459  // loop over the nodes
460  for (int n = 0; n < n_node; n++)
461  {
462  // loop over the degree of freedom
463  for (int k = 0; k < n_value; k++)
464  {
465  // determine local eqn number
466  int local_eqn_number = this->nodal_local_eqn(n, k);
467 
468  // if local equation number is less than zero then nodal dof pinned
469  // then ignore
470  if (local_eqn_number >= 0)
471  {
472  // store dof lookup in temporary pair
473  dof_lookup.first = this->eqn_number(local_eqn_number);
474  dof_lookup.second = k;
475  dof_lookup_list.push_front(dof_lookup);
476  // add to list
477  }
478  }
479  }
480  }
481 
482 
483  /// Access functions for the source function pointer
485  {
486  return Source_fct_pt;
487  }
488 
489  /// Access functions for the source function pointers (const version)
491  {
492  return Source_fct_pt;
493  }
494 
495 
496  protected:
497  /// Compute element residual Vector only (if JFLAG=and/or element
498  /// Jacobian matrix
500  Vector<double>& residuals, DenseMatrix<double>& jacobian, unsigned JFLAG);
501 
502 
503  /// Pointer to source function:
505 
506 
507  /// Array to hold local eqn numbers: Local_eqn[n] (=-1 for BC)
509 
510  private:
511  // number of degrees of freedom of second derivative
512  static const unsigned d2_dof[];
513  };
514 
515 
516  // declares number of degrees of freedom of second derivative
517  template<unsigned DIM>
518  const unsigned BiharmonicEquations<DIM>::d2_dof[3] = {1, 3, 6};
519 
520 
521  //=============================================================================
522  /// biharmonic element class
523  //=============================================================================
524  template<unsigned DIM>
525  class BiharmonicElement : public virtual QHermiteElement<DIM>,
526  public virtual BiharmonicEquations<DIM>
527  {
528  public:
529  /// access function for value, kth dof of node n
530  inline double u(const unsigned& n, const unsigned& k) const
531  {
532  return this->node_pt(n)->value(k);
533  }
534 
535 
536  /// Constructor: Call constructors for QElement and
537  /// Poisson equations
539 
540 
542 
543 
544  /// Required # of `values' (pinned or dofs)
545  /// at node n
546  inline unsigned required_nvalue(const unsigned& n) const
547  {
548  return DIM * 2;
549  }
550 
551 
552  /// Output
553  void output(std::ostream& outfile)
554  {
556  }
557 
558  /// output wrapper
559  void output(std::ostream& outfile, const unsigned& n_plot)
560  {
561  BiharmonicEquations<DIM>::output(outfile, n_plot);
562  }
563 
564  /// C-style output
565  void output(FILE* file_pt)
566  {
568  }
569 
570  /// C_style output at n_plot points
571  void output(FILE* file_pt, const unsigned& n_plot)
572  {
573  BiharmonicEquations<DIM>::output(file_pt, n_plot);
574  }
575 
576 
577  /// analytic solution wrapper
578  void output_fct(std::ostream& outfile,
579  const unsigned& nplot,
581  {
582  BiharmonicEquations<DIM>::output_fct(outfile, nplot, exact_soln_pt);
583  }
584 
585 
586  /// Final override
587  void output_fct(std::ostream& outfile,
588  const unsigned& nplot,
589  const double& time,
591  {
592  BiharmonicEquations<DIM>::output_fct(outfile, nplot, time, exact_soln_pt);
593  }
594 
595 
596  /// computes error
597  void compute_error(std::ostream& outfile,
599  double& error,
600  double& norm)
601  {
603  outfile, exact_soln_pt, error, norm);
604  }
605 
606  /// Call the equations-class overloaded unsteady error calculation
607  void compute_error(std::ostream& outfile,
609  const double& time,
610  double& error,
611  double& norm)
612  {
614  outfile, exact_soln_pt, time, error, norm);
615  }
616  };
617 
618 
619 } // namespace oomph
620 #endif
static char t char * s
Definition: cfortran.h:568
cstr elem_len * i
Definition: cfortran.h:603
biharmonic element class
void output(std::ostream &outfile)
Output.
void compute_error(std::ostream &outfile, FiniteElement::UnsteadyExactSolutionFctPt exact_soln_pt, const double &time, double &error, double &norm)
Call the equations-class overloaded unsteady error calculation.
unsigned required_nvalue(const unsigned &n) const
Required # of ‘values’ (pinned or dofs) at node n.
void output_fct(std::ostream &outfile, const unsigned &nplot, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt)
analytic solution wrapper
void output(FILE *file_pt)
C-style output.
void compute_error(std::ostream &outfile, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt, double &error, double &norm)
computes error
void output(FILE *file_pt, const unsigned &n_plot)
C_style output at n_plot points.
void output_fct(std::ostream &outfile, const unsigned &nplot, const double &time, FiniteElement::UnsteadyExactSolutionFctPt exact_soln_pt)
Final override.
double u(const unsigned &n, const unsigned &k) const
access function for value, kth dof of node n
void output(std::ostream &outfile, const unsigned &n_plot)
output wrapper
BiharmonicElement()
Constructor: Call constructors for QElement and Poisson equations.
Biharmonic Equation Class - contains the equations.
void fill_in_contribution_to_residuals(Vector< double > &residuals)
wrapper function, adds contribution to residual
void output_fluid_velocity(std::ostream &outfile, const unsigned &nplot)
output fluid velocity field
SourceFctPt Source_fct_pt
Pointer to source function:
void output(std::ostream &outfile)
Output at default number of plot points.
void compute_error(std::ostream &outfile, FiniteElement::UnsteadyExactSolutionFctPt exact_soln_pt, const double &time, double &error, double &norm)
Plot the error when compared against a given time-dependent exact solution . Also calculates the norm...
virtual void get_source(const unsigned &ipt, const Vector< double > &x, double &source) const
gets source strength
BiharmonicEquations()
Constructor (must initialise the Source_fct_pt to null)
static const unsigned d2_dof[]
unsigned self_test()
self test wrapper
void output_fct(std::ostream &outfile, const unsigned &nplot, const double &time, FiniteElement::UnsteadyExactSolutionFctPt exact_soln_pt)
Output exact solution specified via function pointer at a given time and at a given number of plot po...
virtual void fill_in_generic_residual_contribution_biharmonic(Vector< double > &residuals, DenseMatrix< double > &jacobian, unsigned JFLAG)
Compute element residual Vector only (if JFLAG=and/or element Jacobian matrix.
void output(std::ostream &outfile, const unsigned &nplot)
output with nplot points
void interpolated_dudx(const Vector< double > &s, Vector< double > &dudx)
output with nplot points
void compute_error(std::ostream &outfile, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt, double &error, double &norm)
computes the error
void output(FILE *file_pt)
C-style output.
void(* SourceFctPt)(const Vector< double > &x, double &u)
source function type definition
void fill_in_contribution_to_jacobian(Vector< double > &residual, DenseMatrix< double > &jacobian)
wrapper function, adds contribution to residual and generic
unsigned ndof_types() const
The number of "DOF types" that degrees of freedom in this element are sub-divided into (for block pre...
double interpolated_u_biharmonic(const Vector< double > &s)
calculates interpolated u at s
void output_fct(std::ostream &outfile, const unsigned &nplot, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt)
output analytic solution
SourceFctPt & source_fct_pt()
Access functions for the source function pointer.
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 get_d2_dof()
return number of second derivate degrees of freedom
virtual double u(const unsigned &n, const unsigned &k) const =0
Access function: Nodal function value at local node n Uses suitably interpolated value for hanging no...
Vector< int > Local_eqn
Array to hold local eqn numbers: Local_eqn[n] (=-1 for BC)
void output(FILE *file_pt, const unsigned &n_plot)
C_style output at n_plot points.
SourceFctPt source_fct_pt() const
Access functions for the source function pointers (const version)
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
Class of matrices containing doubles, and stored as a DenseMatrix<double>, but with solving functiona...
Definition: matrices.h:1271
A general Finite Element class.
Definition: elements.h:1313
virtual double J_eulerian(const Vector< double > &s) const
Return the Jacobian of mapping from local to global coordinates at local position s.
Definition: elements.cc:4103
virtual void output_fct(std::ostream &outfile, const unsigned &n_plot, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt)
Output an exact solution over the element.
Definition: elements.h:3104
Node *& node_pt(const unsigned &n)
Return a pointer to the local node n.
Definition: elements.h:2175
virtual void output(std::ostream &outfile)
Output the element data — typically the values at the nodes in a format suitable for post-processing.
Definition: elements.h:3050
virtual std::string tecplot_zone_string(const unsigned &nplot) const
Return string for tecplot zone header (when plotting nplot points in each "coordinate direction")
Definition: elements.h:3161
virtual unsigned required_nvalue(const unsigned &n) const
Number of values that must be stored at local node n by the element. The default is 0,...
Definition: elements.h:2455
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...
virtual double interpolated_x(const Vector< double > &s, const unsigned &i) const
Return FE interpolated coordinate x[i] at local coordinate s.
Definition: elements.cc:3962
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 nnode() const
Return the number of nodes.
Definition: elements.h:2210
virtual void compute_error(FiniteElement::SteadyExactSolutionFctPt exact_soln_pt, double &error, double &norm)
Calculate the norm of the error and that of the exact solution.
Definition: elements.h:3198
void(* SteadyExactSolutionFctPt)(const Vector< double > &, Vector< double > &)
Function pointer for function that computes vector-valued steady "exact solution" as .
Definition: elements.h:1759
Integral *const & integral_pt() const
Return the pointer to the integration scheme (const version)
Definition: elements.h:1963
virtual void get_s_plot(const unsigned &i, const unsigned &nplot, Vector< double > &s, const bool &shifted_to_interior=false) const
Get cector of local coordinates of plot point i (when plotting nplot points in each "coordinate direc...
Definition: elements.h:3148
virtual unsigned nplot_points(const unsigned &nplot) const
Return total number of plot points (when plotting nplot points in each "coordinate direction")
Definition: elements.h:3186
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
virtual void write_tecplot_zone_footer(std::ostream &outfile, const unsigned &nplot) const
Add tecplot zone "footer" to output stream (when plotting nplot points in each "coordinate direction"...
Definition: elements.h:3174
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
virtual unsigned self_test()
Self-test: Check inversion of element & do self-test for GeneralisedElement. Return 0 if OK.
Definition: elements.cc:4440
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
int local_eqn_number(const unsigned long &ieqn_global) const
Return the local equation number corresponding to the ieqn_global-th global equation number....
Definition: elements.h:726
virtual double knot(const unsigned &i, const unsigned &j) const =0
Return local coordinate s[j] of i-th integration point.
virtual unsigned nweight() const =0
Return the number of integration points of the scheme.
virtual double weight(const unsigned &i) const =0
Return weight of i-th integration point.
double value(const unsigned &i) const
Return i-th value (dofs or pinned) at this node either directly or via hanging node representation....
Definition: nodes.cc:2408
/////////////////////////////////////////////////////////////////// /////////////////////////////////...
A Class for shape functions. In simple cases, the shape functions have only one index that can be tho...
Definition: shape.h:76
//////////////////////////////////////////////////////////////////// ////////////////////////////////...