27 #ifndef OOMPH_MIXED_ORDER_PETROV_GALERKIN_SPACE_TIME_NAVIER_STOKES_ELEMENTS_HEADER
28 #define OOMPH_MIXED_ORDER_PETROV_GALERKIN_SPACE_TIME_NAVIER_STOKES_ELEMENTS_HEADER
32 #include <oomph-lib-config.h>
65 const unsigned& compute_jacobian_flag) = 0;
80 template<
class ELEMENT>
93 const bool& called_from_refineable_constructor =
false)
102 if (!called_from_refineable_constructor)
105 if (element_pt->
dim() == 3)
118 std::ostringstream error_message_stream;
122 <<
"This flux element will not work correctly "
123 <<
"if nodes are hanging!" << std::endl;
127 OOMPH_CURRENT_FUNCTION,
128 OOMPH_EXCEPTION_LOCATION);
147 const unsigned& flag);
154 std::ostringstream error_message;
157 error_message <<
"fill_in_contribution_to_residuals() must not be "
158 <<
"called directly.\nsince it uses the local equation "
159 <<
"numbering of the bulk element\nwhich calls the "
160 <<
"relevant helper function directly." << std::endl;
164 error_message.str(), OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION);
173 std::ostringstream error_message;
176 error_message <<
"fill_in_contribution_to_jacobian() must not be "
177 <<
"called directly.\nsince it uses the local equation "
178 <<
"numbering of the bulk element\nwhich calls the "
179 <<
"relevant helper function directly." << std::endl;
183 error_message.str(), OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION);
196 void output(std::ostream& outfile,
const unsigned& nplot)
214 template<
class ELEMENT>
219 const unsigned& flag)
223 OOMPH_CURRENT_FUNCTION,
224 OOMPH_EXCEPTION_LOCATION);
228 unsigned my_dim = this->dim();
245 unsigned n_intpt = this->integral_pt()->nweight();
251 int local_unknown = 0;
254 ELEMENT* bulk_el_pt =
dynamic_cast<ELEMENT*
>(this->bulk_element_pt());
257 unsigned n_pres = bulk_el_pt->npres_nst();
260 double re = bulk_el_pt->re();
269 for (
unsigned ipt = 0; ipt < n_intpt; ipt++)
272 double w = this->integral_pt()->weight(ipt);
275 for (
unsigned i = 0;
i < my_dim;
i++)
278 s[
i] = this->integral_pt()->knot(ipt,
i);
282 s_bulk = this->local_coordinate_in_bulk(
s);
288 this->outer_unit_normal(ipt, unit_normal);
291 bulk_el_pt->interpolated_u_nst(s_bulk, velocity);
297 for (
unsigned i = 0;
i < my_dim;
i++)
300 flux += velocity[
i] * unit_normal[
i];
311 double interpolated_press = bulk_el_pt->interpolated_p_nst(s_bulk);
314 bulk_el_pt->pshape_nst(s_bulk, psip, testp);
317 double J = this->J_eulerian(
s);
324 for (
unsigned l = 0; l < n_pres; l++)
327 local_eqn = bulk_el_pt->p_local_eqn(l);
333 residuals[local_eqn] -= re * flux * interpolated_press * testp[l] *
W;
339 for (
unsigned l2 = 0; l2 < n_pres; l2++)
342 local_unknown = bulk_el_pt->p_local_eqn(l2);
345 if (local_unknown >= 0)
348 jacobian(local_eqn, local_unknown) -=
349 re * flux * psip[l2] * testp[l] *
W;
410 std::map<
Data*, std::vector<int>>& eqn_number_backup) = 0;
417 const unsigned& face_index) = 0;
433 const unsigned& which_one = 0) = 0;
466 template<
unsigned DIM>
566 DShape& dtestdx)
const = 0;
577 DShape& dtestdx)
const = 0;
601 DShape& dpsidx)
const = 0;
607 DShape& dtestdx)
const = 0;
615 DShape& dppsidx)
const = 0;
623 DShape& dptestdx)
const = 0;
633 DShape& dptestdx)
const = 0;
656 (*Body_force_fct_pt)(time, x, result);
689 for (
unsigned i = 0;
i < DIM;
i++)
698 for (
unsigned j = 0; j < DIM; j++)
701 d_body_force_dx(j,
i) = (body_force_pls[j] - body_force[j]) / eps_fd;
747 double source_pls = 0.0;
753 for (
unsigned i = 0;
i < DIM;
i++)
762 for (
unsigned j = 0; j < DIM; j++)
765 gradient[
i] = (source_pls - source) / eps_fd;
781 const unsigned& flag);
790 const unsigned& flag);
798 double*
const& parameter_pt,
802 const unsigned& flag);
858 Data* reynolds_strouhal_data_pt)
862 if (reynolds_strouhal_data_pt == 0)
865 std::ostringstream error_message_stream;
869 <<
"User supplied the Reynolds x Strouhal number as "
870 <<
"external data\nbut the pointer provided is a "
871 <<
"null pointer!" << std::endl;
875 OOMPH_CURRENT_FUNCTION,
876 OOMPH_EXCEPTION_LOCATION);
900 const double&
re()
const
922 unsigned data_index = 0;
925 unsigned re_st_index = 0;
946 unsigned data_index = 0;
949 unsigned re_st_index = 0;
1078 Shape& test)
const = 0;
1085 double u_nst(
const unsigned& n,
const unsigned&
i)
const
1092 double u_nst(
const unsigned&
t,
const unsigned& n,
const unsigned&
i)
const
1099 throw OomphLibError(
"Space-time elements cannot have history values!",
1100 OOMPH_CURRENT_FUNCTION,
1101 OOMPH_EXCEPTION_LOCATION);
1121 std::ostringstream error_message_stream;
1124 error_message_stream <<
"Input index " <<
i <<
" does not correspond "
1125 <<
"to a velocity component when DIM=" << DIM
1126 <<
"!" << std::endl;
1130 OOMPH_CURRENT_FUNCTION,
1131 OOMPH_EXCEPTION_LOCATION);
1178 const unsigned&
i)
const
1181 unsigned n_node =
nnode();
1187 DShape dpsidx(n_node, DIM + 1);
1191 const unsigned el_dim = this->
dim();
1208 double interpolated_dudt = 0.0;
1214 for (
unsigned l = 0; l < n_node; l++)
1217 interpolated_dudt +=
nodal_value(l, u_nodal_index) * dpsidx(l, DIM);
1221 return interpolated_dudt;
1245 virtual double p_nst(
const unsigned& n_p)
const = 0;
1248 virtual double p_nst(
const unsigned&
t,
const unsigned& n_p)
const = 0;
1251 virtual void fix_pressure(
const unsigned& p_dof,
const double& p_value) = 0;
1330 const unsigned& which_one = 0);
1344 const unsigned& nplot)
const
1353 for (
unsigned iplot = 0; iplot < num_plot_points; iplot++)
1375 std::stringstream error_stream;
1378 error_stream <<
"These Navier Stokes elements only store " << DIM + 1
1380 <<
"but i is currently " <<
i << std::endl;
1384 OOMPH_CURRENT_FUNCTION,
1385 OOMPH_EXCEPTION_LOCATION);
1395 std::ofstream& file_out,
1397 const unsigned& nplot,
1405 std::stringstream error_stream;
1408 error_stream <<
"These Navier Stokes elements only store " << DIM + 1
1409 <<
" fields, but i is currently " <<
i << std::endl;
1413 error_stream.str(), OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION);
1421 double interpolated_t = 0.0;
1430 for (
unsigned iplot = 0; iplot < num_plot_points; iplot++)
1436 for (
unsigned j = 0; j < DIM; j++)
1452 (*exact_soln_pt)(interpolated_t, spatial_coordinates, exact_soln);
1455 file_out << exact_soln[
i] << std::endl;
1481 std::stringstream error_stream;
1484 error_stream <<
"These Navier Stokes elements only store " << DIM + 1
1485 <<
" fields,\nbut i is currently " <<
i << std::endl;
1489 error_stream.str(), OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION);
1502 unsigned n_plot = 5;
1511 void output(std::ostream& outfile,
const unsigned& n_plot);
1519 unsigned n_plot = 5;
1528 void output(FILE* file_pt,
const unsigned& n_plot);
1537 unsigned n_plot = 5;
1547 void full_output(std::ostream& outfile,
const unsigned& n_plot);
1554 const unsigned& nplot,
1567 const unsigned& nplot,
1575 const unsigned& nplot,
1637 unsigned compute_jacobian_flag = 0;
1644 compute_jacobian_flag);
1654 unsigned compute_jacobian_flag = 1;
1661 compute_jacobian_flag);
1673 unsigned compute_matrices_flag = 2;
1677 residuals, jacobian, mass_matrix, compute_matrices_flag);
1687 unsigned compute_jacobian_flag = 0;
1695 compute_jacobian_flag);
1702 double*
const& parameter_pt,
1707 unsigned compute_jacobian_flag = 1;
1715 compute_jacobian_flag);
1722 double*
const& parameter_pt,
1728 unsigned compute_matrices_flag = 2;
1734 dmass_matrix_dparam,
1735 compute_matrices_flag);
1745 unsigned compute_jacobian_flag = 0;
1759 unsigned compute_jacobian_flag = 1;
1763 residuals, jacobian, compute_jacobian_flag);
1769 std::map<
Data*, std::vector<int>>& eqn_number_backup)
1774 for (
unsigned j = 0; j < nint; j++)
1777 if (eqn_number_backup[data_pt].
size() == 0)
1779 unsigned nvalue = data_pt->
nvalue();
1780 eqn_number_backup[data_pt].resize(nvalue);
1781 for (
unsigned i = 0;
i < nvalue;
i++)
1784 eqn_number_backup[data_pt][
i] = data_pt->
eqn_number(
i);
1793 unsigned nnod = this->
nnode();
1794 for (
unsigned j = 0; j < nnod; j++)
1797 if (eqn_number_backup[nod_pt].
size() == 0)
1799 unsigned nvalue = nod_pt->
nvalue();
1800 eqn_number_backup[nod_pt].resize(nvalue);
1801 for (
unsigned i = 0;
i < nvalue;
i++)
1828 if (solid_nod_pt != 0)
1831 if (eqn_number_backup[solid_posn_data_pt].
size() == 0)
1833 unsigned nvalue = solid_posn_data_pt->
nvalue();
1834 eqn_number_backup[solid_posn_data_pt].resize(nvalue);
1835 for (
unsigned i = 0;
i < nvalue;
i++)
1838 eqn_number_backup[solid_posn_data_pt][
i] =
1842 solid_posn_data_pt->
pin(
i);
1855 const unsigned& face_index) = 0;
1862 std::ostream& outfile)
1865 for (
unsigned e = 0;
e < nel;
e++)
1869 outfile <<
"ZONE" << std::endl;
1874 for (
unsigned ipt = 0; ipt < n; ipt++)
1876 for (
unsigned i = 0;
i < DIM;
i++)
1882 for (
unsigned i = 0;
i < DIM + 1;
i++)
1884 outfile << x[
i] <<
" ";
1886 for (
unsigned i = 0;
i < DIM;
i++)
1888 outfile << unit_normal[
i] <<
" ";
1890 outfile << std::endl;
1902 for (
unsigned e = 0;
e < nel;
e++)
1923 unsigned n_node =
nnode();
1932 for (
unsigned i = 0;
i < DIM;
i++)
1941 for (
unsigned l = 0; l < n_node; l++)
1954 unsigned n_node =
nnode();
1966 double interpolated_u = 0.0;
1969 for (
unsigned l = 0; l < n_node; l++)
1972 interpolated_u +=
nodal_value(l, u_nodal_index) * psi[l];
1976 return interpolated_u;
1984 const unsigned&
i)
const
1987 std::ostringstream error_message_stream;
1990 error_message_stream <<
"This interface doesn't make sense in "
1991 <<
"space-time elements!" << std::endl;
1995 OOMPH_CURRENT_FUNCTION,
1996 OOMPH_EXCEPTION_LOCATION);
2011 unsigned n_node =
nnode();
2023 unsigned n_u_dof = 0;
2026 for (
unsigned l = 0; l < n_node; l++)
2032 if (global_eqn >= 0)
2040 du_ddata.resize(n_u_dof, 0.0);
2041 global_eqn_number.resize(n_u_dof, 0);
2046 for (
unsigned l = 0; l < n_node; l++)
2050 if (global_eqn >= 0)
2053 global_eqn_number[count] = global_eqn;
2055 du_ddata[count] = psi[l];
2076 double interpolated_p = 0.0;
2079 for (
unsigned l = 0; l < n_pres; l++)
2081 interpolated_p +=
p_nst(l) * psi[l];
2084 return (interpolated_p);
2101 double interpolated_p = 0.0;
2104 for (
unsigned l = 0; l < n_pres; l++)
2106 interpolated_p +=
p_nst(
t, l) * psi[l];
2109 return (interpolated_p);
2118 data.resize(2 * DIM + 2);
2121 for (
unsigned i = 0;
i < DIM + 1;
i++)
2128 for (
unsigned i = 0;
i < DIM;
i++)
2149 template<
unsigned DIM>
2151 :
public virtual QElement<DIM + 1, 3>,
2188 const unsigned& ipt,
2257 return this->
node_pt(this->Pconv[n_p]);
2269 const unsigned* it = std::find(this->Pconv, this->Pconv + n_p, n);
2272 if (it != this->Pconv + n_p)
2275 return (it - (this->Pconv));
2292 return std::find(this->Pconv, this->Pconv + n_p, n) != this->Pconv + n_p;
2315 return static_cast<int>(DIM);
2338 double p_nst(
const unsigned&
t,
const unsigned& n_p)
const
2350 return static_cast<unsigned>(pow(2.0,
static_cast<int>(DIM + 1)));
2386 std::set<std::pair<Data*, unsigned>>& paired_load_data);
2398 std::set<std::pair<Data*, unsigned>>& paired_pressure_data);
2410 void output(std::ostream& outfile,
const unsigned& nplot)
2426 void output(FILE* file_pt,
const unsigned& nplot)
2436 unsigned ndof_types() const
2438 // Return the number of dof types being used in the mesh
2440 } // End of ndof_types
2453 void get_dof_numbers_for_unknowns(
2454 std::list<std::pair<unsigned long,unsigned> >& dof_lookup_list) const;
2471 unsigned n_node_1d = 3;
2474 double psi_values[n_dim][n_node_1d];
2487 for (
unsigned i = 0;
i < n_node_1d;
i++)
2490 for (
unsigned j = 0; j < n_node_1d; j++)
2493 for (
unsigned k = 0; k < n_node_1d; k++)
2496 psi[index] = psi_values[0][k] * psi_values[1][j] * psi_values[2][
i];
2517 unsigned n_node_1d = 3;
2520 double psi_values[n_dim][n_node_1d];
2521 double dpsi_values[n_dim][n_node_1d];
2537 for (
unsigned i = 0;
i < n_node_1d;
i++)
2540 for (
unsigned j = 0; j < n_node_1d; j++)
2543 for (
unsigned k = 0; k < n_node_1d; k++)
2547 dpsi_values[0][k] * psi_values[1][j] * psi_values[2][
i];
2551 psi_values[0][k] * dpsi_values[1][j] * psi_values[2][
i];
2555 psi_values[0][k] * psi_values[1][j] * dpsi_values[2][
i];
2558 psi[index] = psi_values[0][k] * psi_values[1][j] * psi_values[2][
i];
2579 unsigned n_node_1d = 3;
2582 double test_values[n_dim][n_node_1d];
2583 double dtest_values[n_dim][n_node_1d];
2601 for (
unsigned i = 0;
i < n_node_1d;
i++)
2604 for (
unsigned j = 0; j < n_node_1d; j++)
2607 for (
unsigned k = 0; k < n_node_1d; k++)
2611 dtest_values[0][k] * test_values[1][j] * test_values[2][
i];
2615 test_values[0][k] * dtest_values[1][j] * test_values[2][
i];
2619 test_values[0][k] * test_values[1][j] * dtest_values[2][
i];
2623 test_values[0][k] * test_values[1][j] * test_values[2][
i];
2638 template<
unsigned DIM>
2650 const unsigned el_dim = this->dim();
2656 std::ostringstream error_message_stream;
2659 error_message_stream <<
"Need 3D space-time elements for this to work!"
2664 OOMPH_CURRENT_FUNCTION,
2665 OOMPH_EXCEPTION_LOCATION);
2670 dshape_local_u_nst(
s, psi, dpsidx);
2677 this->local_to_eulerian_mapping(dpsidx, inverse_jacobian);
2680 this->transform_derivatives(inverse_jacobian, dpsidx);
2687 dtest_local_u_nst(
s, test, dtestdx);
2690 this->transform_derivatives(inverse_jacobian, dtestdx);
2701 template<
unsigned DIM>
2703 DIM>::dshape_and_dtest_eulerian_at_knot_nst(
const unsigned& ipt,
2710 const unsigned el_dim = DIM + 1;
2716 for (
unsigned i = 0;
i < el_dim;
i++)
2719 s[
i] = this->integral_pt()->knot(ipt,
i);
2723 return dshape_and_dtest_eulerian_nst(
s, psi, dpsidx, test, dtestdx);
2765 for (
unsigned i = 0;
i < 2;
i++)
2768 for (
unsigned j = 0; j < 2; j++)
2771 for (
unsigned k = 0; k < 2; k++)
2774 ppsi[4 *
i + 2 * j + k] = psi3[
i] * psi2[j] * psi1[k];
2778 dppsidx(4 *
i + 2 * j + k, 0) = psi3[
i] * psi2[j] * dpsi1[k];
2782 dppsidx(4 *
i + 2 * j + k, 1) = psi3[
i] * dpsi2[j] * psi1[k];
2786 dppsidx(4 *
i + 2 * j + k, 2) = dpsi3[
i] * psi2[j] * psi1[k];
2798 dshape_local_u_nst(
s, psi, dpsi);
2804 const double det = local_to_eulerian_mapping(dpsi, inverse_jacobian);
2808 transform_derivatives(inverse_jacobian, dppsidx);
2872 for (
unsigned i = 0;
i < 2;
i++)
2875 for (
unsigned j = 0; j < 2; j++)
2878 for (
unsigned k = 0; k < 2; k++)
2881 ptest[4 *
i + 2 * j + k] = test3[
i] * test2[j] * test1[k];
2885 dptestdx(4 *
i + 2 * j + k, 0) = test3[
i] * test2[j] * dtest1[k];
2889 dptestdx(4 *
i + 2 * j + k, 1) = test3[
i] * dtest2[j] * test1[k];
2893 dptestdx(4 *
i + 2 * j + k, 2) = dtest3[
i] * test2[j] * test1[k];
2905 dshape_local_u_nst(
s, psi, dpsi);
2911 const double det = local_to_eulerian_mapping(dpsi, inverse_jacobian);
2915 transform_derivatives(inverse_jacobian, dptestdx);
2936 dptest_eulerian(
s, ptest, dptestdx);
2939 return this->dpshape_eulerian(
s, ppsi, dppsidx);
2955 const unsigned& ipt,
2965 throw OomphLibError(
"Hasn't been implemented properly yet!",
2966 OOMPH_CURRENT_FUNCTION,
2967 OOMPH_EXCEPTION_LOCATION);
3002 for (
unsigned i = 0;
i < 2;
i++)
3005 for (
unsigned j = 0; j < 2; j++)
3008 for (
unsigned k = 0; k < 2; k++)
3011 psi[4 *
i + 2 * j + k] = psi3[
i] * psi2[j] * psi1[k];
3052 for (
unsigned i = 0;
i < 2;
i++)
3055 for (
unsigned j = 0; j < 2; j++)
3058 for (
unsigned k = 0; k < 2; k++)
3061 test[4 *
i + 2 * j + k] = test3[
i] * test2[j] * test1[k];
3071 template<
unsigned DIM>
3114 template<
class TAYLOR_HOOD_ELEMENT>
3136 if (fld < this->dim() - 1)
3139 unsigned nnod = this->nnode();
3142 for (
unsigned j = 0; j < nnod; j++)
3145 data_values.push_back(std::make_pair(this->node_pt(j), fld));
3153 unsigned Pconv_size = this->dim();
3156 for (
unsigned j = 0; j < Pconv_size; j++)
3159 unsigned vertex_index = this->Pconv[j];
3162 data_values.push_back(
3163 std::make_pair(this->node_pt(vertex_index), fld));
3187 if (fld == this->dim())
3190 return this->node_pt(0)->ntstorage();
3196 return this->node_pt(0)->ntstorage();
3206 return this->node_pt(0)->position_time_stepper_pt()->ntstorage();
3217 unsigned n_dim = this->dim();
3220 unsigned n_node = this->nnode();
3226 this->pshape_nst(
s, psi);
3232 Shape testf(n_node);
3235 DShape dpsifdx(n_node, n_dim);
3238 DShape dtestfdx(n_node, n_dim);
3241 double J = this->dshape_and_dtest_eulerian_nst(
3242 s, psif, dpsifdx, testf, dtestfdx);
3251 Shape testf(n_node);
3254 DShape dpsifdx(n_node, n_dim);
3257 DShape dtestfdx(n_node, n_dim);
3261 this->dshape_and_dtest_eulerian_nst(
s, psi, dpsifdx, testf, dtestfdx);
3272 const unsigned& fld,
3275 unsigned n_dim = this->dim();
3276 unsigned n_node = this->nnode();
3281 return this->interpolated_p_nst(
t,
s);
3287 unsigned u_nodal_index = this->u_index_nst(fld);
3293 this->shape_u_nst(
s, psi);
3296 double interpolated_u = 0.0;
3299 for (
unsigned l = 0; l < n_node; l++)
3301 interpolated_u += this->nodal_value(
t, l, u_nodal_index) * psi[l];
3303 return interpolated_u;
3311 if (fld == this->dim())
3313 return this->npres_nst();
3317 return this->nnode();
3325 if (fld == this->dim())
3327 return this->p_local_eqn(j);
3331 const unsigned u_nodal_index = this->u_index_nst(fld);
3332 return this->nodal_local_eqn(j, u_nodal_index);
3342 template<
class ELEMENT>
3354 template<
class ELEMENT>
A Class for the derivatives of shape functions The class design is essentially the same as Shape,...
A class that represents a collection of data; each Data object may contain many different individual ...
long & eqn_number(const unsigned &i)
Return the equation number of the i-th stored variable.
void pin(const unsigned &i)
Pin the i-th stored variable.
double * value_pt(const unsigned &i) const
Return the pointer to the i-the stored value. Typically this is required when direct access to the st...
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...
unsigned nvalue() const
Return number of values stored in data object (incl pinned ones).
/////////////////////////////////////////////////////////////////////////
FaceElements are elements that coincide with the faces of higher-dimensional "bulk" elements....
int & face_index()
Index of the face (a number that uniquely identifies the face in the element)
void outer_unit_normal(const Vector< double > &s, Vector< double > &unit_normal) const
Compute outer unit normal at the specified local coordinate.
double interpolated_x(const Vector< double > &s, const unsigned &i) const
Return FE interpolated coordinate x[i] at local coordinate s. Overloaded to get information from bulk...
FaceGeometry()
Constructor; empty.
FaceGeometry()
Constructor; empty.
//////////////////////////////////////////////////////////////////// ////////////////////////////////...
A general Finite Element class.
virtual unsigned nplot_points_paraview(const unsigned &nplot) const
Return the number of actual plot points for paraview plot with parameter nplot. Broken virtual; can b...
virtual void local_coordinate_of_node(const unsigned &j, Vector< double > &s) const
Get local coordinates of node j in the element; vector sets its own size (broken virtual)
Node *& node_pt(const unsigned &n)
Return a pointer to the local node n.
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...
virtual void output(std::ostream &outfile)
Output the element data — typically the values at the nodes in a format suitable for post-processing.
virtual void transform_derivatives(const DenseMatrix< double > &inverse_jacobian, DShape &dbasis) const
Convert derivative w.r.t.local coordinates to derivatives w.r.t the coordinates used to assemble the ...
double size() const
Calculate the size of the element (length, area, volume,...) in Eulerian computational coordinates....
virtual double interpolated_x(const Vector< double > &s, const unsigned &i) const
Return FE interpolated coordinate x[i] at local coordinate s.
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.
unsigned dim() const
Return the spatial dimension of the element, i.e. the number of local coordinates required to paramet...
unsigned nnode() const
Return the number of nodes.
void(* SteadyExactSolutionFctPt)(const Vector< double > &, Vector< double > &)
Function pointer for function that computes vector-valued steady "exact solution" as .
Integral *const & integral_pt() const
Return the pointer to the integration scheme (const version)
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...
virtual double local_to_eulerian_mapping(const DShape &dpsids, DenseMatrix< double > &jacobian, DenseMatrix< double > &inverse_jacobian) const
Calculate the mapping from local to Eulerian coordinates, given the derivatives of the shape function...
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....
void(* UnsteadyExactSolutionFctPt)(const double &, const Vector< double > &, Vector< double > &)
Function pointer for function that computes Vector-valued time-dependent function as .
bool has_hanging_nodes() const
Return boolean to indicate if any of the element's nodes are geometrically hanging.
Helper class for elements that impose Robin boundary conditions on pressure advection diffusion probl...
virtual ~FpPressureAdvDiffRobinBCMixedOrderSpaceTimeElementBase()
Empty virtual destructor.
FpPressureAdvDiffRobinBCMixedOrderSpaceTimeElementBase()
Constructor.
virtual void fill_in_generic_residual_contribution_fp_press_adv_diff_robin_bc(Vector< double > &residuals, DenseMatrix< double > &jacobian, const unsigned &compute_jacobian_flag)=0
This function returns the residuals for the traction function. If compute_jacobian_flag=1 (or 0): do ...
//////////////////////////////////////////////////////////////////// ////////////////////////////////...
virtual void fill_in_generic_residual_contribution_fp_press_adv_diff_robin_bc(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 t...
FpPressureAdvDiffRobinBCMixedOrderSpaceTimeElement(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....
void fill_in_contribution_to_residuals(Vector< double > &residuals)
This function returns just the residuals.
~FpPressureAdvDiffRobinBCMixedOrderSpaceTimeElement()
Empty destructor.
void fill_in_contribution_to_jacobian(Vector< double > &residuals, DenseMatrix< double > &jacobian)
This function returns the residuals and the jacobian.
void output(std::ostream &outfile, const unsigned &nplot)
Output function: x,y,[z],u,v,[w],p in tecplot format.
void output(std::ostream &outfile)
Overload the output function.
static double Default_fd_jacobian_step
Double used for the default finite difference step in elemental jacobian calculations.
Data *& external_data_pt(const unsigned &i)
Return a pointer to i-th external data object.
Data *& internal_data_pt(const unsigned &i)
Return a pointer to i-th internal data object.
static DenseMatrix< double > Dummy_matrix
Empty dense matrix used as a dummy argument to combined residual and jacobian functions in the case w...
unsigned ninternal_data() const
Return the number of internal data objects.
unsigned add_external_data(Data *const &data_pt, const bool &fd=true)
Add a (pointer to an) external data object to the element and return its index (i....
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.
//////////////////////////////////////////////////////////////////// ////////////////////////////////...
Nodes are derived from Data, but, in addition, have a definite (Eulerian) position in a space of a gi...
An OomphLibError object which should be thrown when an run-time error is encountered....
Wrapper class for projectable elements. Adds "projectability" to the underlying ELEMENT.
///////////////////////////////////////////////////////////////// ///////////////////////////////////...
unsigned nvalue_of_field(const unsigned &fld)
Return number of values in field fld.
double jacobian_and_shape_of_field(const unsigned &fld, const Vector< double > &s, Shape &psi)
Return the Jacobian of the mapping and the shape functions of field fld at local coordinate s.
int local_equation(const unsigned &fld, const unsigned &j)
Return local equation number of value j in field fld.
unsigned nhistory_values_for_coordinate_projection()
Number of positional history values (Note: count includes current value!)
ProjectableTaylorHoodMixedOrderSpaceTimeElement()
Constructor [this was only required explicitly from gcc 4.5.2 onwards...].
unsigned nfields_for_projection()
Number of fields to be projected: dim+1, corresponding to velocity components and pressure.
Vector< std::pair< Data *, unsigned > > data_values_of_field(const unsigned &fld)
Specify the values associated with field fld. The information is returned in a vector of pairs which ...
unsigned nhistory_values_for_projection(const unsigned &fld)
Number of history values to be stored for fld-th field. Whatever the timestepper has set up for the v...
double get_field(const unsigned &t, const unsigned &fld, const Vector< double > &s)
Return interpolated field fld at local coordinate s, at time level t (t=0: present; t>0: history valu...
/////////////////////////////////////////////////////////////////////// /////////////////////////////...
/////////////////////////////////////////////////////////////////////////// /////////////////////////...
void pshape_nst(const Vector< double > &s, Shape &psi, Shape &test) const
Pressure shape and test functions at local coordinte s.
void ptest_nst(const Vector< double > &s, Shape &psi) const
Pressure test functions at local coordinate s.
unsigned npres_nst() const
Return number of pressure values.
double dpshape_and_dptest_eulerian_nst(const Vector< double > &s, Shape &ppsi, DShape &dppsidx, Shape &ptest, DShape &dptestdx) const
Pressure shape and test functions and their derivs w.r.t. to global coords at local coordinate s (tak...
static const unsigned Initial_Nvalue[]
Static array of ints to hold number of variables at node.
double dpshape_eulerian(const Vector< double > &s, Shape &ppsi, DShape &dppsidx) const
Pressure shape functions and their derivs w.r.t. to global coords at local coordinate s (taken from g...
static const unsigned Pconv[]
Static array of ints to hold conversion from pressure node numbers to actual node numbers.
void shape_u_nst(const Vector< double > &s, Shape &psi) const
DRAIG: Fill in later...
double p_nst(const unsigned &n_p) const
Access function for the pressure values at local pressure node n_p (const version)
void fix_pressure(const unsigned &p_dof, const double &p_value)
Pin p_dof-th pressure dof and set it to value specified by p_value.
double dshape_and_dtest_eulerian_at_knot_nst(const unsigned &ipt, Shape &psi, DShape &dpsidx, RankFourTensor< double > &d_dpsidx_dX, Shape &test, DShape &dtestdx, RankFourTensor< double > &d_dtestdx_dX, DenseMatrix< double > &djacobian_dX) const
Shape/test functions and derivs w.r.t. to global coords at integration point ipt; return Jacobian of ...
void output(FILE *file_pt, const unsigned &nplot)
Redirect output to NavierStokesMixedOrderEquations output.
bool is_pressure_node(const unsigned &n) const
Helper function which indicates whether or not the provided node is a pressure node.
void pshape_nst(const Vector< double > &s, Shape &psi) const
Pressure shape functions at local coordinate s.
void output(std::ostream &outfile)
Redirect output to NavierStokesMixedOrderEquations output.
double dshape_and_dtest_eulerian_nst(const Vector< double > &s, Shape &psi, DShape &dpsidx, Shape &test, DShape &dtestdx) const
Velocity shape and test functions and their derivs w.r.t. to global coords at local coordinate s (tak...
void identify_pressure_data(std::set< std::pair< Data *, unsigned >> &paired_pressure_data)
Add to the set paired_pressure_data pairs containing.
void identify_load_data(std::set< std::pair< Data *, unsigned >> &paired_load_data)
Add to the set paired_load_data pairs containing.
void dtest_local_u_nst(const Vector< double > &s, Shape &test, DShape &dtestdx) const
DRAIG: Fill in later...
~QTaylorHoodMixedOrderSpaceTimeElement()
virtual unsigned required_nvalue(const unsigned &n) const
Number of values (pinned or dofs) required at node n. Can be overwritten for hanging node version.
void build_fp_press_adv_diff_robin_bc_element(const unsigned &face_index)
Build FaceElements that apply the Robin boundary condition to the pressure advection diffusion proble...
QTaylorHoodMixedOrderSpaceTimeElement()
Constructor, no internal data points.
int p_local_eqn(const unsigned &n) const
Return the local equation numbers for the pressure values.
virtual int p_nodal_index_nst() const
Set the value at which the pressure is stored in the nodes.
void output(FILE *file_pt)
Redirect output to NavierStokesMixedOrderEquations output.
double dptest_eulerian(const Vector< double > &s, Shape &ptest, DShape &dptestdx) const
Pressure test functions and their derivs w.r.t. to global coords at local coordinate s (taken from ge...
double p_nst(const unsigned &t, const unsigned &n_p) const
Access function for the pressure values at local pressure node n_p (const version)
double dshape_and_dtest_eulerian_at_knot_nst(const unsigned &ipt, Shape &psi, DShape &dpsidx, Shape &test, DShape &dtestdx) const
Velocity shape and test functions and their derivs w.r.t. to global coords at local coordinate s (tak...
Node * pressure_node_pt(const unsigned &n_p)
Pointer to n_p-th pressure node.
void dshape_local_u_nst(const Vector< double > &s, Shape &psi, DShape &dpsidx) const
DRAIG: Fill in later...
int element_node_index_to_pressure_node_index(const unsigned &n)
Helper function which deduce the pressure node index associated with elemental node n (we're assuming...
void output(std::ostream &outfile, const unsigned &nplot)
Redirect output to NavierStokesMixedOrderEquations output.
////////////////////////////////////////////////////////////////// //////////////////////////////////...
////////////////////////////////////////////////////////////////// //////////////////////////////////...
RefineableElements are FiniteElements that may be subdivided into children to provide a better local ...
A Class for shape functions. In simple cases, the shape functions have only one index that can be tho...
A Class for nodes that deform elastically (i.e. position is an unknown in the problem)....
Data *const & variable_position_pt() const
Pointer to variable_position data (const version)
//////////////////////////////////////////////////////////////////// ////////////////////////////////...
void pin_all_non_pressure_dofs(std::map< Data *, std::vector< int >> &eqn_number_backup)
Pin all non-pressure dofs and backup eqn numbers.
virtual void dinterpolated_u_nst_ddata(const Vector< double > &s, const unsigned &i, Vector< double > &du_ddata, Vector< unsigned > &global_eqn_number)
Compute the derivatives of the i-th component of velocity at point s with respect to all data that ca...
void fill_in_contribution_to_hessian_vector_products(Vector< double > const &Y, DenseMatrix< double > const &C, DenseMatrix< double > &product)
Compute the hessian tensor vector products required to perform continuation of bifurcations analytica...
virtual double dpshape_and_dptest_eulerian_nst(const Vector< double > &s, Shape &ppsi, DShape &dppsidx, Shape &ptest, DShape &dptestdx) const =0
Compute the pressure shape and test functions and derivatives w.r.t. global coords at local coordinat...
virtual double dptest_eulerian(const Vector< double > &s, Shape &ptest, DShape &dptestdx) const =0
Pressure test functions and their derivs w.r.t. to global coords at local coordinate s (taken from ge...
double get_du_dt(const unsigned &n, const unsigned &i) const
i-th component of du/dt at local node n. Uses suitably interpolated value for hanging nodes....
virtual double get_source_nst(const double &time, const unsigned &ipt, const Vector< double > &x)
Calculate the source fct at a given time and Eulerian position.
double *& viscosity_ratio_pt()
Pointer to Viscosity Ratio.
virtual void fill_in_generic_dresidual_contribution_nst(double *const ¶meter_pt, Vector< double > &dres_dparam, DenseMatrix< double > &djac_dparam, DenseMatrix< double > &dmass_matrix_dparam, const unsigned &flag)
Compute the derivatives of the residuals for the Navier-Stokes equations with respect to a parameter ...
void get_vorticity(const Vector< double > &s, double &vorticity) const
Compute the vorticity vector at local coordinate s.
void full_output(std::ostream &outfile)
Full output function: x,y,[z],u,v,[w],p,du/dt,dv/dt,[dw/dt],dissipation in tecplot format....
virtual void dshape_local_u_nst(const Vector< double > &s, Shape &psi, DShape &dpsidx) const =0
DRAIG: Fill in later...
virtual double p_nst(const unsigned &t, const unsigned &n_p) const =0
Pressure at local pressure "node" n_p at time level t.
void output(FILE *file_pt)
C-style output function: x,y,[z],u,v,[w],p in tecplot format. The default number of plot points is fi...
void fill_in_contribution_to_jacobian(Vector< double > &residuals, DenseMatrix< double > &jacobian)
Compute the element's residual Vector and the jacobian matrix Virtual function can be overloaded by h...
virtual double dpshape_eulerian(const Vector< double > &s, Shape &ppsi, DShape &dppsidx) const =0
Pressure shape functions and their derivs w.r.t. to global coords at local coordinate s (taken from g...
double *& re_pt()
Pointer to Reynolds number.
NavierStokesMixedOrderBodyForceFctPt Body_force_fct_pt
Pointer to body force function.
double * re_st_pt() const
Pointer to Strouhal parameter (const. version)
void get_load(const Vector< double > &s, const Vector< double > &N, Vector< double > &load)
This implements a pure virtual function defined in the FSIFluidElement class. The function computes t...
void enable_ALE()
(Re-)enable ALE, i.e. take possible mesh motion into account when evaluating the time-derivative....
void delete_pressure_advection_diffusion_robin_elements()
Delete the FaceElements that apply the Robin boundary condition to the pressure advection diffusion p...
double(* NavierStokesMixedOrderPressureAdvDiffSourceFctPt)(const Vector< double > &x)
Function pointer to source function fct(x) for the pressure advection diffusion equation (only used d...
virtual double dshape_and_dtest_eulerian_nst(const Vector< double > &s, Shape &psi, DShape &dpsidx, Shape &test, DShape &dtestdx) const =0
Compute the shape functions and derivatives w.r.t. global coords at local coordinate s....
Vector< double > *& g_pt()
Pointer to Vector of gravitational components.
double interpolated_du_dt_nst(const Vector< double > &s, const unsigned &i) const
Return FE representation of function value du_i/dt(s) at local coordinate s.
double * ReSt_pt
Pointer to global Reynolds number x Strouhal number (=Womersley)
double *& re_invfr_pt()
Pointer to global inverse Froude number.
static int Pressure_not_stored_at_node
Static "magic" number that indicates that the pressure is not stored at a node.
unsigned n_u_nst() const
Return the number of velocity components (used in FluidInterfaceElements)
virtual void get_dresidual_dnodal_coordinates(RankThreeTensor< double > &dresidual_dnodal_coordinates)
Compute derivatives of elemental residual vector with respect to nodal coordinates....
double * Viscosity_Ratio_pt
Pointer to the viscosity ratio (relative to the viscosity used in the definition of the Reynolds numb...
virtual double dshape_and_dtest_eulerian_at_knot_nst(const unsigned &ipt, Shape &psi, DShape &dpsidx, Shape &test, DShape &dtestdx) const =0
Compute the shape functions and derivatives w.r.t. global coords at ipt-th integration point Return J...
void interpolated_u_nst(const Vector< double > &s, Vector< double > &velocity) const
Compute vector of FE interpolated velocity u at local coordinate s.
void fill_in_contribution_to_residuals(Vector< double > &residuals)
Compute the element's residual Vector.
double interpolated_p_nst(const unsigned &t, const Vector< double > &s) const
Return FE interpolated pressure at local coordinate s at time level t.
double(* NavierStokesMixedOrderSourceFctPt)(const double &time, const Vector< double > &x)
Function pointer to source function fct(t,x) (x is a Vector!)
const Vector< double > & g() const
Vector of gravitational components.
void fill_in_contribution_to_djacobian_and_dmass_matrix_dparameter(double *const ¶meter_pt, Vector< double > &dres_dparam, DenseMatrix< double > &djac_dparam, DenseMatrix< double > &dmass_matrix_dparam)
Add the element's contribution to its residuals vector, jacobian matrix and mass matrix.
double d_kin_energy_dt() const
Get integral of time derivative of kinetic energy over element.
const double & viscosity_ratio() const
Viscosity ratio for element: Element's viscosity relative to the viscosity used in the definition of ...
double u_nst(const unsigned &n, const unsigned &i) const
Velocity i at local node n. Uses suitably interpolated value for hanging nodes. The use of u_index_ns...
bool ReynoldsStrouhal_is_stored_as_external_data
Boolean to indicate whether or not the Strouhal value is stored as external data (if it's also an unk...
double * Re_pt
Pointer to global Reynolds number.
void scalar_value_fct_paraview(std::ofstream &file_out, const unsigned &i, const unsigned &nplot, const double &time, FiniteElement::UnsteadyExactSolutionFctPt exact_soln_pt) const
Write values of the i_field-th scalar field at the plot points. Needs to be implemented for each new ...
SpaceTimeNavierStokesMixedOrderEquations()
Constructor: NULL the body force and source function and make sure the ALE terms are included by defa...
double interpolated_u_nst(const Vector< double > &s, const unsigned &i) const
Return FE interpolated velocity u[i] at local coordinate s.
void fill_in_pressure_advection_diffusion_residuals(Vector< double > &residuals)
Compute the residuals for the associated pressure advection diffusion problem. Used by the Fp precond...
virtual void get_body_force_gradient_nst(const double &time, const unsigned &ipt, const Vector< double > &s, const Vector< double > &x, DenseMatrix< double > &d_body_force_dx)
Get gradient of body force term at (Eulerian) position x. This function is virtual to allow overloadi...
std::string scalar_name_paraview(const unsigned &i) const
Name of the i-th scalar field. Default implementation returns V1 for the first one,...
NavierStokesMixedOrderPressureAdvDiffSourceFctPt & source_fct_for_pressure_adv_diff()
Access function for the source-function pointer for pressure advection diffusion (used for validation...
virtual void pshape_nst(const Vector< double > &s, Shape &psi) const =0
Compute the pressure shape functions at local coordinate s.
bool is_reynolds_strouhal_stored_as_external_data() const
Are we storing the Strouhal number as external data?
static double Default_Physical_Ratio_Value
Static default value for the physical ratios (all are initialised to one)
double du_dt_nst(const unsigned &n, const unsigned &i) const
i-th component of du/dt at local node n. Uses suitably interpolated value for hanging nodes.
virtual void get_source_gradient_nst(const double &time, const unsigned &ipt, const Vector< double > &x, Vector< double > &gradient)
Get gradient of source term at (Eulerian) position x. This function is virtual to allow overloading i...
virtual void build_fp_press_adv_diff_robin_bc_element(const unsigned &face_index)=0
Build FaceElements that apply the Robin boundary condition to the pressure advection diffusion proble...
virtual void pshape_nst(const Vector< double > &s, Shape &psi, Shape &test) const =0
Compute the pressure shape and test functions at local coordinate s.
virtual void dtest_local_u_nst(const Vector< double > &s, Shape &test, DShape &dtestdx) const =0
DRAIG: Fill in later...
void output_veloc(std::ostream &outfile, const unsigned &nplot, const unsigned &t)
Output function: x,y,t,u,v in tecplot format. Use n_plot points in each coordinate direction at times...
bool ALE_is_disabled
Boolean flag to indicate if ALE formulation is disabled when time-derivatives are computed....
double *& density_ratio_pt()
Pointer to Density ratio.
static double Default_Physical_Constant_Value
Static default value for the physical constants (all initialised to zero)
void output_fct(std::ostream &outfile, const unsigned &nplot, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt)
Output exact solution specified via function pointer at a given number of plot points....
virtual int p_nodal_index_nst() const
Return the index at which the pressure is stored if it is stored at the nodes. If not stored at the n...
virtual void fill_in_generic_pressure_advection_diffusion_contribution_nst(Vector< double > &residuals, DenseMatrix< double > &jacobian, const unsigned &flag)
Compute the residuals for the associated pressure advection diffusion problem. Used by the Fp precond...
void get_pressure_and_velocity_mass_matrix_diagonal(Vector< double > &press_mass_diag, Vector< double > &veloc_mass_diag, const unsigned &which_one=0)
Compute the diagonal of the velocity/pressure mass matrices. If which one=0, both are computed,...
virtual void fill_in_generic_residual_contribution_nst(Vector< double > &residuals, DenseMatrix< double > &jacobian, DenseMatrix< double > &mass_matrix, const unsigned &flag)
Compute the residuals for the Navier-Stokes equations. Flag=1 (or 0): do (or don't) compute the Jacob...
double *& re_st_pt()
Pointer to ReSt number (can only assign to private member data)
virtual void fix_pressure(const unsigned &p_dof, const double &p_value)=0
Pin p_dof-th pressure dof and set it to value specified by p_value.
int Pinned_fp_pressure_eqn
Global eqn number of pressure dof that's pinned in pressure advection diffusion problem (defaults to ...
unsigned nscalar_paraview() const
Number of scalars/fields output by this element. Reimplements broken virtual function in base class.
double * ReInvFr_pt
Pointer to global Reynolds number x inverse Froude number (= Bond number / Capillary number)
const double & re() const
Reynolds number.
void scalar_value_paraview(std::ofstream &file_out, const unsigned &i, const unsigned &nplot) const
Write values of the i-th scalar field at the plot points. Needs to be implemented for each new specif...
virtual double p_nst(const unsigned &n_p) const =0
Pressure at local pressure "node" n_p Uses suitably interpolated value for hanging nodes.
double u_nst(const unsigned &t, const unsigned &n, const unsigned &i) const
Velocity i at local node n at timestep t (t=0: present; t>0: previous). Uses suitably interpolated va...
virtual unsigned npres_nst() const =0
Function to return number of pressure degrees of freedom.
double kin_energy() const
Get integral of kinetic energy over element.
void compute_norm(Vector< double > &norm)
Compute the vector norm of the FEM solution.
void strain_rate(const Vector< double > &s, DenseMatrix< double > &strain_rate) const
Strain-rate tensor: 1/2 (du_i/dx_j+du_j/dx_i)
virtual double interpolated_p_nst(const Vector< double > &s) const
Return FE interpolated pressure at local coordinate s.
const double & density_ratio() const
Density ratio for element: Element's density relative to the viscosity used in the definition of the ...
static Vector< double > Gamma
Vector to decide whether the stress-divergence form is used or not N.B. This needs to be public so th...
void(* NavierStokesMixedOrderBodyForceFctPt)(const double &time, const Vector< double > &x, Vector< double > &body_force)
Function pointer to body force function fct(t,x,f(x)) x is a Vector!
NavierStokesMixedOrderPressureAdvDiffSourceFctPt source_fct_for_pressure_adv_diff() const
Access function for the source-function pointer for pressure advection diffusion (used for validation...
void disable_ALE()
Disable ALE, i.e. assert the mesh is not moving – you do this at your own risk!
void output_vorticity(std::ostream &outfile, const unsigned &nplot)
Output function: x,y,t,omega in tecplot format. nplot points in each coordinate direction.
Vector< double > * G_pt
Pointer to global gravity Vector.
double dissipation() const
Return integral of dissipation over element.
NavierStokesMixedOrderPressureAdvDiffSourceFctPt Press_adv_diff_source_fct_pt
Pointer to source function pressure advection diffusion equation (only to be used during validation)
NavierStokesMixedOrderSourceFctPt Source_fct_pt
Pointer to volumetric source function.
const double & re_st() const
ReSt parameter (const. version)
NavierStokesMixedOrderBodyForceFctPt & body_force_fct_pt()
Access function for the body-force pointer.
double * Density_Ratio_pt
Pointer to the density ratio (relative to the density used in the definition of the Reynolds number)
void get_vorticity(const Vector< double > &s, Vector< double > &vorticity) const
Compute the vorticity vector at local coordinate s.
void output_pressure_advection_diffusion_robin_elements(std::ostream &outfile)
Output the FaceElements that apply the Robin boundary condition to the pressure advection diffusion p...
virtual void get_body_force_nst(const double &time, const unsigned &ipt, const Vector< double > &s, const Vector< double > &x, Vector< double > &result)
Calculate the body force at a given time and local and/or Eulerian position. This function is virtual...
NavierStokesMixedOrderBodyForceFctPt body_force_fct_pt() const
Access function for the body-force pointer. Const version.
void fill_in_contribution_to_dresiduals_dparameter(double *const ¶meter_pt, Vector< double > &dres_dparam)
Compute the element's residual Vector (differentiated w.r.t. a parameter)
NavierStokesMixedOrderSourceFctPt source_fct_pt() const
Access function for the source-function pointer. Const version.
virtual double dshape_and_dtest_eulerian_at_knot_nst(const unsigned &ipt, Shape &psi, DShape &dpsidx, RankFourTensor< double > &d_dpsidx_dX, Shape &test, DShape &dtestdx, RankFourTensor< double > &d_dtestdx_dX, DenseMatrix< double > &djacobian_dX) const =0
Shape/test functions and derivs w.r.t. to global coords at integration point ipt; return Jacobian of ...
double interpolated_u_nst(const unsigned &t, const Vector< double > &s, const unsigned &i) const
Return FE interpolated velocity u[i] at local coordinate s time level, t. Purposely broken for space-...
void store_reynolds_strouhal_as_external_data(Data *reynolds_strouhal_data_pt)
Function that tells us whether the period is stored as external data.
void output(std::ostream &outfile)
Output function: x,y,t,u,v,p in tecplot format. The default number of plot points is five.
void fill_in_contribution_to_jacobian_and_mass_matrix(Vector< double > &residuals, DenseMatrix< double > &jacobian, DenseMatrix< double > &mass_matrix)
Add the element's contribution to its residuals vector, jacobian matrix and mass matrix.
void get_traction(const Vector< double > &s, const Vector< double > &N, Vector< double > &traction)
Compute traction (on the viscous scale) exerted onto the fluid at local coordinate s....
static Vector< double > Default_Gravity_vector
Static default value for the gravity vector.
double pressure_integral() const
Integral of pressure over element.
const double & re_invfr() const
Global inverse Froude number.
void fill_in_pressure_advection_diffusion_jacobian(Vector< double > &residuals, DenseMatrix< double > &jacobian)
Compute the residuals and Jacobian for the associated pressure advection diffusion problem....
int & pinned_fp_pressure_eqn()
Global eqn number of pressure dof that's pinned in pressure adv diff problem.
NavierStokesMixedOrderSourceFctPt & source_fct_pt()
Access function for the source-function pointer.
virtual void shape_u_nst(const Vector< double > &s, Shape &psi) const =0
DRAIG: Fill in later...
void fill_in_contribution_to_djacobian_dparameter(double *const ¶meter_pt, Vector< double > &dres_dparam, DenseMatrix< double > &djac_dparam)
Compute the element's residual Vector and the jacobian matrix Virtual function can be overloaded by h...
virtual unsigned u_index_nst(const unsigned &i) const
Return the index at which the i-th unknown velocity component is stored. The default value,...
void compute_error(std::ostream &outfile, FiniteElement::UnsteadyExactSolutionFctPt exact_soln_pt, const double &time, double &error, double &norm)
Validate against exact solution at given time Solution is provided via function pointer....
Vector< FpPressureAdvDiffRobinBCMixedOrderSpaceTimeElementBase * > Pressure_advection_diffusion_robin_element_pt
Storage for FaceElements that apply Robin BC for pressure adv diff equation used in Fp preconditioner...
void point_output_data(const Vector< double > &s, Vector< double > &data)
Output solution in data vector at local cordinates s: x,y,z,u,v,p.
//////////////////////////////////////////////////////////////////// ////////////////////////////////...
virtual void pin_all_non_pressure_dofs(std::map< Data *, std::vector< int >> &eqn_number_backup)=0
Pin all non-pressure dofs and backup eqn numbers of all Data.
virtual ~TemplateFreeSpaceTimeNavierStokesMixedOrderEquationsBase()
Virtual destructor (empty)
virtual int & pinned_fp_pressure_eqn()=0
Global eqn number of pressure dof that's pinned in pressure adv diff problem.
virtual int p_nodal_index_nst() const =0
Return the index at which the pressure is stored if it is stored at the nodes. If not stored at the n...
virtual void fill_in_pressure_advection_diffusion_jacobian(Vector< double > &residuals, DenseMatrix< double > &jacobian)=0
Compute the residuals and Jacobian for the associated pressure advection diffusion problem....
virtual void get_pressure_and_velocity_mass_matrix_diagonal(Vector< double > &press_mass_diag, Vector< double > &veloc_mass_diag, const unsigned &which_one=0)=0
Compute the diagonal of the velocity/pressure mass matrices. If which one=0, both are computed,...
virtual void fill_in_pressure_advection_diffusion_residuals(Vector< double > &residuals)=0
Compute the residuals for the associated pressure advection diffusion problem. Used by the Fp precond...
virtual void build_fp_press_adv_diff_robin_bc_element(const unsigned &face_index)=0
Build FaceElements that apply the Robin boundary condition to the pressure advection diffusion proble...
TemplateFreeSpaceTimeNavierStokesMixedOrderEquationsBase()
Constructor (empty)
virtual int p_local_eqn(const unsigned &n) const =0
Access function for the local equation number information for the pressure. p_local_eqn[n] = local eq...
virtual void delete_pressure_advection_diffusion_robin_elements()=0
Delete the FaceElements that apply the Robin boundary condition to the pressure advection diffusion p...
void initialise(const _Tp &__value)
Iterate over all values and set to the desired value.
std::string string(const unsigned &i)
Return the i-th string or "" if the relevant string hasn't been defined.
void dshape< 3 >(const double &s, double *DPsi)
Derivatives of 1D shape functions specialised to quadratic order (3 Nodes)
void dshape< 2 >(const double &s, double *DPsi)
Derivatives of 1D shape functions specialised to linear order (2 Nodes)
void shape< 3 >(const double &s, double *Psi)
1D shape functions specialised to quadratic order (3 Nodes)
void shape< 2 >(const double &s, double *Psi)
1D shape functions specialised to linear order (2 Nodes)
std::string to_string(T object, unsigned float_precision=8)
Conversion function that should work for anything with operator<< defined (at least all basic types).
//////////////////////////////////////////////////////////////////// ////////////////////////////////...