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// Version 1.0; svn revision $LastChangedRevision: 1097 $
7 //LIC//
8 //LIC// $LastChangedDate: 2015-12-17 11:53:17 +0000 (Thu, 17 Dec 2015) $
9 //LIC//
10 //LIC// Copyright (C) 2006-2016 Matthias Heil and Andrew Hazel
11 //LIC//
12 //LIC// This library is free software; you can redistribute it and/or
13 //LIC// modify it under the terms of the GNU Lesser General Public
14 //LIC// License as published by the Free Software Foundation; either
15 //LIC// version 2.1 of the License, or (at your option) any later version.
16 //LIC//
17 //LIC// This library is distributed in the hope that it will be useful,
18 //LIC// but WITHOUT ANY WARRANTY; without even the implied warranty of
19 //LIC// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
20 //LIC// Lesser General Public License for more details.
21 //LIC//
22 //LIC// You should have received a copy of the GNU Lesser General Public
23 //LIC// License along with this library; if not, write to the Free Software
24 //LIC// Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
25 //LIC// 02110-1301 USA.
26 //LIC//
27 //LIC// The authors may be contacted at oomph-lib@maths.man.ac.uk.
28 //LIC//
29 //LIC//====================================================================
30 //Non-inline functions for FoepplvonKarman elements
31 
33 
34 #include <iostream>
35 
36 namespace oomph
37 {
38 
39 
40 //======================================================================
41 /// Set the data for the number of Variables at each node - 8
42 //======================================================================
43 template<unsigned NNODE_1D>
44 const unsigned QFoepplvonKarmanElement<NNODE_1D>::Initial_Nvalue = 8;
45 
46 // Foeppl von Karman equations static data
47 
48 /// Default value physical constants
50 
51 //======================================================================
52 /// Compute contribution to element residual Vector
53 ///
54 /// Pure version without hanging nodes
55 //======================================================================
58 {
59  //Find out how many nodes there are
60  const unsigned n_node = nnode();
61 
62  //Set up memory for the shape and test functions
63  Shape psi(n_node), test(n_node);
64  DShape dpsidx(n_node,2), dtestdx(n_node,2);
65 
66  //Indices at which the unknowns are stored
67  const unsigned w_nodal_index = nodal_index_fvk(0);
68  const unsigned laplacian_w_nodal_index = nodal_index_fvk(1);
69  const unsigned phi_nodal_index = nodal_index_fvk(2);
70  const unsigned laplacian_phi_nodal_index = nodal_index_fvk(3);
71  const unsigned smooth_dwdx_nodal_index = nodal_index_fvk(4);
72  const unsigned smooth_dwdy_nodal_index = nodal_index_fvk(5);
73  const unsigned smooth_dphidx_nodal_index = nodal_index_fvk(6);
74  const unsigned smooth_dphidy_nodal_index = nodal_index_fvk(7);
75 
76  //Set the value of n_intpt
77  const unsigned n_intpt = integral_pt()->nweight();
78 
79  //Integers to store the local equation numbers
80  int local_eqn=0;
81 
82  //Loop over the integration points
83  for(unsigned ipt=0;ipt<n_intpt;ipt++)
84  {
85  //Get the integral weight
86  double w = integral_pt()->weight(ipt);
87 
88  //Call the derivatives of the shape and test functions
89  double J = dshape_and_dtest_eulerian_at_knot_fvk(ipt,psi,dpsidx,
90  test,dtestdx);
91 
92  //Premultiply the weights and the Jacobian
93  double W = w*J;
94 
95  //Allocate and initialise to zero storage for the interpolated values
97 
98  double interpolated_w = 0;
99  double interpolated_laplacian_w = 0;
100  double interpolated_phi = 0;
101  double interpolated_laplacian_phi = 0;
102 
103  Vector<double> interpolated_dwdx(2,0.0);
104  Vector<double> interpolated_dlaplacian_wdx(2,0.0);
105  Vector<double> interpolated_dphidx(2,0.0);
106  Vector<double> interpolated_dlaplacian_phidx(2,0.0);
107 
108  Vector<double> interpolated_smooth_dwdx(2,0.0);
109  Vector<double> interpolated_smooth_dphidx(2,0.0);
110  double interpolated_continuous_d2wdx2 = 0;
111  double interpolated_continuous_d2wdy2 = 0;
112  double interpolated_continuous_d2phidx2 = 0;
113  double interpolated_continuous_d2phidy2 = 0;
114  double interpolated_continuous_d2wdxdy = 0;
115  double interpolated_continuous_d2phidxdy = 0;
116 
117  //Calculate function values and derivatives:
118  //-----------------------------------------
120  // Loop over nodes
121  for(unsigned l=0;l<n_node;l++)
122  {
123  //Get the nodal values
124  nodal_value[0] = raw_nodal_value(l,w_nodal_index);
125  nodal_value[1] = raw_nodal_value(l,laplacian_w_nodal_index);
126 
128  {
129  nodal_value[2] = raw_nodal_value(l,phi_nodal_index);
130  nodal_value[3] = raw_nodal_value(l,laplacian_phi_nodal_index);
131  nodal_value[4] = raw_nodal_value(l,smooth_dwdx_nodal_index);
132  nodal_value[5] = raw_nodal_value(l,smooth_dwdy_nodal_index);
133  nodal_value[6] = raw_nodal_value(l,smooth_dphidx_nodal_index);
134  nodal_value[7] = raw_nodal_value(l,smooth_dphidy_nodal_index);
135  }
136 
137  //Add contributions from current node/shape function
138  interpolated_w += nodal_value[0]*psi(l);
139  interpolated_laplacian_w += nodal_value[1]*psi(l);
140 
142  {
143  interpolated_phi += nodal_value[2]*psi(l);
144  interpolated_laplacian_phi += nodal_value[3]*psi(l);
145 
146  interpolated_smooth_dwdx[0] += nodal_value[4]*psi(l);
147  interpolated_smooth_dwdx[1] += nodal_value[5]*psi(l);
148  interpolated_smooth_dphidx[0] += nodal_value[6]*psi(l);
149  interpolated_smooth_dphidx[1] += nodal_value[7]*psi(l);
150 
151  interpolated_continuous_d2wdx2 += nodal_value[4]*dpsidx(l,0);
152  interpolated_continuous_d2wdy2 += nodal_value[5]*dpsidx(l,1);
153  interpolated_continuous_d2phidx2 += nodal_value[6]*dpsidx(l,0);
154  interpolated_continuous_d2phidy2 += nodal_value[7]*dpsidx(l,1);
155  //mjr CHECK THESE
156  interpolated_continuous_d2wdxdy += 0.5*(nodal_value[4]*dpsidx(l,1)
157  + nodal_value[5]*dpsidx(l,0));
158  interpolated_continuous_d2phidxdy += 0.5*(nodal_value[6]*dpsidx(l,1)
159  + nodal_value[7]*dpsidx(l,0));
160  }
161 
162  // Loop over directions
163  for(unsigned j=0;j<2;j++)
164  {
165  interpolated_x[j] += raw_nodal_position(l,j)*psi(l);
166  interpolated_dwdx[j] += nodal_value[0]*dpsidx(l,j);
167  interpolated_dlaplacian_wdx[j] += nodal_value[1]*dpsidx(l,j);
168 
170  {
171  interpolated_dphidx[j] += nodal_value[2]*dpsidx(l,j);
172  interpolated_dlaplacian_phidx[j] += nodal_value[3]*dpsidx(l,j);
173  }
174  }
175  }
176 
177  //Get pressure function
178  //-------------------
179  double pressure;
180  get_pressure_fvk(ipt,interpolated_x,pressure);
181 
182  double airy_forcing;
183  get_airy_forcing_fvk(ipt,interpolated_x,airy_forcing);
184 
185  // Assemble residuals and Jacobian
186  //--------------------------------
187 
188  // Loop over the test functions
189  for(unsigned l=0;l<n_node;l++)
190  {
191  //Get the local equation
192  local_eqn = nodal_local_eqn(l,w_nodal_index);
193 
194  // IF it's not a boundary condition
195  if(local_eqn >= 0)
196  {
197  residuals[local_eqn] += pressure*test(l)*W;
198 
199  // Reduced order biharmonic operator
200  for(unsigned k=0;k<2;k++)
201  {
202  residuals[local_eqn] += interpolated_dlaplacian_wdx[k]*dtestdx(l,k)*W;
203  }
204 
206  {
207  // Monge-Ampere part
208  residuals[local_eqn] +=
209  eta()*
210  (interpolated_continuous_d2wdx2*interpolated_continuous_d2phidy2
211  + interpolated_continuous_d2wdy2*interpolated_continuous_d2phidx2
212  - 2*interpolated_continuous_d2wdxdy*interpolated_continuous_d2phidxdy)
213  *test(l)*W;
214  }
215  }
216 
217  //Get the local equation
218  local_eqn = nodal_local_eqn(l,laplacian_w_nodal_index);
219 
220  // IF it's not a boundary condition
221  if(local_eqn >= 0)
222  {
223  // The coupled Poisson equations for the biharmonic operator
224  residuals[local_eqn] += interpolated_laplacian_w*test(l)*W;
225 
226  for(unsigned k=0;k<2;k++)
227  {
228  residuals[local_eqn] += interpolated_dwdx[k]*dtestdx(l,k)*W;
229  }
230  }
231 
232  //Get the local equation
233  local_eqn = nodal_local_eqn(l,phi_nodal_index);
234 
235  // IF it's not a boundary condition
236  if(local_eqn >= 0)
237  {
238  residuals[local_eqn] += airy_forcing*test(l)*W;
239 
240  // Reduced order biharmonic operator
241  for(unsigned k=0;k<2;k++)
242  {
243  residuals[local_eqn] +=
244  interpolated_dlaplacian_phidx[k]*dtestdx(l,k)*W;
245  }
246 
247  // Monge-Ampere part
248  residuals[local_eqn] -=
249  (interpolated_continuous_d2wdx2*interpolated_continuous_d2wdy2
250  - interpolated_continuous_d2wdxdy*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)*W;
292  }
293 
294  local_eqn = nodal_local_eqn(l,smooth_dphidy_nodal_index);
295 
296  if(local_eqn >= 0)
297  {
298  residuals[local_eqn] +=
299  (interpolated_dphidx[1] - interpolated_smooth_dphidx[1])*test(l)*W;
300  }
301  }
302 
303  } // End of loop over integration points
304 
305 
306  // Finally: My contribution to the volume constraint equation
307  // (if any). Note this must call get_bounded_volume since the
308  // definition of the bounded volume can be overloaded in derived
309  // elements.
311  {
312  local_eqn=external_local_eqn(
314  if (local_eqn>=0)
315  {
316  residuals[local_eqn]+=get_bounded_volume();
317  }
318  }
319 
320 }
321 
322 /*
323 void FoepplvonKarmanEquations::fill_in_contribution_to_jacobian(Vector<double> &residuals,
324  DenseMatrix<double> &jacobian)
325 {
326  //Add the contribution to the residuals
327  FoepplvonKarmanEquations::fill_in_contribution_to_residuals(residuals);
328  //Allocate storage for the full residuals (residuals of entire element)
329  unsigned n_dof = ndof();
330  Vector<double> full_residuals(n_dof);
331  //Get the residuals for the entire element
332  FoepplvonKarmanEquations::get_residuals(full_residuals);
333  //Calculate the contributions from the internal dofs
334  //(finite-difference the lot by default)
335  fill_in_jacobian_from_internal_by_fd(full_residuals,jacobian,true);
336  //Calculate the contributions from the external dofs
337  //(finite-difference the lot by default)
338  fill_in_jacobian_from_external_by_fd(full_residuals,jacobian,true);
339  //Calculate the contributions from the nodal dofs
340  fill_in_jacobian_from_nodal_by_fd(full_residuals,jacobian);
341 }
342 */
343 
344 
345 //======================================================================
346 /// Self-test: Return 0 for OK
347 //======================================================================
349 {
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 //======================================================================
373 /// Compute in-plane stresses
374 //======================================================================
376  double& sigma_xx,
377  double& sigma_yy,
378  double& sigma_xy)
379 {
380 
381  // No in plane stresses if linear bending
383  {
384  sigma_xx=0.0;
385  sigma_yy=0.0;
386  sigma_xy=0.0;
387  return;
388  }
389 
390  // Get shape fcts and derivs
391  unsigned n_dim=this->dim();
392  unsigned n_node=this->nnode();
393  Shape psi(n_node);
394  DShape dpsidx(n_node,n_dim);
395  dshape_eulerian(s,psi,dpsidx);
396  double interpolated_continuous_d2phidx2 = 0;
397  double interpolated_continuous_d2phidy2 = 0;
398  double interpolated_continuous_d2phidxdy = 0;
399 
400  const unsigned smooth_dphidx_nodal_index = nodal_index_fvk(6);
401  const unsigned smooth_dphidy_nodal_index = nodal_index_fvk(7);
402 
403  // Loop over nodes
404  for(unsigned l=0;l<n_node;l++)
405  {
406  interpolated_continuous_d2phidx2 +=
407  raw_nodal_value(l,smooth_dphidx_nodal_index)*dpsidx(l,0);
408  interpolated_continuous_d2phidy2 +=
409  raw_nodal_value(l,smooth_dphidy_nodal_index)*dpsidx(l,1);
410  interpolated_continuous_d2phidxdy += 0.5*(
411  raw_nodal_value(l,smooth_dphidx_nodal_index)*dpsidx(l,1)+
412  raw_nodal_value(l,smooth_dphidy_nodal_index)*dpsidx(l,0));
413  }
414 
415  sigma_xx=interpolated_continuous_d2phidy2;
416  sigma_yy=interpolated_continuous_d2phidx2;
417  sigma_xy=-interpolated_continuous_d2phidxdy;
418 
419 }
420 
421 
422 
423 //======================================================================
424 /// Output function:
425 ///
426 /// x,y,w
427 ///
428 /// nplot points in each coordinate direction
429 //======================================================================
430 void FoepplvonKarmanEquations::output(std::ostream &outfile,
431  const unsigned &nplot)
432 {
433 
434  //Vector of local coordinates
435  Vector<double> s(2);
436 
437  // Tecplot header info
438  outfile << tecplot_zone_string(nplot);
439 
440  // Loop over plot points
441  unsigned num_plot_points=nplot_points(nplot);
442  for (unsigned iplot=0;iplot<num_plot_points;iplot++)
443  {
444 
445  // Get local coordinates of plot point
446  get_s_plot(iplot,nplot,s);
447 
448  for(unsigned i=0;i<2;i++)
449  {
450  outfile << interpolated_x(s,i) << " ";
451  }
452  outfile << interpolated_w_fvk(s) << std::endl;
453 
454  }
455 
456  // Write tecplot footer (e.g. FE connectivity lists)
457  write_tecplot_zone_footer(outfile,nplot);
458 
459 }
460 
461 
462 //======================================================================
463 /// C-style output function:
464 ///
465 /// x,y,w
466 ///
467 /// nplot points in each coordinate direction
468 //======================================================================
470  const unsigned &nplot)
471 {
472  //Vector of local coordinates
473  Vector<double> s(2);
474 
475  // Tecplot header info
476  fprintf(file_pt,"%s",tecplot_zone_string(nplot).c_str());
477 
478  // Loop over plot points
479  unsigned num_plot_points=nplot_points(nplot);
480  for (unsigned iplot=0;iplot<num_plot_points;iplot++)
481  {
482  // Get local coordinates of plot point
483  get_s_plot(iplot,nplot,s);
484 
485  for(unsigned i=0;i<2;i++)
486  {
487  fprintf(file_pt,"%g ",interpolated_x(s,i));
488  }
489  fprintf(file_pt,"%g \n",interpolated_w_fvk(s));
490  }
491 
492  // Write tecplot footer (e.g. FE connectivity lists)
493  write_tecplot_zone_footer(file_pt,nplot);
494 }
495 
496 
497 
498 //======================================================================
499  /// Output exact solution
500  ///
501  /// Solution is provided via function pointer.
502  /// Plot at a given number of plot points.
503  ///
504  /// x,y,w_exact
505 //======================================================================
506 void FoepplvonKarmanEquations::output_fct(std::ostream &outfile,
507  const unsigned &nplot,
509 {
510  //Vector of local coordinates
511  Vector<double> s(2);
512 
513  // Vector for coordintes
514  Vector<double> x(2);
515 
516  // Tecplot header info
517  outfile << tecplot_zone_string(nplot);
518 
519  // Exact solution Vector (here a scalar)
520  //Vector<double> exact_soln(1);
521  Vector<double> exact_soln(2);
522 
523  // Loop over plot points
524  unsigned num_plot_points=nplot_points(nplot);
525  for (unsigned iplot=0;iplot<num_plot_points;iplot++)
526  {
527 
528  // Get local coordinates of plot point
529  get_s_plot(iplot,nplot,s);
530 
531  // Get x position as Vector
532  interpolated_x(s,x);
533 
534  // Get exact solution at this point
535  (*exact_soln_pt)(x,exact_soln);
536 
537  //Output x,y,w_exact
538  for(unsigned i=0;i<2;i++)
539  {
540  outfile << x[i] << " ";
541  }
542  outfile << exact_soln[0] << std::endl;
543  }
544 
545  // Write tecplot footer (e.g. FE connectivity lists)
546  write_tecplot_zone_footer(outfile,nplot);
547 }
548 
549 
550 
551 
552 //======================================================================
553  /// Validate against exact solution
554  ///
555  /// Solution is provided via function pointer.
556  /// Plot error at a given number of plot points.
557  ///
558 //======================================================================
559 void FoepplvonKarmanEquations::compute_error(std::ostream &outfile,
561  double& error, double& norm)
562 {
563 
564  // Initialise
565  error=0.0;
566  norm=0.0;
567 
568  //Vector of local coordinates
569  Vector<double> s(2);
570 
571  // Vector for coordintes
572  Vector<double> x(2);
573 
574  //Find out how many nodes there are in the element
575  unsigned n_node = nnode();
576 
577  Shape psi(n_node);
578 
579  //Set the value of n_intpt
580  unsigned n_intpt = integral_pt()->nweight();
581 
582  // Tecplot
583  outfile << "ZONE" << std::endl;
584 
585  // Exact solution Vector (here a scalar)
586  //Vector<double> exact_soln(1);
587  Vector<double> exact_soln(2);
588 
589  //Loop over the integration points
590  for(unsigned ipt=0;ipt<n_intpt;ipt++)
591  {
592 
593  //Assign values of s
594  for(unsigned i=0;i<2;i++)
595  {
596  s[i] = integral_pt()->knot(ipt,i);
597  }
598 
599  //Get the integral weight
600  double w = integral_pt()->weight(ipt);
601 
602  // Get jacobian of mapping
603  double J=J_eulerian(s);
604 
605  //Premultiply the weights and the Jacobian
606  double W = w*J;
607 
608  // Get x position as Vector
609  interpolated_x(s,x);
610 
611  // Get FE function value
612  double w_fe=interpolated_w_fvk(s);
613 
614  // Get exact solution at this point
615  (*exact_soln_pt)(x,exact_soln);
616 
617  //Output x,y,error
618  for(unsigned i=0;i<2;i++)
619  {
620  outfile << x[i] << " ";
621  }
622  outfile << exact_soln[0] << " " << exact_soln[0]-w_fe << std::endl;
623 
624  // Add to error and norm
625  norm+=exact_soln[0]*exact_soln[0]*W;
626  error+=(exact_soln[0]-w_fe)*(exact_soln[0]-w_fe)*W;
627 
628  }
629 }
630 
631 
632 //====================================================================
633 // Force build of templates
634 //====================================================================
635 template class QFoepplvonKarmanElement<2>;
636 template class QFoepplvonKarmanElement<3>;
637 template class QFoepplvonKarmanElement<4>;
638 
639 }
640 
void compute_error(std::ostream &outfile, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt, double &error, double &norm)
Get error against and norm of exact solution.
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.
virtual void get_s_plot(const unsigned &i, const unsigned &nplot, Vector< double > &s) const
Get cector of local coordinates of plot point i (when plotting nplot points in each "coordinate direc...
Definition: elements.h:2938
static double Default_Physical_Constant_Value
Default value for physical constants.
bool Linear_bending_model
Flag which stores whether we are using a linear, pure bending model instead of the full non-linear Fo...
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:3225
cstr elem_len * i
Definition: cfortran.h:607
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:2962
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:2976
virtual unsigned self_test()
Self-test: Check inversion of element & do self-test for GeneralisedElement. Return 0 if OK...
Definition: elements.cc:4247
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:2458
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:3981
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. NOTE: Moved to cc file because of a possible compiler bug in gcc (yes, really!). The move to the cc file avoids inlining which appears to cause problems (only) when compiled with gcc and -O3; offensive "illegal read" is in optimised-out section of code and data that is allegedly illegal is readily readable (by other means) just before this function is called so I can't really see how we could possibly be responsible for this...
Definition: elements.cc:1657
virtual double weight(const unsigned &i) const =0
Return weight of i-th integration point.
int Volume_constraint_pressure_external_data_index
Index of the external Data object that represents the volume constraint pressure (initialised to -1 i...
void interpolated_stress(const Vector< double > &s, double &sigma_xx, double &sigma_yy, double &sigma_xy)
Compute in-plane stresses.
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:3841
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...
static char t char * s
Definition: cfortran.h:572
unsigned self_test()
Self-test: Return 0 for OK.
virtual double knot(const unsigned &i, const unsigned &j) const =0
Return local coordinate s[j] of i-th integration point.
Integral *const & integral_pt() const
Return the pointer to the integration scheme (const version)
Definition: elements.h:1900
void(* SteadyExactSolutionFctPt)(const Vector< double > &, Vector< double > &)
Function pointer for function that computes vector-valued steady "exact solution" as ...
Definition: elements.h:1720
void fill_in_contribution_to_residuals(Vector< double > &residuals)
Fill in the residuals with this element's contribution.
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 ...
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:1383
void output(std::ostream &outfile)
Output with default number of plot points.
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 ...
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:2446
virtual unsigned nweight() const =0
Return the number of integration points of the scheme.
unsigned dim() const
Return the spatial dimension of the element, i.e. the number of local coordinates required to paramet...
Definition: elements.h:2470
unsigned nnode() const
Return the number of nodes.
Definition: elements.h:2134
virtual double get_bounded_volume() const
Return the integral of the displacement over the current element, effectively calculating its contrib...
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:2949
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...
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:313
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...