mixed_order_petrov_galerkin_space_time_fluid_traction_elements.h
Go to the documentation of this file.
1// LIC// ====================================================================
2// LIC// This file forms part of oomph-lib, the object-oriented,
3// LIC// multi-physics finite-element library, available
4// LIC// at http://www.oomph-lib.org.
5// LIC//
6// LIC// Copyright (C) 2006-2022 Matthias Heil and Andrew Hazel
7// LIC//
8// LIC// This library is free software; you can redistribute it and/or
9// LIC// modify it under the terms of the GNU Lesser General Public
10// LIC// License as published by the Free Software Foundation; either
11// LIC// version 2.1 of the License, or (at your option) any later version.
12// LIC//
13// LIC// This library is distributed in the hope that it will be useful,
14// LIC// but WITHOUT ANY WARRANTY; without even the implied warranty of
15// LIC// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
16// LIC// Lesser General Public License for more details.
17// LIC//
18// LIC// You should have received a copy of the GNU Lesser General Public
19// LIC// License along with this library; if not, write to the Free Software
20// LIC// Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
21// LIC// 02110-1301 USA.
22// LIC//
23// LIC// The authors may be contacted at oomph-lib@maths.man.ac.uk.
24// LIC//
25// LIC//====================================================================
26// Header file for elements that are used to integrate fluid tractions
27// This includes the guts (i.e. equations) because we want to inline them
28// for faster operation, although it slows down the compilation!
29#ifndef OOMPH_MIXED_ORDER_PETROV_GALERKIN_SPACE_TIME_FLUID_TRACTION_ELEMENTS_HEADER
30#define OOMPH_MIXED_ORDER_PETROV_GALERKIN_SPACE_TIME_FLUID_TRACTION_ELEMENTS_HEADER
31
32// Config header generated by autoconfig
33#ifdef HAVE_CONFIG_H
34#include <oomph-lib-config.h>
35#endif
36
37// OOMPH-LIB headers
38#include "generic/Qelements.h"
39#include "generic/Telements.h"
40
41namespace oomph
42{
43 //======================================================================
44 /// A class for elements that allow the imposition of an applied traction
45 /// to the Navier--Stokes equations
46 /// The geometrical information can be read from the FaceGeometry<ELEMENT>
47 /// class and and thus, we can be generic enough without the need to have
48 /// a separate equations class
49 //======================================================================
50 template<class ELEMENT>
52 : public virtual FaceGeometry<ELEMENT>,
53 public virtual FaceElement
54 {
55 private:
56 /// Pointer to an imposed traction function
57 void (*Traction_fct_pt)(const double& time,
58 const Vector<double>& x,
59 const Vector<double>& n,
60 Vector<double>& result);
61
62 public:
63 /// Constructor, which takes a "bulk" element and the value of the index
64 /// and its limit
66 FiniteElement* const& element_pt,
67 const int& face_index,
68 const bool& called_from_refineable_constructor = false)
69 : FaceGeometry<ELEMENT>(), FaceElement()
70 {
71 // Attach the geometrical information to the element. N.B. This function
72 // also assigns nbulk_value from the required_nvalue of the bulk element
73 element_pt->build_face_element(face_index, this);
74
75#ifdef PARANOID
76 {
77 // Check that the element is not a refineable 3D element
78 if (!called_from_refineable_constructor)
79 {
80 // If it's 3D
81 if (element_pt->dim() == 3)
82 {
83 // Is it refineable
84 RefineableElement* ref_el_pt =
85 dynamic_cast<RefineableElement*>(element_pt);
86
87 // If the upcast was successful
88 if (ref_el_pt != 0)
89 {
90 // If there are any hanging nodes
91 if (this->has_hanging_nodes())
92 {
93 // Create an output stream
94 std::ostringstream error_message_stream;
95
96 // Create an error message
97 error_message_stream
98 << "This flux element will not work correctly "
99 << "if nodes are hanging" << std::endl;
100
101 // Throw the error message
102 throw OomphLibError(error_message_stream.str(),
103 OOMPH_CURRENT_FUNCTION,
104 OOMPH_EXCEPTION_LOCATION);
105 }
106 } // if (ref_el_pt!=0)
107 } // if (element_pt->dim()==3)
108 } // if (!called_from_refineable_constructor)
109 }
110#endif
111
112 // Set the body force function pointer to zero
113 Traction_fct_pt = 0;
114
115 // Set the dimension from the dimension of the first node
116 Dim = this->node_pt(0)->ndim();
117 } // End of NavierStokesMixedOrderSpaceTimeTractionElement
118
119
120 /// Destructor should not delete anything
122
123
124 /// Access function for the imposed traction pointer
125 void (*&traction_fct_pt())(const double& t,
126 const Vector<double>& x,
127 const Vector<double>& n,
128 Vector<double>& result)
129 {
130 // Return the function pointer
131 return Traction_fct_pt;
132 } // End of traction_fct_pt
133
134
135 /// This function returns just the residuals
137 {
138 // Call the generic residuals function with flag set to 0
139 // using a dummy matrix argument
142 }
143
144
145 /// This function returns the residuals and the jacobian
147 DenseMatrix<double>& jacobian)
148 {
149 // Call the generic routine with the flag set to 1
151 residuals, jacobian, 1);
152 }
153
154
155 /// Overload the output function
156 void output(std::ostream& outfile)
157 {
158 FiniteElement::output(outfile);
159 }
160
161
162 /// Output function: x,y,[z],u,v,[w],p in tecplot format
163 void output(std::ostream& outfile, const unsigned& nplot)
164 {
165 FiniteElement::output(outfile, nplot);
166 }
167
168 protected:
169 /// This function returns the residuals for the traction function.
170 /// flag=1(or 0): do (or don't) compute the Jacobian as well.
172 Vector<double>& residuals,
173 DenseMatrix<double>& jacobian,
174 const unsigned& flag);
175
176
177 /// The "global" intrinsic coordinate of the element when
178 /// viewed as part of a geometric object should be given by
179 /// the FaceSpaceTimeElement representation, by default
180 double zeta_nodal(const unsigned& n,
181 const unsigned& k,
182 const unsigned& i) const
183 {
184 return FaceElement::zeta_nodal(n, k, i);
185 } // End of zeta_nodal
186
187
188 /// Access function that returns the local equation numbers
189 /// for velocity components.
190 /// u_local_eqn(n,i)=local equation number or < 0 if pinned.
191 /// The default is to asssume that n is the local node number
192 /// and the i-th velocity component is the i-th unknown stored at the node.
193 virtual inline int u_local_eqn(const unsigned& n, const unsigned& i)
194 {
195 return nodal_local_eqn(n, i);
196 }
197
198 /// Function to compute the shape and test functions and to return
199 /// the Jacobian of mapping
200 inline double shape_and_test_at_knot(const unsigned& ipt,
201 Shape& psi,
202 Shape& test) const
203 {
204 // Number of nodes in each direction
205 const unsigned nnode_1d = 3;
206
207 // Find the dimension of the element
208 const unsigned el_dim = dim();
209
210 // Sanity check
211 if (el_dim != 2)
212 {
213 // Throw an error
214 throw OomphLibError(
215 "Can only deal with 2D space-time traction elements!",
216 OOMPH_CURRENT_FUNCTION,
217 OOMPH_EXCEPTION_LOCATION);
218 }
219
220 // Storage for the local coordinates of the integration point
221 Vector<double> s(el_dim);
222
223 // Set the local coordinate
224 for (unsigned i = 0; i < el_dim; i++)
225 {
226 // Get the i-th local coordinate associated with the ipt-th knot point
227 s[i] = integral_pt()->knot(ipt, i);
228 }
229
230 //--------------------------
231 // Call the shape functions:
232 //--------------------------
233 //--------start_of_shape------------------------------------------------
234 // Local storage
235 double psi_values[2][nnode_1d];
236
237 // Call the OneDimensional shape function (for x or y) with quadratic
238 // interpolation
239 OneDimLagrange::shape<3>(s[0], psi_values[0]);
240
241 // Call the OneDimensional shape function. NOTE: The second local
242 // coordinate will always be time for a space-time element!
244
245 // Index for the shape functions
246 unsigned index = 0;
247
248 // Loop over the first local coordinate direction
249 for (unsigned j = 0; j < nnode_1d; j++)
250 {
251 // Loop over the first local coordinate direction
252 for (unsigned i = 0; i < nnode_1d; i++)
253 {
254 // Multiply the two 1D functions together to get the 2D function
255 psi[index] = psi_values[0][i] * psi_values[1][j];
256
257 // Increment the index
258 index++;
259 }
260 } // for (unsigned j=0;j<nnode_1d;j++)
261 //--------end_of_shape--------------------------------------------------
262
263 //-------------------------
264 // Call the test functions:
265 //-------------------------
266 //--------start_of_shape------------------------------------------------
267 // Local storage
268 double test_values[2][nnode_1d];
269
270 // Call the OneDimensional shape function (for x or y) with quadratic
271 // interpolation
272 OneDimLagrange::shape<3>(s[0], test_values[0]);
273
274 // Call the OneDimensional shape function. NOTE: The second local
275 // coordinate will always be time for a space-time element!
277
278 // Reset the value of index to loop over the test functions
279 index = 0;
280
281 // Loop over the first local coordinate direction
282 for (unsigned j = 0; j < nnode_1d; j++)
283 {
284 // Loop over the first local coordinate direction
285 for (unsigned i = 0; i < nnode_1d; i++)
286 {
287 // Multiply the two 1D functions together to get the 2D function
288 test[index] = test_values[0][i] * test_values[1][j];
289
290 // Increment the index
291 index++;
292 }
293 } // for (unsigned j=0;j<nnode_1d;j++)
294 //--------end_of_shape--------------------------------------------------
295
296 // Return the value of the jacobian
297 return J_eulerian_at_knot(ipt);
298 } // End of shape_and_test_at_knot
299
300
301 /// Function to calculate the traction applied to the fluid
302 void get_traction(const double& time,
303 const Vector<double>& x,
304 const Vector<double>& n,
305 Vector<double>& result)
306 {
307 // If the function pointer is zero return zero
308 if (Traction_fct_pt == 0)
309 {
310 // Loop over the (spatial) dimensions
311 for (unsigned i = 0; i < Dim - 1; i++)
312 {
313 // Set the i-th traction component to zero
314 result[i] = 0.0;
315 }
316 }
317 // If the traction function pointer has been set
318 else
319 {
320 // Call the function pointer
321 (*Traction_fct_pt)(time, x, n, result);
322 } // if (Traction_fct_pt==0)
323 } // End of get_traction
324
325
326 /// The highest dimension of the problem
327 unsigned Dim;
328 };
329
330 /// ////////////////////////////////////////////////////////////////////
331 /// ////////////////////////////////////////////////////////////////////
332 /// ////////////////////////////////////////////////////////////////////
333
334 //=========================================================================
335 /// Function that returns the residuals for the imposed traction
336 /// Navier-Stokes equations
337 //=========================================================================
338 template<class ELEMENT>
341 Vector<double>& residuals,
342 DenseMatrix<double>& jacobian,
343 const unsigned& flag)
344 {
345 // Find out how many nodes there are
346 unsigned n_node = nnode();
347
348 // Set up memory for the shape functions
349 Shape psif(n_node);
350
351 // Set up memory for the test functions
352 Shape testf(n_node);
353
354 // Number of integration points
355 unsigned n_intpt = integral_pt()->nweight();
356
357 // Integers to store local equation numbers
358 int local_eqn = 0;
359
360 // Loop over the integration points
361 for (unsigned ipt = 0; ipt < n_intpt; ipt++)
362 {
363 // Get the integral weight
364 double w = integral_pt()->weight(ipt);
365
366 // Find the shape and test functions and return the Jacobian
367 // of the mapping
368 double J = shape_and_test_at_knot(ipt, psif, testf);
369
370 // Premultiply the weights and the Jacobian
371 double W = w * J;
372
373 // Get continuous time
374 double interpolated_t = 0.0;
375
376 // Need to find position to feed into traction function
377 Vector<double> interpolated_x(Dim - 1, 0.0);
378
379 // Allocate space for the outer unit normal
380 Vector<double> interpolated_spacetime_n(Dim, 0.0);
381
382 // Allocate space for the outer unit normal (in the given time slice)
383 Vector<double> interpolated_n(Dim - 1, 0.0);
384
385 // Allocate space for the traction vector
386 Vector<double> traction(Dim - 1, 0.0);
387
388 // Calculate velocities and derivatives
389 for (unsigned l = 0; l < n_node; l++)
390 {
391 // Update the interpolated time value
392 interpolated_t += raw_nodal_position(l, Dim - 1) * psif(l);
393
394 // Loop over velocity components
395 for (unsigned i = 0; i < Dim - 1; i++)
396 {
397 // Update the i-th interpolated coordinate value
398 interpolated_x[i] += nodal_position(l, i) * psif(l);
399 }
400 } // for (unsigned l=0;l<n_node;l++)
401
402 // Calculate the outer unit normal
403 outer_unit_normal(ipt, interpolated_spacetime_n);
404
405#ifdef PARANOID
406 // Sanity check: Make sure the normal component in the time-direction is
407 // zero. If it isn't then something has gone wrong...
408 if (std::fabs(interpolated_spacetime_n[Dim - 1]) > 1.0e-15)
409 {
410 // Create an output stream
411 std::ostringstream error_message_stream;
412
413 // Create the error message
414 error_message_stream
415 << "The component of the outer unit normal in the "
416 << "time-direction\nshould be zero but the actual "
417 << "value is: " << interpolated_spacetime_n[Dim - 1] << std::endl;
418
419 // Throw the error message
420 throw OomphLibError(error_message_stream.str(),
421 OOMPH_CURRENT_FUNCTION,
422 OOMPH_EXCEPTION_LOCATION);
423 }
424#endif
425
426 // Loop over spatial coordinates
427 for (unsigned i = 0; i < Dim - 1; i++)
428 {
429 // Copy the i-th normal component over
430 interpolated_n[i] = interpolated_spacetime_n[i];
431 }
432
433 // Get the user-defined traction
434 get_traction(interpolated_t, interpolated_x, interpolated_n, traction);
435
436 // Now add to the appropriate equations:
437
438 // Loop over the test functions
439 for (unsigned l = 0; l < n_node; l++)
440 {
441 // Loop over the velocity components
442 for (unsigned i = 0; i < Dim - 1; i++)
443 {
444 // Get the local equation number associated with this unknown and node
445 local_eqn = u_local_eqn(l, i);
446
447 // If it's not a boundary condition
448 if (local_eqn >= 0)
449 {
450 // Add the user-defined traction terms
451 residuals[local_eqn] += traction[i] * testf[l] * W;
452
453 // Assuming the traction DOES NOT depend upon velocities
454 // or pressures, the jacobian is always zero, so no jacobian
455 // terms are required
456 }
457 } // for (unsigned i=0;i<Dim;i++)
458 } // for (unsigned l=0;l<n_node;l++)
459 } // for (unsigned ipt=0;ipt<n_intpt;ipt++)
460 } // End of fill_in_generic_residual_contribution_fluid_traction
461
462 /// //////////////////////////////////////////////////////////////////////
463 /// //////////////////////////////////////////////////////////////////////
464 /// //////////////////////////////////////////////////////////////////////
465
466 //======================================================================
467 /// A class for elements that allow the imposition of an applied traction
468 /// to the Navier--Stokes equations
469 /// The geometrical information can be read from the FaceGeometry<ELEMENT>
470 /// class and and thus, we can be generic enough without the need to have
471 /// a separate equations class.
472 ///
473 /// THIS IS THE REFINEABLE VERSION.
474 //======================================================================
475 template<class ELEMENT>
477 : public virtual NavierStokesMixedOrderSpaceTimeTractionElement<ELEMENT>,
479 {
480 public:
481 /// Constructor, which takes a "bulk" element and the face index
483 FiniteElement* const& element_pt, const int& face_index)
484 : // we're calling this from the constructor of the refineable version.
486 element_pt, face_index, true)
487 {
488 }
489
490
491 /// Destructor should not delete anything
493
494
495 /// Number of continuously interpolated values are the
496 /// same as those in the bulk element.
498 {
499 return dynamic_cast<ELEMENT*>(this->bulk_element_pt())
501 }
502
503 /// This function returns just the residuals
505 {
506 // Call the generic residuals function using a dummy matrix argument
509 }
510
511 /// This function returns the residuals and the Jacobian
513 DenseMatrix<double>& jacobian)
514 {
515 // Call the generic routine
517 residuals, jacobian, 1);
518 }
519
520
521 protected:
522 /// This function returns the residuals for the traction function.
523 /// flag=1(or 0): do (or don't) compute the Jacobian as well.
525 Vector<double>& residuals,
526 DenseMatrix<double>& jacobian,
527 const unsigned& flag);
528 };
529
530 /// ////////////////////////////////////////////////////////////////////
531 /// ////////////////////////////////////////////////////////////////////
532 /// ////////////////////////////////////////////////////////////////////
533
534 //============================================================================
535 /// Function that returns the residuals for the imposed traction Navier_Stokes
536 /// equations
537 //============================================================================
538 template<class ELEMENT>
541 Vector<double>& residuals,
542 DenseMatrix<double>& jacobian,
543 const unsigned& flag)
544 {
545 // Get the indices at which the velocity components are stored
546 unsigned u_nodal_index[this->Dim - 1];
547
548 // Loop over the velocity components
549 for (unsigned i = 0; i < this->Dim - 1; i++)
550 {
551 // Get the i-th velocity component index
552 u_nodal_index[i] =
553 dynamic_cast<ELEMENT*>(this->bulk_element_pt())->u_index_nst(i);
554 }
555
556 // Find out how many nodes there are
557 unsigned n_node = nnode();
558
559 // Set up memory for the shape functions
560 Shape psif(n_node);
561
562 // Set up memory for the test functions
563 Shape testf(n_node);
564
565 // Number of integration points
566 unsigned n_intpt = integral_pt()->nweight();
567
568 // Integers to store local equation numbers
569 int local_eqn = 0;
570
571 // Loop over the integration points
572 for (unsigned ipt = 0; ipt < n_intpt; ipt++)
573 {
574 // Get the integral weight
575 double w = integral_pt()->weight(ipt);
576
577 // Find the shape and test functions and return the Jacobian
578 // of the mapping
579 double J = this->shape_and_test_at_knot(ipt, psif, testf);
580
581 // Premultiply the weights and the Jacobian
582 double W = w * J;
583
584 // Get continuous time
585 double interpolated_t = 0.0;
586
587 // Need to find position to feed into traction function
588 Vector<double> interpolated_x(this->Dim - 1, 0.0);
589
590 // Allocate space for the outer unit normal
591 Vector<double> interpolated_spacetime_n(this->Dim, 0.0);
592
593 // Allocate space for the outer unit normal (in the given time slice)
594 Vector<double> interpolated_n(this->Dim - 1, 0.0);
595
596 // Allocate space for the traction vector
597 Vector<double> traction(this->Dim - 1, 0.0);
598
599 // Calculate velocities and derivatives
600 for (unsigned l = 0; l < n_node; l++)
601 {
602 // Update the interpolated time value
603 interpolated_t += raw_nodal_position(l, this->Dim - 1) * psif(l);
604
605 // Loop over velocity components
606 for (unsigned i = 0; i < this->Dim - 1; i++)
607 {
608 // Update the i-th interpolated coordinate value
609 interpolated_x[i] += this->nodal_position(l, i) * psif(l);
610 }
611 } // for (unsigned l=0;l<n_node;l++)
612
613 // Calculate the outer unit normal
614 this->outer_unit_normal(ipt, interpolated_spacetime_n);
615
616#ifdef PARANOID
617 // Sanity check: Make sure the normal component in the time-direction is
618 // zero. If it isn't then something has gone wrong...
619 if (std::fabs(interpolated_n[this->Dim - 1]) > 1.0e-15)
620 {
621 // Create an output stream
622 std::ostringstream error_message_stream;
623
624 // Create the error message
625 error_message_stream << "The component of the outer unit normal in the "
626 << "time-direction\nshould be zero but the actual "
627 << "value is: " << interpolated_n[this->Dim - 1]
628 << std::endl;
629
630 // Throw the error message
631 throw OomphLibError(error_message_stream.str(),
632 OOMPH_CURRENT_FUNCTION,
633 OOMPH_EXCEPTION_LOCATION);
634 }
635#endif
636
637 // Loop over spatial coordinates
638 for (unsigned i = 0; i < this->Dim - 1; i++)
639 {
640 // Copy the i-th normal component over
641 interpolated_n[i] = interpolated_spacetime_n[i];
642 }
643
644 // Get the user-defined traction
645 this->get_traction(
646 interpolated_t, interpolated_x, interpolated_n, traction);
647
648 // Now add to the appropriate equations:
649
650 // Number of master nodes and storage for the weight of the shape function
651 unsigned n_master = 1;
652
653 // The weight of the shape function
654 double hang_weight = 1.0;
655
656 // Pointer to hang info object
657 HangInfo* hang_info_pt = 0;
658
659 // Loop over the nodes for the test functions/equations
660 //-----------------------------------------------------
661 for (unsigned l = 0; l < n_node; l++)
662 {
663 // Local boolean to indicate whether the node is hanging
664 bool is_node_hanging = this->node_pt(l)->is_hanging();
665
666 // If the node is hanging
667 if (is_node_hanging)
668 {
669 // Get the HangInfo pointer associated with this hanging node
670 hang_info_pt = this->node_pt(l)->hanging_pt();
671
672 // Read out number of master nodes from hanging data
673 n_master = hang_info_pt->nmaster();
674 }
675 // Otherwise the node is its own master
676 else
677 {
678 // If it's not hanging then it is its (one and only) master node
679 n_master = 1;
680 }
681
682 // Loop over the master nodes
683 for (unsigned m = 0; m < n_master; m++)
684 {
685 // Loop over velocity components for equations
686 for (unsigned i = 0; i < this->Dim - 1; i++)
687 {
688 // Get the equation number:
689
690 // If the node is hanging
691 if (is_node_hanging)
692 {
693 // Get the equation number from the master node
694 local_eqn = this->local_hang_eqn(hang_info_pt->master_node_pt(m),
695 u_nodal_index[i]);
696
697 // Get the hang weight from the master node
698 hang_weight = hang_info_pt->master_weight(m);
699 }
700 // If the node is not hanging
701 else
702 {
703 // Local equation number
704 local_eqn = this->nodal_local_eqn(l, u_nodal_index[i]);
705
706 // Node contributes with full weight
707 hang_weight = 1.0;
708 }
709
710 // If it's not a boundary condition...
711 if (local_eqn >= 0)
712 {
713 // Add the user-defined traction terms
714 residuals[local_eqn] += traction[i] * testf[l] * W * hang_weight;
715
716 // Assuming the the traction DOES NOT depend upon velocities
717 // or pressures, the jacobian is always zero, so no jacobian
718 // terms are required.
719 }
720 } // for (unsigned i=0;i<this->Dim-1;i++)
721 } // for (unsigned m=0;m<n_master;m++)
722 } // for (unsigned l=0;l<n_node;l++)
723 } // for (unsigned ipt=0;ipt<n_intpt;ipt++)
724 } // End of refineable_fill_in_generic_residual_contribution_fluid_traction
725} // End of namespace oomph
726
727#endif
static char t char * s
Definition: cfortran.h:568
cstr elem_len * i
Definition: cfortran.h:603
char t
Definition: cfortran.h:568
FaceElements are elements that coincide with the faces of higher-dimensional "bulk" elements....
Definition: elements.h:4338
int & face_index()
Index of the face (a number that uniquely identifies the face in the element)
Definition: elements.h:4626
double zeta_nodal(const unsigned &n, const unsigned &k, const unsigned &i) const
In a FaceElement, the "global" intrinsic coordinate of the element along the boundary,...
Definition: elements.h:4497
FiniteElement *& bulk_element_pt()
Pointer to higher-dimensional "bulk" element.
Definition: elements.h:4735
double J_eulerian_at_knot(const unsigned &ipt) const
Return the Jacobian of the mapping from local to global coordinates at the ipt-th integration point O...
Definition: elements.cc:5328
//////////////////////////////////////////////////////////////////// ////////////////////////////////...
Definition: elements.h:4998
A general Finite Element class.
Definition: elements.h:1313
Integral *const & integral_pt() const
Return the pointer to the integration scheme (const version)
Definition: elements.h:1963
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
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
Node *& node_pt(const unsigned &n)
Return a pointer to the local node n.
Definition: elements.h:2175
virtual void build_face_element(const int &face_index, FaceElement *face_element_pt)
Function for building a lower dimensional FaceElement on the specified face of the FiniteElement....
Definition: elements.cc:5132
virtual unsigned nnode_1d() const
Return the number of nodes along one edge of the element Default is to return zero — must be overload...
Definition: elements.h:2218
bool has_hanging_nodes() const
Return boolean to indicate if any of the element's nodes are geometrically hanging.
Definition: elements.h:2470
static DenseMatrix< double > Dummy_matrix
Empty dense matrix used as a dummy argument to combined residual and jacobian functions in the case w...
Definition: elements.h:227
Class that contains data for hanging nodes.
Definition: nodes.h:742
Node *const & master_node_pt(const unsigned &i) const
Return a pointer to the i-th master node.
Definition: nodes.h:791
unsigned nmaster() const
Return the number of master nodes.
Definition: nodes.h:785
double const & master_weight(const unsigned &i) const
Return weight for dofs on i-th master node.
Definition: nodes.h:808
virtual double knot(const unsigned &i, const unsigned &j) const =0
Return local coordinate s[j] of i-th integration point.
A class for elements that allow the imposition of an applied traction to the Navier–Stokes equations ...
virtual int u_local_eqn(const unsigned &n, const unsigned &i)
Access function that returns the local equation numbers for velocity components. u_local_eqn(n,...
NavierStokesMixedOrderSpaceTimeTractionElement(FiniteElement *const &element_pt, const int &face_index, const bool &called_from_refineable_constructor=false)
Constructor, which takes a "bulk" element and the value of the index and its limit.
double shape_and_test_at_knot(const unsigned &ipt, Shape &psi, Shape &test) const
Function to compute the shape and test functions and to return the Jacobian of mapping.
void fill_in_generic_residual_contribution_fluid_traction(Vector< double > &residuals, DenseMatrix< double > &jacobian, const unsigned &flag)
This function returns the residuals for the traction function. flag=1(or 0): do (or don't) compute th...
double zeta_nodal(const unsigned &n, const unsigned &k, const unsigned &i) const
The "global" intrinsic coordinate of the element when viewed as part of a geometric object should be ...
void fill_in_contribution_to_residuals(Vector< double > &residuals)
This function returns just the residuals.
void fill_in_contribution_to_jacobian(Vector< double > &residuals, DenseMatrix< double > &jacobian)
This function returns the residuals and the jacobian.
void get_traction(const double &time, const Vector< double > &x, const Vector< double > &n, Vector< double > &result)
Function to calculate the traction applied to the fluid.
void output(std::ostream &outfile, const unsigned &nplot)
Output function: x,y,[z],u,v,[w],p in tecplot format.
void(*&)(const double &t, const Vector< double > &x, const Vector< double > &n, Vector< double > &result) traction_fct_pt()
Access function for the imposed traction pointer.
void(* Traction_fct_pt)(const double &time, const Vector< double > &x, const Vector< double > &n, Vector< double > &result)
Pointer to an imposed traction function.
unsigned ndim() const
Return (Eulerian) spatial dimension of the node.
Definition: nodes.h:1054
/////////////////////////////////////////////////////////////////// /////////////////////////////////...
An OomphLibError object which should be thrown when an run-time error is encountered....
RefineableElements are FiniteElements that may be subdivided into children to provide a better local ...
////////////////////////////////////////////////////////////////////// //////////////////////////////...
void fill_in_contribution_to_residuals(Vector< double > &residuals)
This function returns just the residuals.
RefineableNavierStokesMixedOrderSpaceTimeTractionElement(FiniteElement *const &element_pt, const int &face_index)
Constructor, which takes a "bulk" element and the face index.
void fill_in_contribution_to_jacobian(Vector< double > &residuals, DenseMatrix< double > &jacobian)
This function returns the residuals and the Jacobian.
void refineable_fill_in_generic_residual_contribution_fluid_traction(Vector< double > &residuals, DenseMatrix< double > &jacobian, const unsigned &flag)
This function returns the residuals for the traction function. flag=1(or 0): do (or don't) compute th...
unsigned ncont_interpolated_values() const
Number of continuously interpolated values are the same as those in the bulk element.
A Class for shape functions. In simple cases, the shape functions have only one index that can be tho...
Definition: shape.h:76
unsigned Dim
Dimension of zeta tuples (set by get_dim_helper) – needed because we store the scalar coordinates in ...
Definition: multi_domain.cc:60
void shape< 3 >(const double &s, double *Psi)
1D shape functions specialised to quadratic order (3 Nodes)
Definition: shape.h:967
void shape< 3 >(const double &s, double *Psi)
1D shape functions specialised to quadratic order (3 Nodes)
Definition: shape.h:1124
void shape< 3 >(const double &s, double *Psi)
1D shape functions specialised to quadratic order (3 Nodes)
Definition: shape.h:635
//////////////////////////////////////////////////////////////////// ////////////////////////////////...