Taxisym_poroelasticity_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// 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 #ifndef OOMPH_TAXISYM_POROELASTICITY_ELEMENTS_HEADER
31 #define OOMPH_TAXISYM_POROELASTICITY_ELEMENTS_HEADER
32 
33 // Config header generated by autoconfig
34 #ifdef HAVE_CONFIG_H
35  #include <oomph-lib-config.h>
36 #endif
37 
39 #include "../generic/Telements.h"
40 
41 namespace oomph
42 {
43 
44  ///=================================================================
45  /// Element which solves the Darcy/linear elasticity equations
46  /// using TElements
47  /// Geometrically the element is always a six noded triangle.
48  /// We use the mid-side nodes to store edge-based flux degrees of
49  /// freedom and internal data for the discontinuous pressure
50  /// and internal flux dofs.
51  ///=================================================================
52  template<unsigned ORDER>
54  public TElement<2,3>,
56  {
57  private:
58 
59  /// The number of values stored at each node
60  static const unsigned Initial_Nvalue[];
61 
62  /// Face index associated with edge flux degree of freedom
63  static const unsigned Face_index_of_edge_flux[];
64 
65  /// \short Conversion scheme from an edge degree of freedom to the node
66  /// it's stored at
67  static const unsigned Q_edge_conv[];
68 
69  /// The points along each edge where the fluxes are taken to be
70  static const double Flux_interpolation_point[];
71 
72  /// The internal data index where the internal q degrees of freedom are stored
74 
75  /// The internal data index where the p degrees of freedom are stored
77 
78  /// \short Unit normal signs associated with each edge to ensure inter-element
79  /// continuity of the flux
80  std::vector<short> Sign_edge;
81 
82  public:
83 
84  /// Constructor
86 
87  /// Destructor
89 
90  /// Number of values required at node n
91  unsigned required_nvalue(const unsigned &n) const
92  {
93  return Initial_Nvalue[n];
94  }
95 
96  /// Return the face index associated with specified edge
97  unsigned face_index_of_edge(const unsigned& j) const
98  {
99  return (j+2)%3;
100  }
101 
102  /// \short Compute the face element coordinates of the nth flux interpolation
103  /// point along specified edge
105  const unsigned &n,
106  Vector<double> &s)
107  const
108  {
109  // Get the location of the n-th flux interpolation point along
110  // the edge in terms of the distance along the edge itself
111  Vector<double> flux_interpolation_point=
113 
114  // Convert the edge number to the number of the mid-edge node along that
115  // edge
116  unsigned node_number=Q_edge_conv[edge];
117 
118  // The edge basis functions are defined in a clockwise manner, so we have
119  // to effectively "flip" some coordinates
120  switch(node_number)
121  {
122  case 3:
123  s[0]=flux_interpolation_point[0];
124  break;
125  case 4:
126  s[0]=1.0-flux_interpolation_point[0];
127  break;
128  case 5:
129  s[0]=flux_interpolation_point[0];
130  break;
131  }
132  }
133 
134  /// Return the face index associated with j-th edge flux degree of freedom
135  unsigned face_index_of_q_edge_basis_fct(const unsigned& j) const
136  {
137  return Face_index_of_edge_flux[j];
138  }
139 
140  /// \short Return the nodal index of the j-th solid displacement unknown
141  /// [0: r; 1: z]
142  unsigned u_index_axisym_poroelasticity(const unsigned &j) const
143  {
144 #ifdef RANGE_CHECKING
145  if(j >= 2)
146  {
147  std::ostringstream error_message;
148  error_message << "Range Error: j " << j
149  << " is not in the range (0,1)";
150  throw OomphLibError(error_message.str(),
151  OOMPH_CURRENT_FUNCTION,
152  OOMPH_EXCEPTION_LOCATION);
153  }
154 #endif
155  return j;
156  }
157 
158 
159  /// Return the equation number of the j-th edge (flux) degree of freedom
160  int q_edge_local_eqn(const unsigned &j) const
161  {
162 #ifdef RANGE_CHECKING
163  if(j >= nq_basis_edge())
164  {
165  std::ostringstream error_message;
166  error_message << "Range Error: j " << j
167  << " is not in the range (0,"
168  << nq_basis_edge()-1 <<")";
169  throw OomphLibError(error_message.str(),
170  OOMPH_CURRENT_FUNCTION,
171  OOMPH_EXCEPTION_LOCATION);
172  }
173 #endif
174  return this->nodal_local_eqn(q_edge_node_number(j),q_edge_index(j));
175  }
176 
177  /// Return the equation number of the j-th internal degree of freedom
178  int q_internal_local_eqn(const unsigned &j) const
179  {
180 #ifdef RANGE_CHECKING
181  if(j >= nq_basis_internal())
182  {
183  std::ostringstream error_message;
184  error_message << "Range Error: j " << j
185  << " is not in the range (0,"
186  << nq_basis_internal()-1 <<")";
187  throw OomphLibError(error_message.str(),
188  OOMPH_CURRENT_FUNCTION,
189  OOMPH_EXCEPTION_LOCATION);
190  }
191 #endif
193  }
194 
195  /// \short Return vector of pointers to the Data objects that store the
196  /// edge flux values
198  {
199  // It's the mid-side nodes:
200  Vector<Data*> data_pt(3);
201  data_pt[0]=node_pt(3);
202  data_pt[1]=node_pt(4);
203  data_pt[2]=node_pt(5);
204  return data_pt;
205  }
206 
207  /// Return pointer to the Data object that stores the internal flux values
209  {
211  }
212 
213  /// \short Return the index of the internal data where the q_internal degrees
214  /// of freedom are stored
215  unsigned q_internal_index() const
216  {
217  return Q_internal_data_index;
218  }
219 
220  /// Return the nodal index at which the jth edge unknown is stored
221  unsigned q_edge_index(const unsigned &j) const
222  {
223 #ifdef RANGE_CHECKING
224  if(j >= (nq_basis_edge()))
225  {
226  std::ostringstream error_message;
227  error_message << "Range Error: j " << j
228  << " is not in the range (0,"
229  << nq_basis_edge()-1 <<")";
230  throw OomphLibError(error_message.str(),
231  OOMPH_CURRENT_FUNCTION,
232  OOMPH_EXCEPTION_LOCATION);
233  }
234 #endif
235  return j%(ORDER+1)+2;
236  }
237 
238  /// Return the number of the node where the jth edge unknown is stored
239  unsigned q_edge_node_number(const unsigned &j) const
240  {
241 #ifdef RANGE_CHECKING
242  if(j >= (nq_basis_edge()))
243  {
244  std::ostringstream error_message;
245  error_message << "Range Error: j " << j
246  << " is not in the range (0,"
247  << nq_basis_edge()-1 <<")";
248  throw OomphLibError(error_message.str(),
249  OOMPH_CURRENT_FUNCTION,
250  OOMPH_EXCEPTION_LOCATION);
251  }
252 #endif
253  return Q_edge_conv[j/(ORDER+1)];
254  }
255 
256 
257  /// Get pointer to node that stores the edge flux dofs for specified edge
258  Node* edge_flux_node_pt(const unsigned& edge)
259  {
260  return node_pt(Q_edge_conv[edge]);
261  }
262 
263  /// Return the values of the j-th edge (flux) degree of freedom
264  double q_edge(const unsigned &j) const
265  {
266 #ifdef RANGE_CHECKING
267  if(j >= (nq_basis_edge()))
268  {
269  std::ostringstream error_message;
270  error_message << "Range Error: j " << j
271  << " is not in the range (0,"
272  << nq_basis_edge()-1 <<")";
273  throw OomphLibError(error_message.str(),
274  OOMPH_CURRENT_FUNCTION,
275  OOMPH_EXCEPTION_LOCATION);
276  }
277 #endif
279  }
280 
281  /// \short Return the values of the j-th edge (flux) degree of
282  /// freedom at time history level t
283  double q_edge(const unsigned &t, const unsigned &j) const
284  {
285 #ifdef RANGE_CHECKING
286  if(j >= (nq_basis_edge()))
287  {
288  std::ostringstream error_message;
289  error_message << "Range Error: j " << j
290  << " is not in the range (0,"
291  << nq_basis_edge()-1 <<")";
292  throw OomphLibError(error_message.str(),
293  OOMPH_CURRENT_FUNCTION,
294  OOMPH_EXCEPTION_LOCATION);
295  }
296 #endif
298  }
299 
300  /// Return the values of the internal degree of freedom
301  double q_internal(const unsigned &j) const
302  {
303 #ifdef RANGE_CHECKING
304  if(j >= nq_basis_internal())
305  {
306  std::ostringstream error_message;
307  error_message << "Range Error: j " << j
308  << " is not in the range (0,"
309  << nq_basis_internal()-1 <<")";
310  throw OomphLibError(error_message.str(),
311  OOMPH_CURRENT_FUNCTION,
312  OOMPH_EXCEPTION_LOCATION);
313  }
314 #endif
315  return this->internal_data_pt(q_internal_index())->value(j);
316  }
317 
318  /// \short Return the value of the j-th internal degree of freedom at
319  /// time history level t
320  double q_internal(const unsigned &t,const unsigned &j) const
321  {
322 #ifdef RANGE_CHECKING
323  if(j >= nq_basis_internal())
324  {
325  std::ostringstream error_message;
326  error_message << "Range Error: j " << j
327  << " is not in the range (0,"
328  << nq_basis_internal()-1 <<")";
329  throw OomphLibError(error_message.str(),
330  OOMPH_CURRENT_FUNCTION,
331  OOMPH_EXCEPTION_LOCATION);
332  }
333 #endif
334  return this->internal_data_pt(q_internal_index())->value(t,j);
335  }
336 
337  /// Pin the j-th edge (flux) degree of freedom and set it to specified value
338  void pin_q_edge_value(const unsigned &j, const double& value)
339  {
342  }
343 
344  /// Set the values of the j-th edge (flux) degree of freedom
345  void set_q_edge(const unsigned &j, const double& value)
346  {
348  }
349 
350  /// \short Set the values of the j-th edge (flux) degree of freedom at
351  /// time history level t
352  void set_q_edge(const unsigned &t, const unsigned &j, const double& value)
353  {
355  }
356 
357 
358  /// \short Set the values of the j-th internal degree of freedom
359  void set_q_internal(const unsigned &j, const double& value)
360  {
361  this->internal_data_pt(q_internal_index())->set_value(j,value);
362  }
363 
364  /// \short Set the values of the j-th internal degree of freedom at
365  /// time history level t
366  void set_q_internal(const unsigned &t, const unsigned &j, const double& value)
367  {
368  this->internal_data_pt(q_internal_index())->set_value(t,j,value);
369  }
370 
371  /// Return the number of edge basis functions for flux q
372  unsigned nq_basis_edge() const;
373 
374  /// Return the number of internal basis functions for flux q
375  unsigned nq_basis_internal() const;
376 
377  /// Returns the local form of the q basis at local coordinate s
378  void get_q_basis_local(const Vector<double> &s,
379  Shape &q_basis) const;
380 
381  /// Returns the local form of the q basis and dbasis/ds at local coordinate s
383  Shape &div_q_basis_ds) const;
384 
385  /// \short Returns the number of flux_interpolation points along each
386  /// edge of the element
387  unsigned nedge_flux_interpolation_point() const;
388 
389  /// \short Returns the local coordinate of the jth flux_interpolation point
390  /// along specified edge
392  const unsigned &j) const
393  {
394 #ifdef RANGE_CHECKING
395  if(edge >= 3)
396  {
397  std::ostringstream error_message;
398  error_message << "Range Error: edge " << edge
399  << " is not in the range (0,2)";
400  throw OomphLibError(error_message.str(),
401  OOMPH_CURRENT_FUNCTION,
402  OOMPH_EXCEPTION_LOCATION);
403  }
405  {
406  std::ostringstream error_message;
407  error_message << "Range Error: j " << j
408  << " is not in the range (0,"
409  << nedge_flux_interpolation_point()-1 <<")";
410  throw OomphLibError(error_message.str(),
411  OOMPH_CURRENT_FUNCTION,
412  OOMPH_EXCEPTION_LOCATION);
413  }
414 #endif
415  Vector<double> coord(1);
416  coord[0]=(1.0-sign_edge(edge))/2.0+sign_edge(edge)*
418  return coord;
419  }
420 
421  /// \short Compute the global coordinates of the jth flux_interpolation
422  /// point along specified edge
423  void edge_flux_interpolation_point_global(const unsigned &edge,
424  const unsigned &j,
425  Vector<double> &x) const
426  {
427 #ifdef RANGE_CHECKING
428  if(edge >= 3)
429  {
430  std::ostringstream error_message;
431  error_message << "Range Error: edge " << edge
432  << " is not in the range (0,2)";
433  throw OomphLibError(error_message.str(),
434  OOMPH_CURRENT_FUNCTION,
435  OOMPH_EXCEPTION_LOCATION);
436  }
438  {
439  std::ostringstream error_message;
440  error_message << "Range Error: j " << j
441  << " is not in the range (0,"
442  << nedge_flux_interpolation_point()-1 <<")";
443  throw OomphLibError(error_message.str(),
444  OOMPH_CURRENT_FUNCTION,
445  OOMPH_EXCEPTION_LOCATION);
446  }
447 #endif
448 
449  // Get the location of the n-th flux_interpolation point along the
450  // edge in terms of the distance along the edge itself
451  Vector<double> flux_interpolation_point=
453 
454  // Convert the edge number to the number of the mid-edge node along that
455  // edge
456  unsigned node_number=Q_edge_conv[edge];
457 
458  // Storage for the local coords of the flux_interpolation point
459  Vector<double> s_flux_interpolation(2,0);
460 
461  // The edge basis functions are defined in a clockwise manner, so we have
462  // to effectively "flip" the coordinates along edges 0 and 1 to match this
463  switch(node_number)
464  {
465  case 3:
466  s_flux_interpolation[0]=1.0-flux_interpolation_point[0];
467  s_flux_interpolation[1]=flux_interpolation_point[0];
468  break;
469  case 4:
470  s_flux_interpolation[0]=0.0;
471  s_flux_interpolation[1]=1.0-flux_interpolation_point[0];
472  break;
473  case 5:
474  s_flux_interpolation[0]=flux_interpolation_point[0];
475  s_flux_interpolation[1]=0.0;
476  break;
477  }
478 
479  // Calculate the global coordinates from the local ones
480  interpolated_x(s_flux_interpolation,x);
481  }
482 
483  /// Pin the jth internal q value and set it to specified value
484  void pin_q_internal_value(const unsigned &j, const double& q)
485  {
486 #ifdef RANGE_CHECKING
487  if(j >= nq_basis_internal())
488  {
489  std::ostringstream error_message;
490  error_message << "Range Error: j " << j
491  << " is not in the range (0,"
492  << nq_basis_internal()-1 <<")";
493  throw OomphLibError(error_message.str(),
494  OOMPH_CURRENT_FUNCTION,
495  OOMPH_EXCEPTION_LOCATION);
496  }
497 #endif
500  }
501 
502  /// Return the equation number of the j-th pressure degree of freedom
503  int p_local_eqn(const unsigned &j) const
504  {
505 #ifdef RANGE_CHECKING
506  if(j >= np_basis())
507  {
508  std::ostringstream error_message;
509  error_message << "Range Error: j " << j
510  << " is not in the range (0,"
511  << np_basis()-1 <<")";
512  throw OomphLibError(error_message.str(),
513  OOMPH_CURRENT_FUNCTION,
514  OOMPH_EXCEPTION_LOCATION);
515  }
516 #endif
518  }
519 
520  /// Return the jth pressure value
521  double p_value(const unsigned &j) const
522  {
523 #ifdef RANGE_CHECKING
524  if(j >= np_basis())
525  {
526  std::ostringstream error_message;
527  error_message << "Range Error: j " << j
528  << " is not in the range (0,"
529  << np_basis()-1 <<")";
530  throw OomphLibError(error_message.str(),
531  OOMPH_CURRENT_FUNCTION,
532  OOMPH_EXCEPTION_LOCATION);
533  }
534 #endif
535  return this->internal_data_pt(P_internal_data_index)->value(j);
536  }
537 
538  /// Return the total number of pressure basis functions
539  unsigned np_basis() const;
540 
541  /// Return the pressure basis
542  void get_p_basis(const Vector<double> &s,
543  Shape &p_basis) const;
544 
545  /// Pin the jth pressure value and set to specified value
546  void pin_p_value(const unsigned &j, const double &p)
547  {
548 #ifdef RANGE_CHECKING
549  if(j >= np_basis())
550  {
551  std::ostringstream error_message;
552  error_message << "Range Error: j " << j
553  << " is not in the range (0,"
554  << np_basis()-1 <<")";
555  throw OomphLibError(error_message.str(),
556  OOMPH_CURRENT_FUNCTION,
557  OOMPH_EXCEPTION_LOCATION);
558  }
559 #endif
562  }
563 
564  /// Return pointer to the Data object that stores the pressure values
565  Data* p_data_pt() const
566  {
568  }
569 
570  /// Set the jth pressure value
571  void set_p_value(const unsigned &j, const double& value)
572  {
574  }
575 
576  /// Scale the edge basis to allow arbitrary edge mappings
577  void scale_basis(Shape& basis) const
578  {
579  // Storage for the lengths of the edges of the element
580  Vector<double> length(3,0.0);
581 
582  // Temporary storage for the vertex positions
583  double x0,y0,x1,y1;
584 
585  // loop over the edges of the element and calculate their lengths (in x-y
586  // space)
587  for(unsigned i=0;i<3;i++)
588  {
589  x0=this->node_pt(i)->x(0);
590  y0=this->node_pt(i)->x(1);
591  x1=this->node_pt((i+1)%3)->x(0);
592  y1=this->node_pt((i+1)%3)->x(1);
593 
594  length[i]=std::sqrt(std::pow(y1-y0,2)+std::pow(x1-x0,2));
595  }
596 
597  // lengths of the sides of the reference element (in the same order as the
598  // basis functions)
599  const double ref_length[3]={std::sqrt(2.0),1,1};
600 
601  // get the number of basis functions associated with the edges
602  unsigned n_q_basis_edge=nq_basis_edge();
603 
604  // rescale the edge basis functions to allow arbitrary edge mappings from
605  // element to ref. element
606  const unsigned n_index2 = basis.nindex2();
607  for(unsigned i=0;i<n_index2;i++)
608  {
609  for(unsigned l=0;l<n_q_basis_edge;l++)
610  {
611  basis(l,i)*=(length[l/(ORDER+1)]/ref_length[l/(ORDER+1)]);
612  }
613  }
614  }
615 
616  /// Accessor for the unit normal sign of edge n (const version)
617  const short &sign_edge(const unsigned &n) const
618  {
619  return Sign_edge[n];
620  }
621 
622  /// Accessor for the unit normal sign of edge n
623  short &sign_edge(const unsigned &n)
624  {
625  return Sign_edge[n];
626  }
627 
628  /// Output with default number of plot points
629  void output(std::ostream &outfile)
630  {
632  }
633 
634  /// \short Output FE representation of soln: x,y,u1,u2,div_q,p at
635  /// Nplot^DIM plot points
636  void output(std::ostream &outfile, const unsigned &Nplot)
637  {
639  }
640 
641 
642  /// \short Number of vertex nodes in the element
643  unsigned nvertex_node() const
644  {return TElement<2,3>::nvertex_node();}
645 
646  /// \short Pointer to the j-th vertex node in the element
647  Node* vertex_node_pt(const unsigned& j) const
648  {return TElement<2,3>::vertex_node_pt(j);}
649 
650 /// Recovery order for Z2 error estimator
651  unsigned nrecovery_order()
652  {
653  return 2; // need to experiment with this...
654  }
655 
656  protected:
657 
658  /// \short Returns the geometric basis, and the u, p and divergence basis
659  /// functions and test functions at local coordinate s
661  Shape &psi,
662  DShape &dpsi,
663  Shape &u_basis,
664  Shape &u_test,
665  DShape &du_basis_dx,
666  DShape &du_test_dx,
667  Shape &q_basis,
668  Shape &q_test,
669  Shape &p_basis,
670  Shape &p_test,
671  Shape &div_q_basis_ds,
672  Shape &div_q_test_ds) const
673  {
674  const unsigned n_q_basis = this->nq_basis();
675 
676  Shape q_basis_local(n_q_basis,2);
677  this->get_q_basis_local(s,q_basis_local);
678  this->get_p_basis(s,p_basis);
679  this->get_div_q_basis_local(s,div_q_basis_ds);
680 
681  double J = this->transform_basis(s,q_basis_local,psi,dpsi,q_basis);
682 
683  // u_basis consists of the normal Lagrangian shape functions
684  u_basis=psi;
685  du_basis_dx=dpsi;
686 
687  u_test=psi;
688  du_test_dx=dpsi;
689 
690  q_test=q_basis;
691  p_test=p_basis;
692  div_q_test_ds=div_q_basis_ds;
693 
694  return J;
695  }
696 
697  /// \short Returns the geometric basis, and the u, p and divergence basis
698  /// functions and test functions at integration point ipt
699  double shape_basis_test_local_at_knot(const unsigned &ipt,
700  Shape &psi,
701  DShape &dpsi,
702  Shape &u_basis,
703  Shape &u_test,
704  DShape &du_basis_dx,
705  DShape &du_test_dx,
706  Shape &q_basis,
707  Shape &q_test,
708  Shape &p_basis,
709  Shape &p_test,
710  Shape &div_q_basis_ds,
711  Shape &div_q_test_ds) const
712  {
713  Vector<double> s(2);
714  for(unsigned i=0;i<2;i++) { s[i] = this->integral_pt()->knot(ipt,i); }
715 
716  return shape_basis_test_local(
717  s,psi,dpsi,
718  u_basis,u_test,
719  du_basis_dx,du_test_dx,
720  q_basis,q_test,
721  p_basis,p_test,
722  div_q_basis_ds,div_q_test_ds);
723  }
724  };
725 
726 
727 
728  //================================================================
729  /// Face geometry for TAxisymmetricPoroelasticityElement<0>
730  //================================================================
731  template<>
733  public virtual TElement<1,3>
734  {
735  public:
736 
737  /// Constructor: Call constructor of base
738  FaceGeometry() : TElement<1,3>() {}
739  };
740 
741 
742  //================================================================
743  /// Face geometry for TAxisymmetricPoroelasticityElement<1>
744  //================================================================
745  template<>
747  public virtual TElement<1,3>
748  {
749  public:
750 
751  /// Constructor: Call constructor of base class
752  FaceGeometry() : TElement<1,3>() {}
753  };
754 
755 
756 
757 } // End of oomph namespace
758 
759 #endif
760 
void set_p_value(const unsigned &j, const double &value)
Set the jth pressure value.
void pin_q_internal_value(const unsigned &j, const double &q)
Pin the jth internal q value and set it to specified value.
int q_edge_local_eqn(const unsigned &j) const
Return the equation number of the j-th edge (flux) degree of freedom.
int internal_local_eqn(const unsigned &i, const unsigned &j) const
Return the local equation number corresponding to the j-th value stored at the i-th internal data...
Definition: elements.h:268
static const unsigned Face_index_of_edge_flux[]
Face index associated with edge flux degree of freedom.
static const double Flux_interpolation_point[]
The points along each edge where the fluxes are taken to be.
double shape_basis_test_local_at_knot(const unsigned &ipt, Shape &psi, DShape &dpsi, Shape &u_basis, Shape &u_test, DShape &du_basis_dx, DShape &du_test_dx, Shape &q_basis, Shape &q_test, Shape &p_basis, Shape &p_test, Shape &div_q_basis_ds, Shape &div_q_test_ds) const
Returns the geometric basis, and the u, p and divergence basis functions and test functions at integr...
double q_edge(const unsigned &j) const
Return the values of the j-th edge (flux) degree of freedom.
unsigned u_index_axisym_poroelasticity(const unsigned &j) const
Return the nodal index of the j-th solid displacement unknown [0: r; 1: z].
double transform_basis(const Vector< double > &s, const Shape &q_basis_local, Shape &psi, DShape &dpsi, Shape &q_basis) const
Performs a div-conserving transformation of the vector basis functions from the reference element to ...
unsigned nedge_flux_interpolation_point() const
Returns the number of flux_interpolation points along each edge of the element.
void scale_basis(Shape &basis) const
Scale the edge basis to allow arbitrary edge mappings.
void get_q_basis_local(const Vector< double > &s, Shape &q_basis) const
Returns the local form of the q basis at local coordinate s.
cstr elem_len * i
Definition: cfortran.h:607
void get_p_basis(const Vector< double > &s, Shape &p_basis) const
Return the pressure basis.
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
unsigned np_basis() const
Return the total number of pressure basis functions.
char t
Definition: cfortran.h:572
Nodes are derived from Data, but, in addition, have a definite (Eulerian) position in a space of a gi...
Definition: nodes.h:852
short & sign_edge(const unsigned &n)
Accessor for the unit normal sign of edge n.
unsigned q_edge_node_number(const unsigned &j) const
Return the number of the node where the jth edge unknown is stored.
void output(std::ostream &outfile)
Output with default number of plot points.
Node * vertex_node_pt(const unsigned &j) const
Pointer to the j-th vertex node in the element.
unsigned nrecovery_order()
Recovery order for Z2 error estimator.
const short & sign_edge(const unsigned &n) const
Accessor for the unit normal sign of edge n (const version)
Data * p_data_pt() const
Return pointer to the Data object that stores the pressure values.
void pin(const unsigned &i)
Pin the i-th stored variable.
Definition: nodes.h:383
unsigned q_edge_index(const unsigned &j) const
Return the nodal index at which the jth edge unknown is stored.
unsigned Q_internal_data_index
The internal data index where the internal q degrees of freedom are stored.
Data * q_internal_data_pt() const
Return pointer to the Data object that stores the internal flux values.
Vector< double > edge_flux_interpolation_point(const unsigned &edge, const unsigned &j) const
Returns the local coordinate of the jth flux_interpolation point along specified edge.
void pin_q_edge_value(const unsigned &j, const double &value)
Pin the j-th edge (flux) degree of freedom and set it to specified value.
void set_q_edge(const unsigned &j, const double &value)
Set the values of the j-th edge (flux) degree of freedom.
void set_q_internal(const unsigned &j, const double &value)
Set the values of the j-th internal degree of freedom.
double & x(const unsigned &i)
Return the i-th nodal coordinate.
Definition: nodes.h:995
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
unsigned nvertex_node() const
Number of vertex nodes in the element.
double shape_basis_test_local(const Vector< double > &s, Shape &psi, DShape &dpsi, Shape &u_basis, Shape &u_test, DShape &du_basis_dx, DShape &du_test_dx, Shape &q_basis, Shape &q_test, Shape &p_basis, Shape &p_test, Shape &div_q_basis_ds, Shape &div_q_test_ds) const
Returns the geometric basis, and the u, p and divergence basis functions and test functions at local ...
void get_div_q_basis_local(const Vector< double > &s, Shape &div_q_basis_ds) const
Returns the local form of the q basis and dbasis/ds at local coordinate s.
unsigned q_internal_index() const
Return the index of the internal data where the q_internal degrees of freedom are stored...
unsigned nq_basis_internal() const
Return the number of internal basis functions for flux q.
static char t char * s
Definition: cfortran.h:572
void output(std::ostream &outfile, const unsigned &Nplot)
Output FE representation of soln: x,y,u1,u2,div_q,p at Nplot^DIM plot points.
virtual double knot(const unsigned &i, const unsigned &j) const =0
Return local coordinate s[j] of i-th integration point.
void set_value(const unsigned &i, const double &value_)
Set the i-th stored data value to specified value. The only reason that we require an explicit set fu...
Definition: nodes.h:267
Data *& internal_data_pt(const unsigned &i)
Return a pointer to i-th internal data object.
Definition: elements.h:623
Integral *const & integral_pt() const
Return the pointer to the integration scheme (const version)
Definition: elements.h:1900
A class that represents a collection of data; each Data object may contain many different individual ...
Definition: nodes.h:89
void edge_flux_interpolation_point_global(const unsigned &edge, const unsigned &j, Vector< double > &x) const
Compute the global coordinates of the jth flux_interpolation point along specified edge...
double q_internal(const unsigned &j) const
Return the values of the internal degree of freedom.
Node *& node_pt(const unsigned &n)
Return a pointer to the local node n.
Definition: elements.h:2097
void set_q_internal(const unsigned &t, const unsigned &j, const double &value)
Set the values of the j-th internal degree of freedom at time history level t.
Node * edge_flux_node_pt(const unsigned &edge)
Get pointer to node that stores the edge flux dofs for specified edge.
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
Vector< Data * > q_edge_data_pt() const
Return vector of pointers to the Data objects that store the edge flux values.
std::vector< short > Sign_edge
Unit normal signs associated with each edge to ensure inter-element continuity of the flux...
unsigned P_internal_data_index
The internal data index where the p degrees of freedom are stored.
void pin_p_value(const unsigned &j, const double &p)
Pin the jth pressure value and set to specified value.
unsigned nq_basis_edge() const
Return the number of edge basis functions for flux q.
virtual unsigned nq_basis() const
Return the total number of computational basis functions for q.
double q_internal(const unsigned &t, const unsigned &j) const
Return the value of the j-th internal degree of freedom at time history level t.
unsigned face_index_of_edge(const unsigned &j) const
Return the face index associated with specified edge.
double p_value(const unsigned &j) const
Return the jth pressure value.
void set_q_edge(const unsigned &t, const unsigned &j, const double &value)
Set the values of the j-th edge (flux) degree of freedom at time history level t. ...
Class implementing the generic maths of the axisym poroelasticity equations: axisym linear elasticity...
double q_edge(const unsigned &t, const unsigned &j) const
Return the values of the j-th edge (flux) degree of freedom at time history level t...
int q_internal_local_eqn(const unsigned &j) const
Return the equation number of the j-th internal degree of freedom.
int p_local_eqn(const unsigned &j) const
Return the equation number of the j-th pressure degree of freedom.
double value(const unsigned &i) const
Return i-th stored value. This function is not virtual so that it can be inlined. This means that if ...
Definition: nodes.h:291
unsigned face_index_of_q_edge_basis_fct(const unsigned &j) const
Return the face index associated with j-th edge flux degree of freedom.
void output(std::ostream &outfile)
Output with default number of plot points.
unsigned required_nvalue(const unsigned &n) const
Number of values required at node n.
void face_local_coordinate_of_flux_interpolation_point(const unsigned &edge, const unsigned &n, Vector< double > &s) const
Compute the face element coordinates of the nth flux interpolation point along specified edge...
unsigned nindex2() const
Return the range of index 2 of the shape function object.
Definition: shape.h:265