foeppl_von_karman_elements.cc
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 // Non-inline functions for FoepplvonKarman elements
27 
29 
30 #include <iostream>
31 
32 namespace oomph
33 {
34  //======================================================================
35  /// Set the data for the number of Variables at each node - 8
36  //======================================================================
37  template<unsigned NNODE_1D>
39 
40  // Foeppl von Karman equations static data
41 
42  /// Default value physical constants
44 
45  //======================================================================
46  /// Compute contribution to element residual Vector
47  ///
48  /// Pure version without hanging nodes
49  //======================================================================
51  Vector<double>& residuals)
52  {
53  // Find out how many nodes there are
54  const unsigned n_node = nnode();
55 
56  // Set up memory for the shape and test functions
57  Shape psi(n_node), test(n_node);
58  DShape dpsidx(n_node, 2), dtestdx(n_node, 2);
59 
60  // Indices at which the unknowns are stored
61  const unsigned w_nodal_index = nodal_index_fvk(0);
62  const unsigned laplacian_w_nodal_index = nodal_index_fvk(1);
63  const unsigned phi_nodal_index = nodal_index_fvk(2);
64  const unsigned laplacian_phi_nodal_index = nodal_index_fvk(3);
65  const unsigned smooth_dwdx_nodal_index = nodal_index_fvk(4);
66  const unsigned smooth_dwdy_nodal_index = nodal_index_fvk(5);
67  const unsigned smooth_dphidx_nodal_index = nodal_index_fvk(6);
68  const unsigned smooth_dphidy_nodal_index = nodal_index_fvk(7);
69 
70  // Set the value of n_intpt
71  const unsigned n_intpt = integral_pt()->nweight();
72 
73  // Integers to store the local equation numbers
74  int local_eqn = 0;
75 
76  // Loop over the integration points
77  for (unsigned ipt = 0; ipt < n_intpt; ipt++)
78  {
79  // Get the integral weight
80  double w = integral_pt()->weight(ipt);
81 
82  // Call the derivatives of the shape and test functions
83  double J =
84  dshape_and_dtest_eulerian_at_knot_fvk(ipt, psi, dpsidx, test, dtestdx);
85 
86  // Premultiply the weights and the Jacobian
87  double W = w * J;
88 
89  // Allocate and initialise to zero storage for the interpolated values
91 
92  double interpolated_w = 0;
93  double interpolated_laplacian_w = 0;
94  double interpolated_phi = 0;
95  double interpolated_laplacian_phi = 0;
96 
97  Vector<double> interpolated_dwdx(2, 0.0);
98  Vector<double> interpolated_dlaplacian_wdx(2, 0.0);
99  Vector<double> interpolated_dphidx(2, 0.0);
100  Vector<double> interpolated_dlaplacian_phidx(2, 0.0);
101 
102  Vector<double> interpolated_smooth_dwdx(2, 0.0);
103  Vector<double> interpolated_smooth_dphidx(2, 0.0);
104  double interpolated_continuous_d2wdx2 = 0;
105  double interpolated_continuous_d2wdy2 = 0;
106  double interpolated_continuous_d2phidx2 = 0;
107  double interpolated_continuous_d2phidy2 = 0;
108  double interpolated_continuous_d2wdxdy = 0;
109  double interpolated_continuous_d2phidxdy = 0;
110 
111  // Calculate function values and derivatives:
112  //-----------------------------------------
113  Vector<double> nodal_value(8, 0.0);
114  // Loop over nodes
115  for (unsigned l = 0; l < n_node; l++)
116  {
117  // Get the nodal values
118  nodal_value[0] = raw_nodal_value(l, w_nodal_index);
119  nodal_value[1] = raw_nodal_value(l, laplacian_w_nodal_index);
120 
122  {
123  nodal_value[2] = raw_nodal_value(l, phi_nodal_index);
124  nodal_value[3] = raw_nodal_value(l, laplacian_phi_nodal_index);
125  nodal_value[4] = raw_nodal_value(l, smooth_dwdx_nodal_index);
126  nodal_value[5] = raw_nodal_value(l, smooth_dwdy_nodal_index);
127  nodal_value[6] = raw_nodal_value(l, smooth_dphidx_nodal_index);
128  nodal_value[7] = raw_nodal_value(l, smooth_dphidy_nodal_index);
129  }
130 
131  // Add contributions from current node/shape function
132  interpolated_w += nodal_value[0] * psi(l);
133  interpolated_laplacian_w += nodal_value[1] * psi(l);
134 
136  {
137  interpolated_phi += nodal_value[2] * psi(l);
138  interpolated_laplacian_phi += nodal_value[3] * psi(l);
139 
140  interpolated_smooth_dwdx[0] += nodal_value[4] * psi(l);
141  interpolated_smooth_dwdx[1] += nodal_value[5] * psi(l);
142  interpolated_smooth_dphidx[0] += nodal_value[6] * psi(l);
143  interpolated_smooth_dphidx[1] += nodal_value[7] * psi(l);
144 
145  interpolated_continuous_d2wdx2 += nodal_value[4] * dpsidx(l, 0);
146  interpolated_continuous_d2wdy2 += nodal_value[5] * dpsidx(l, 1);
147  interpolated_continuous_d2phidx2 += nodal_value[6] * dpsidx(l, 0);
148  interpolated_continuous_d2phidy2 += nodal_value[7] * dpsidx(l, 1);
149  // mjr CHECK THESE
150  interpolated_continuous_d2wdxdy +=
151  0.5 *
152  (nodal_value[4] * dpsidx(l, 1) + nodal_value[5] * dpsidx(l, 0));
153  interpolated_continuous_d2phidxdy +=
154  0.5 *
155  (nodal_value[6] * dpsidx(l, 1) + nodal_value[7] * dpsidx(l, 0));
156  }
157 
158  // Loop over directions
159  for (unsigned j = 0; j < 2; j++)
160  {
161  interpolated_x[j] += raw_nodal_position(l, j) * psi(l);
162  interpolated_dwdx[j] += nodal_value[0] * dpsidx(l, j);
163  interpolated_dlaplacian_wdx[j] += nodal_value[1] * dpsidx(l, j);
164 
166  {
167  interpolated_dphidx[j] += nodal_value[2] * dpsidx(l, j);
168  interpolated_dlaplacian_phidx[j] += nodal_value[3] * dpsidx(l, j);
169  }
170  }
171  }
172 
173  // Get pressure function
174  //-------------------
175  double pressure;
176  get_pressure_fvk(ipt, interpolated_x, pressure);
177 
178  double airy_forcing;
179  get_airy_forcing_fvk(ipt, interpolated_x, airy_forcing);
180 
181  // Assemble residuals and Jacobian
182  //--------------------------------
183 
184  // Loop over the test functions
185  for (unsigned l = 0; l < n_node; l++)
186  {
187  // Get the local equation
188  local_eqn = nodal_local_eqn(l, w_nodal_index);
189 
190  // IF it's not a boundary condition
191  if (local_eqn >= 0)
192  {
193  residuals[local_eqn] += pressure * test(l) * W;
194 
195  // Reduced order biharmonic operator
196  for (unsigned k = 0; k < 2; k++)
197  {
198  residuals[local_eqn] +=
199  interpolated_dlaplacian_wdx[k] * dtestdx(l, k) * W;
200  }
201 
203  {
204  // Monge-Ampere part
205  residuals[local_eqn] += eta() *
206  (interpolated_continuous_d2wdx2 *
207  interpolated_continuous_d2phidy2 +
208  interpolated_continuous_d2wdy2 *
209  interpolated_continuous_d2phidx2 -
210  2.0 * interpolated_continuous_d2wdxdy *
211  interpolated_continuous_d2phidxdy) *
212  test(l) * W;
213  }
214  }
215 
216  // Get the local equation
217  local_eqn = nodal_local_eqn(l, laplacian_w_nodal_index);
218 
219  // IF it's not a boundary condition
220  if (local_eqn >= 0)
221  {
222  // The coupled Poisson equations for the biharmonic operator
223  residuals[local_eqn] += interpolated_laplacian_w * test(l) * W;
224 
225  for (unsigned k = 0; k < 2; k++)
226  {
227  residuals[local_eqn] += interpolated_dwdx[k] * dtestdx(l, k) * W;
228  }
229  }
230 
231  // Get the local equation
232  local_eqn = nodal_local_eqn(l, phi_nodal_index);
233 
234  // IF it's not a boundary condition
235  if (local_eqn >= 0)
236  {
237  residuals[local_eqn] += airy_forcing * test(l) * W;
238 
239  // Reduced order biharmonic operator
240  for (unsigned k = 0; k < 2; k++)
241  {
242  residuals[local_eqn] +=
243  interpolated_dlaplacian_phidx[k] * dtestdx(l, k) * W;
244  }
245 
246  // Monge-Ampere part
247  residuals[local_eqn] -=
248  (interpolated_continuous_d2wdx2 * interpolated_continuous_d2wdy2 -
249  interpolated_continuous_d2wdxdy *
250  interpolated_continuous_d2wdxdy) *
251  test(l) * W;
252  }
253 
254  // Get the local equation
255  local_eqn = nodal_local_eqn(l, laplacian_phi_nodal_index);
256 
257  // IF it's not a boundary condition
258  if (local_eqn >= 0)
259  {
260  // The coupled Poisson equations for the biharmonic operator
261  residuals[local_eqn] += interpolated_laplacian_phi * test(l) * W;
262 
263  for (unsigned k = 0; k < 2; k++)
264  {
265  residuals[local_eqn] += interpolated_dphidx[k] * dtestdx(l, k) * W;
266  }
267  }
268 
269  // Residuals for the smooth derivatives
270  local_eqn = nodal_local_eqn(l, smooth_dwdx_nodal_index);
271 
272  if (local_eqn >= 0)
273  {
274  residuals[local_eqn] +=
275  (interpolated_dwdx[0] - interpolated_smooth_dwdx[0]) * test(l) * W;
276  }
277 
278  local_eqn = nodal_local_eqn(l, smooth_dwdy_nodal_index);
279 
280  if (local_eqn >= 0)
281  {
282  residuals[local_eqn] +=
283  (interpolated_dwdx[1] - interpolated_smooth_dwdx[1]) * test(l) * W;
284  }
285 
286  local_eqn = nodal_local_eqn(l, smooth_dphidx_nodal_index);
287 
288  if (local_eqn >= 0)
289  {
290  residuals[local_eqn] +=
291  (interpolated_dphidx[0] - interpolated_smooth_dphidx[0]) * test(l) *
292  W;
293  }
294 
295  local_eqn = nodal_local_eqn(l, smooth_dphidy_nodal_index);
296 
297  if (local_eqn >= 0)
298  {
299  residuals[local_eqn] +=
300  (interpolated_dphidx[1] - interpolated_smooth_dphidx[1]) * test(l) *
301  W;
302  }
303  }
304 
305  } // End of loop over integration points
306 
307 
308  // Finally: My contribution to the volume constraint equation
309  // (if any). Note this must call get_bounded_volume since the
310  // definition of the bounded volume can be overloaded in derived
311  // elements.
313  {
314  local_eqn =
316  if (local_eqn >= 0)
317  {
318  residuals[local_eqn] += get_bounded_volume();
319  }
320  }
321  }
322 
323  /*
324  void FoepplvonKarmanEquations::fill_in_contribution_to_jacobian(Vector<double>
325  &residuals, DenseMatrix<double> &jacobian)
326  {
327  //Add the contribution to the residuals
328  FoepplvonKarmanEquations::fill_in_contribution_to_residuals(residuals);
329  //Allocate storage for the full residuals (residuals of entire element)
330  unsigned n_dof = ndof();
331  Vector<double> full_residuals(n_dof);
332  //Get the residuals for the entire element
333  FoepplvonKarmanEquations::get_residuals(full_residuals);
334  //Calculate the contributions from the internal dofs
335  //(finite-difference the lot by default)
336  fill_in_jacobian_from_internal_by_fd(full_residuals,jacobian,true);
337  //Calculate the contributions from the external dofs
338  //(finite-difference the lot by default)
339  fill_in_jacobian_from_external_by_fd(full_residuals,jacobian,true);
340  //Calculate the contributions from the nodal dofs
341  fill_in_jacobian_from_nodal_by_fd(full_residuals,jacobian);
342  }
343  */
344 
345 
346  //======================================================================
347  /// Self-test: Return 0 for OK
348  //======================================================================
350  {
351  bool passed = true;
352 
353  // Check lower-level stuff
354  if (FiniteElement::self_test() != 0)
355  {
356  passed = false;
357  }
358 
359  // Return verdict
360  if (passed)
361  {
362  return 0;
363  }
364  else
365  {
366  return 1;
367  }
368  }
369 
370 
371  //======================================================================
372  /// Compute in-plane stresses
373  //======================================================================
375  double& sigma_xx,
376  double& sigma_yy,
377  double& sigma_xy)
378  {
379  // No in plane stresses if linear bending
381  {
382  sigma_xx = 0.0;
383  sigma_yy = 0.0;
384  sigma_xy = 0.0;
385  return;
386  }
387 
388  // Get shape fcts and derivs
389  unsigned n_dim = this->dim();
390  unsigned n_node = this->nnode();
391  Shape psi(n_node);
392  DShape dpsidx(n_node, n_dim);
393  dshape_eulerian(s, psi, dpsidx);
394  double interpolated_continuous_d2phidx2 = 0;
395  double interpolated_continuous_d2phidy2 = 0;
396  double interpolated_continuous_d2phidxdy = 0;
397 
398  const unsigned smooth_dphidx_nodal_index = nodal_index_fvk(6);
399  const unsigned smooth_dphidy_nodal_index = nodal_index_fvk(7);
400 
401  // Loop over nodes
402  for (unsigned l = 0; l < n_node; l++)
403  {
404  interpolated_continuous_d2phidx2 +=
405  raw_nodal_value(l, smooth_dphidx_nodal_index) * dpsidx(l, 0);
406  interpolated_continuous_d2phidy2 +=
407  raw_nodal_value(l, smooth_dphidy_nodal_index) * dpsidx(l, 1);
408  interpolated_continuous_d2phidxdy +=
409  0.5 * (raw_nodal_value(l, smooth_dphidx_nodal_index) * dpsidx(l, 1) +
410  raw_nodal_value(l, smooth_dphidy_nodal_index) * dpsidx(l, 0));
411  }
412 
413  sigma_xx = interpolated_continuous_d2phidy2;
414  sigma_yy = interpolated_continuous_d2phidx2;
415  sigma_xy = -interpolated_continuous_d2phidxdy;
416  }
417 
418 
419  //======================================================================
420  /// Output function:
421  ///
422  /// x,y,w
423  ///
424  /// nplot points in each coordinate direction
425  //======================================================================
426  void FoepplvonKarmanEquations::output(std::ostream& outfile,
427  const unsigned& nplot)
428  {
429  // Vector of local coordinates
430  Vector<double> s(2);
431 
432  // Tecplot header info
433  outfile << tecplot_zone_string(nplot);
434 
435  // Loop over plot points
436  unsigned num_plot_points = nplot_points(nplot);
437  for (unsigned iplot = 0; iplot < num_plot_points; iplot++)
438  {
439  // Get local coordinates of plot point
440  get_s_plot(iplot, nplot, s);
441 
442  for (unsigned i = 0; i < 2; i++)
443  {
444  outfile << interpolated_x(s, i) << " ";
445  }
446  outfile << interpolated_w_fvk(s) << std::endl;
447  }
448 
449  // Write tecplot footer (e.g. FE connectivity lists)
450  write_tecplot_zone_footer(outfile, nplot);
451  }
452 
453 
454  //======================================================================
455  /// C-style output function:
456  ///
457  /// x,y,w
458  ///
459  /// nplot points in each coordinate direction
460  //======================================================================
461  void FoepplvonKarmanEquations::output(FILE* file_pt, const unsigned& nplot)
462  {
463  // Vector of local coordinates
464  Vector<double> s(2);
465 
466  // Tecplot header info
467  fprintf(file_pt, "%s", tecplot_zone_string(nplot).c_str());
468 
469  // Loop over plot points
470  unsigned num_plot_points = nplot_points(nplot);
471  for (unsigned iplot = 0; iplot < num_plot_points; iplot++)
472  {
473  // Get local coordinates of plot point
474  get_s_plot(iplot, nplot, s);
475 
476  for (unsigned i = 0; i < 2; i++)
477  {
478  fprintf(file_pt, "%g ", interpolated_x(s, i));
479  }
480  fprintf(file_pt, "%g \n", interpolated_w_fvk(s));
481  }
482 
483  // Write tecplot footer (e.g. FE connectivity lists)
484  write_tecplot_zone_footer(file_pt, nplot);
485  }
486 
487 
488  //======================================================================
489  /// Output exact solution
490  ///
491  /// Solution is provided via function pointer.
492  /// Plot at a given number of plot points.
493  ///
494  /// x,y,w_exact
495  //======================================================================
497  std::ostream& outfile,
498  const unsigned& nplot,
500  {
501  // Vector of local coordinates
502  Vector<double> s(2);
503 
504  // Vector for coordintes
505  Vector<double> x(2);
506 
507  // Tecplot header info
508  outfile << tecplot_zone_string(nplot);
509 
510  // Exact solution Vector (here a scalar)
511  // Vector<double> exact_soln(1);
512  Vector<double> exact_soln(2);
513 
514  // Loop over plot points
515  unsigned num_plot_points = nplot_points(nplot);
516  for (unsigned iplot = 0; iplot < num_plot_points; iplot++)
517  {
518  // Get local coordinates of plot point
519  get_s_plot(iplot, nplot, s);
520 
521  // Get x position as Vector
522  interpolated_x(s, x);
523 
524  // Get exact solution at this point
525  (*exact_soln_pt)(x, exact_soln);
526 
527  // Output x,y,w_exact
528  for (unsigned i = 0; i < 2; i++)
529  {
530  outfile << x[i] << " ";
531  }
532  outfile << exact_soln[0] << std::endl;
533  }
534 
535  // Write tecplot footer (e.g. FE connectivity lists)
536  write_tecplot_zone_footer(outfile, nplot);
537  }
538 
539 
540  //======================================================================
541  /// Validate against exact solution
542  ///
543  /// Solution is provided via function pointer.
544  /// Plot error at a given number of plot points.
545  ///
546  //======================================================================
548  std::ostream& outfile,
550  double& error,
551  double& norm)
552  {
553  // Initialise
554  error = 0.0;
555  norm = 0.0;
556 
557  // Vector of local coordinates
558  Vector<double> s(2);
559 
560  // Vector for coordintes
561  Vector<double> x(2);
562 
563  // Find out how many nodes there are in the element
564  unsigned n_node = nnode();
565 
566  Shape psi(n_node);
567 
568  // Set the value of n_intpt
569  unsigned n_intpt = integral_pt()->nweight();
570 
571  // Tecplot
572  outfile << "ZONE" << std::endl;
573 
574  // Exact solution Vector (here a scalar)
575  // Vector<double> exact_soln(1);
576  Vector<double> exact_soln(2);
577 
578  // Loop over the integration points
579  for (unsigned ipt = 0; ipt < n_intpt; ipt++)
580  {
581  // Assign values of s
582  for (unsigned i = 0; i < 2; i++)
583  {
584  s[i] = integral_pt()->knot(ipt, i);
585  }
586 
587  // Get the integral weight
588  double w = integral_pt()->weight(ipt);
589 
590  // Get jacobian of mapping
591  double J = J_eulerian(s);
592 
593  // Premultiply the weights and the Jacobian
594  double W = w * J;
595 
596  // Get x position as Vector
597  interpolated_x(s, x);
598 
599  // Get FE function value
600  double w_fe = interpolated_w_fvk(s);
601 
602  // Get exact solution at this point
603  (*exact_soln_pt)(x, exact_soln);
604 
605  // Output x,y,error
606  for (unsigned i = 0; i < 2; i++)
607  {
608  outfile << x[i] << " ";
609  }
610  outfile << exact_soln[0] << " " << exact_soln[0] - w_fe << std::endl;
611 
612  // Add to error and norm
613  norm += exact_soln[0] * exact_soln[0] * W;
614  error += (exact_soln[0] - w_fe) * (exact_soln[0] - w_fe) * W;
615  }
616  }
617 
618 
619  //====================================================================
620  // Force build of templates
621  //====================================================================
622  template class QFoepplvonKarmanElement<2>;
623  template class QFoepplvonKarmanElement<3>;
624  template class QFoepplvonKarmanElement<4>;
625 
626 } // namespace oomph
static char t char * s
Definition: cfortran.h:568
cstr elem_len * i
Definition: cfortran.h:603
A Class for the derivatives of shape functions The class design is essentially the same as Shape,...
Definition: shape.h:278
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
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 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 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 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
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
double raw_nodal_value(const unsigned &n, const unsigned &i) const
Return the i-th value stored at local node n but do NOT take hanging nodes into account.
Definition: elements.h:2576
double raw_nodal_position(const unsigned &n, const unsigned &i) const
Return the i-th coordinate at local node n. Do not use the hanging node representation....
Definition: elements.cc:1686
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
virtual unsigned self_test()
Self-test: Check inversion of element & do self-test for GeneralisedElement. Return 0 if OK.
Definition: elements.cc:4440
void compute_error(std::ostream &outfile, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt, double &error, double &norm)
Get error against and norm of exact solution.
double interpolated_w_fvk(const Vector< double > &s, unsigned index=0) const
Return FE representation of function value w_fvk(s) at local coordinate s (by default - if index > 0,...
void output_fct(std::ostream &outfile, const unsigned &n_plot, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt)
Output exact soln: x,y,w_exact at n_plot^DIM plot points.
void output(std::ostream &outfile)
Output with default number of plot points.
void fill_in_contribution_to_residuals(Vector< double > &residuals)
Fill in the residuals with this element's contribution.
virtual void get_airy_forcing_fvk(const unsigned &ipt, const Vector< double > &x, double &airy_forcing) const
Get Airy forcing term at (Eulerian) position x. This function is virtual to allow overloading in mult...
int Volume_constraint_pressure_external_data_index
Index of the external Data object that represents the volume constraint pressure (initialised to -1 i...
virtual double get_bounded_volume() const
Return the integral of the displacement over the current element, effectively calculating its contrib...
bool Linear_bending_model
Flag which stores whether we are using a linear, pure bending model instead of the full non-linear Fo...
virtual void get_pressure_fvk(const unsigned &ipt, const Vector< double > &x, double &pressure) const
Get pressure term at (Eulerian) position x. This function is virtual to allow overloading in multi-ph...
virtual unsigned nodal_index_fvk(const unsigned &i=0) const
Return the index at which the i-th unknown value is stored. The default value, i, is appropriate for ...
unsigned self_test()
Self-test: Return 0 for OK.
virtual double dshape_and_dtest_eulerian_at_knot_fvk(const unsigned &ipt, Shape &psi, DShape &dpsidx, Shape &test, DShape &dtestdx) const =0
Shape/test functions and derivs w.r.t. to global coords at integration point ipt; return Jacobian of ...
static double Default_Physical_Constant_Value
Default value for physical constants.
void interpolated_stress(const Vector< double > &s, double &sigma_xx, double &sigma_yy, double &sigma_xy)
Compute in-plane stresses.
int external_local_eqn(const unsigned &i, const unsigned &j)
Return the local equation number corresponding to the j-th value stored at the i-th external data.
Definition: elements.h:311
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.
//////////////////////////////////////////////////////////////////////// ////////////////////////////...
static const unsigned Initial_Nvalue
Static int that holds the number of variables at nodes: always the same.
A Class for shape functions. In simple cases, the shape functions have only one index that can be tho...
Definition: shape.h:76
//////////////////////////////////////////////////////////////////// ////////////////////////////////...