Public Member Functions | Private Member Functions | Private Attributes | Static Private Attributes | List of all members
oomph::BuoyantQCrouzeixRaviartElement< DIM > Class Template Reference

////////////////////////////////////////////////////////////////////// ////////////////////////////////////////////////////////////////////// More...

#include <boussinesq_elements.h>

Inheritance diagram for oomph::BuoyantQCrouzeixRaviartElement< DIM >:

Public Member Functions

 BuoyantQCrouzeixRaviartElement ()
 Constructor: call the underlying constructors and initialise the pointer to the Rayleigh number to point to the default value of 0.0. More...
 
void unfix_pressure (const unsigned &p_dof)
 Unpin p_dof-th pressure dof. More...
 
unsigned required_nvalue (const unsigned &n) const
 The required number of values stored at the nodes is the sum of the required values of the two single-physics elements. Note that this step is generic for any multi-physics element of this type. More...
 
const double & ra () const
 Access function for the Rayleigh number (const version) More...
 
double *& ra_pt ()
 Access function for the pointer to the Rayleigh number. More...
 
void disable_ALE ()
 Final override for disable ALE. More...
 
void enable_ALE ()
 Final override for enable ALE. More...
 
unsigned nscalar_paraview () const
 Number of scalars/fields output by this element. Broken virtual. Needs to be implemented for each new specific element type. Temporary dummy. More...
 
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. Broken virtual. Needs to be implemented for each new specific element type. Temporary dummy. More...
 
std::string scalar_name_paraview (const unsigned &i) const
 Name of the i-th scalar field. Default implementation returns V1 for the first one, V2 for the second etc. Can (should!) be overloaded with more meaningful names. More...
 
void output (std::ostream &outfile)
 Overload the standard output function with the broken default. More...
 
void output (std::ostream &outfile, const unsigned &nplot)
 Output function: Output x, y, u, v, p, theta at Nplot^DIM plot points. More...
 
void output (FILE *file_pt)
 C-style output function: Broken default. More...
 
void output (FILE *file_pt, const unsigned &n_plot)
 C-style output function: Broken default. More...
 
void output_fct (std::ostream &outfile, const unsigned &Nplot, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt)
 Output function for an exact solution: Broken default. More...
 
void output_fct (std::ostream &outfile, const unsigned &Nplot, const double &time, FiniteElement::UnsteadyExactSolutionFctPt exact_soln_pt)
 Output function for a time-dependent exact solution: Broken default. More...
 
unsigned u_index_adv_diff () const
 Overload the index at which the temperature variable is stored. We choose to store it after the fluid velocities. More...
 
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. Plot at a given number of plot points and compute L2 error and L2 norm of velocity solution over element Call the broken default. More...
 
void compute_error (std::ostream &outfile, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt, double &error, double &norm)
 Validate against exact solution. Solution is provided via function pointer. Plot at a given number of plot points and compute L2 error and L2 norm of velocity solution over element Call the broken default. More...
 
void get_wind_adv_diff (const unsigned &ipt, const Vector< double > &s, const Vector< double > &x, Vector< double > &wind) const
 Overload the wind function in the advection-diffusion equations. This provides the coupling from the Navier–Stokes equations to the advection-diffusion equations because the wind is the fluid velocity. More...
 
void get_body_force_nst (const double &time, const unsigned &ipt, const Vector< double > &s, const Vector< double > &x, Vector< double > &result)
 Overload the body force in the Navier-Stokes equations This provides the coupling from the advection-diffusion equations to the Navier–Stokes equations, the body force is the temperature multiplied by the Rayleigh number acting in the direction opposite to gravity. More...
 
void fill_in_contribution_to_residuals (Vector< double > &residuals)
 Calculate the element's contribution to the residual vector. Recall that fill_in_* functions MUST NOT initialise the entries in the vector to zero. This allows us to call the fill_in_* functions of the constituent single-physics elements sequentially, without wiping out any previously computed entries. More...
 
void fill_in_contribution_to_jacobian (Vector< double > &residuals, DenseMatrix< double > &jacobian)
 Compute the element's residual vector and the Jacobian matrix. Jacobian is computed by finite-differencing. More...
 
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. More...
 
void fill_in_off_diagonal_jacobian_blocks_by_fd (Vector< double > &residuals, DenseMatrix< double > &jacobian)
 Helper function to get the off-diagonal blocks of the Jacobian matrix by finite differences. More...
 
void fill_in_contribution_to_jacobian (Vector< double > &residuals, DenseMatrix< double > &jacobian)
 Compute the element's residual Vector and the Jacobian matrix. Use finite-differencing only for the off-diagonal blocks. More...
 
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. More...
 
void fill_in_off_diagonal_jacobian_blocks_analytic (Vector< double > &residuals, DenseMatrix< double > &jacobian)
 Helper function to get the off-diagonal blocks of the Jacobian matrix analytically. More...
 
void fill_in_contribution_to_jacobian (Vector< double > &residuals, DenseMatrix< double > &jacobian)
 Compute the element's residual Vector and the Jacobian matrix. Use analytic expressions for the off-diagonal blocks. More...
 
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. More...
 

Private Member Functions

double Default_Physical_Constant_Value
 Set the default physical value to be zero. More...
 
double Default_Physical_Constant_Value
 

Private Attributes

double * Ra_pt
 Pointer to a private data member, the Rayleigh number. More...
 

Static Private Attributes

static double Default_Physical_Constant_Value
 The static default value of the Rayleigh number. More...
 

Detailed Description

template<unsigned DIM>
class oomph::BuoyantQCrouzeixRaviartElement< DIM >

////////////////////////////////////////////////////////////////////// //////////////////////////////////////////////////////////////////////

A class that solves the Boussinesq approximation of the Navier–Stokes and energy equations by coupling two pre-existing classes. The QAdvectionDiffusionElement with bi-quadratic interpolation for the scalar variable (temperature) and QCrouzeixRaviartElement which solves the Navier–Stokes equations using bi-quadratic interpolation for the velocities and a discontinuous bi-linear interpolation for the pressure. Note that we are free to choose the order in which we store the variables at the nodes. In this case we choose to store the variables in the order fluid velocities followed by temperature. We must, therefore, overload the function AdvectionDiffusionEquations<DIM>::u_index_adv_diff() to indicate that the temperature is stored at the DIM-th position not the 0-th. We do not need to overload the corresponding function in the NavierStokesEquations<DIM> class because the velocities are stored first.

Definition at line 65 of file boussinesq_elements.h.

Constructor & Destructor Documentation

◆ BuoyantQCrouzeixRaviartElement()

template<unsigned DIM>
oomph::BuoyantQCrouzeixRaviartElement< DIM >::BuoyantQCrouzeixRaviartElement ( )
inline

Constructor: call the underlying constructors and initialise the pointer to the Rayleigh number to point to the default value of 0.0.

Definition at line 80 of file boussinesq_elements.h.

References oomph::BuoyantQCrouzeixRaviartElement< DIM >::Default_Physical_Constant_Value, and oomph::BuoyantQCrouzeixRaviartElement< DIM >::Ra_pt.

Member Function Documentation

◆ compute_error() [1/2]

template<unsigned DIM>
void oomph::BuoyantQCrouzeixRaviartElement< DIM >::compute_error ( std::ostream &  outfile,
FiniteElement::SteadyExactSolutionFctPt  exact_soln_pt,
double &  error,
double &  norm 
)
inline

Validate against exact solution. Solution is provided via function pointer. Plot at a given number of plot points and compute L2 error and L2 norm of velocity solution over element Call the broken default.

Definition at line 272 of file boussinesq_elements.h.

◆ compute_error() [2/2]

template<unsigned DIM>
void oomph::BuoyantQCrouzeixRaviartElement< DIM >::compute_error ( std::ostream &  outfile,
FiniteElement::UnsteadyExactSolutionFctPt  exact_soln_pt,
const double &  time,
double &  error,
double &  norm 
)
inline

Validate against exact solution at given time Solution is provided via function pointer. Plot at a given number of plot points and compute L2 error and L2 norm of velocity solution over element Call the broken default.

Definition at line 258 of file boussinesq_elements.h.

◆ Default_Physical_Constant_Value() [1/2]

double oomph::BuoyantQCrouzeixRaviartElement< 2 >::Default_Physical_Constant_Value
private

Set the default physical value to be zero.

Definition at line 744 of file boussinesq_elements.h.

◆ Default_Physical_Constant_Value() [2/2]

double oomph::BuoyantQCrouzeixRaviartElement< 3 >::Default_Physical_Constant_Value
private

Definition at line 747 of file boussinesq_elements.h.

◆ disable_ALE()

template<unsigned DIM>
void oomph::BuoyantQCrouzeixRaviartElement< DIM >::disable_ALE ( )
inline

Final override for disable ALE.

Definition at line 114 of file boussinesq_elements.h.

◆ enable_ALE()

template<unsigned DIM>
void oomph::BuoyantQCrouzeixRaviartElement< DIM >::enable_ALE ( )
inline

Final override for enable ALE.

Definition at line 122 of file boussinesq_elements.h.

◆ fill_in_contribution_to_jacobian() [1/3]

template<unsigned DIM>
void oomph::BuoyantQCrouzeixRaviartElement< DIM >::fill_in_contribution_to_jacobian ( Vector< double > &  residuals,
DenseMatrix< double > &  jacobian 
)
inline

Compute the element's residual vector and the Jacobian matrix. Jacobian is computed by finite-differencing.

Definition at line 341 of file boussinesq_elements.h.

◆ fill_in_contribution_to_jacobian() [2/3]

template<unsigned DIM>
void oomph::BuoyantQCrouzeixRaviartElement< DIM >::fill_in_contribution_to_jacobian ( Vector< double > &  residuals,
DenseMatrix< double > &  jacobian 
)
inline

Compute the element's residual Vector and the Jacobian matrix. Use finite-differencing only for the off-diagonal blocks.

Definition at line 516 of file boussinesq_elements.h.

References oomph::BuoyantQCrouzeixRaviartElement< DIM >::fill_in_off_diagonal_jacobian_blocks_by_fd().

◆ fill_in_contribution_to_jacobian() [3/3]

template<unsigned DIM>
void oomph::BuoyantQCrouzeixRaviartElement< DIM >::fill_in_contribution_to_jacobian ( Vector< double > &  residuals,
DenseMatrix< double > &  jacobian 
)
inline

Compute the element's residual Vector and the Jacobian matrix. Use analytic expressions for the off-diagonal blocks.

Definition at line 696 of file boussinesq_elements.h.

References oomph::BuoyantQCrouzeixRaviartElement< DIM >::fill_in_off_diagonal_jacobian_blocks_analytic().

◆ fill_in_contribution_to_jacobian_and_mass_matrix() [1/3]

template<unsigned DIM>
void oomph::BuoyantQCrouzeixRaviartElement< DIM >::fill_in_contribution_to_jacobian_and_mass_matrix ( Vector< double > &  residuals,
DenseMatrix< double > &  jacobian,
DenseMatrix< double > &  mass_matrix 
)
inline

Add the element's contribution to its residuals vector, jacobian matrix and mass matrix.

Definition at line 350 of file boussinesq_elements.h.

Referenced by oomph::BuoyantQCrouzeixRaviartElement< DIM >::fill_in_contribution_to_jacobian_and_mass_matrix().

◆ fill_in_contribution_to_jacobian_and_mass_matrix() [2/3]

template<unsigned DIM>
void oomph::BuoyantQCrouzeixRaviartElement< DIM >::fill_in_contribution_to_jacobian_and_mass_matrix ( Vector< double > &  residuals,
DenseMatrix< double > &  jacobian,
DenseMatrix< double > &  mass_matrix 
)
inline

◆ fill_in_contribution_to_jacobian_and_mass_matrix() [3/3]

template<unsigned DIM>
void oomph::BuoyantQCrouzeixRaviartElement< DIM >::fill_in_contribution_to_jacobian_and_mass_matrix ( Vector< double > &  residuals,
DenseMatrix< double > &  jacobian,
DenseMatrix< double > &  mass_matrix 
)
inline

◆ fill_in_contribution_to_residuals()

template<unsigned DIM>
void oomph::BuoyantQCrouzeixRaviartElement< DIM >::fill_in_contribution_to_residuals ( Vector< double > &  residuals)
inline

Calculate the element's contribution to the residual vector. Recall that fill_in_* functions MUST NOT initialise the entries in the vector to zero. This allows us to call the fill_in_* functions of the constituent single-physics elements sequentially, without wiping out any previously computed entries.

Definition at line 323 of file boussinesq_elements.h.

◆ fill_in_off_diagonal_jacobian_blocks_analytic()

template<unsigned DIM>
void oomph::BuoyantQCrouzeixRaviartElement< DIM >::fill_in_off_diagonal_jacobian_blocks_analytic ( Vector< double > &  residuals,
DenseMatrix< double > &  jacobian 
)
inline

◆ fill_in_off_diagonal_jacobian_blocks_by_fd()

template<unsigned DIM>
void oomph::BuoyantQCrouzeixRaviartElement< DIM >::fill_in_off_diagonal_jacobian_blocks_by_fd ( Vector< double > &  residuals,
DenseMatrix< double > &  jacobian 
)
inline

◆ get_body_force_nst()

template<unsigned DIM>
void oomph::BuoyantQCrouzeixRaviartElement< DIM >::get_body_force_nst ( const double &  time,
const unsigned &  ipt,
const Vector< double > &  s,
const Vector< double > &  x,
Vector< double > &  result 
)
inline

Overload the body force in the Navier-Stokes equations This provides the coupling from the advection-diffusion equations to the Navier–Stokes equations, the body force is the temperature multiplied by the Rayleigh number acting in the direction opposite to gravity.

Definition at line 298 of file boussinesq_elements.h.

References oomph::BuoyantQCrouzeixRaviartElement< DIM >::ra().

◆ get_wind_adv_diff()

template<unsigned DIM>
void oomph::BuoyantQCrouzeixRaviartElement< DIM >::get_wind_adv_diff ( const unsigned &  ipt,
const Vector< double > &  s,
const Vector< double > &  x,
Vector< double > &  wind 
) const
inline

Overload the wind function in the advection-diffusion equations. This provides the coupling from the Navier–Stokes equations to the advection-diffusion equations because the wind is the fluid velocity.

Definition at line 283 of file boussinesq_elements.h.

◆ nscalar_paraview()

template<unsigned DIM>
unsigned oomph::BuoyantQCrouzeixRaviartElement< DIM >::nscalar_paraview ( ) const
inline

Number of scalars/fields output by this element. Broken virtual. Needs to be implemented for each new specific element type. Temporary dummy.

Definition at line 132 of file boussinesq_elements.h.

◆ output() [1/4]

template<unsigned DIM>
void oomph::BuoyantQCrouzeixRaviartElement< DIM >::output ( FILE *  file_pt)
inline

C-style output function: Broken default.

Definition at line 216 of file boussinesq_elements.h.

◆ output() [2/4]

template<unsigned DIM>
void oomph::BuoyantQCrouzeixRaviartElement< DIM >::output ( FILE *  file_pt,
const unsigned &  n_plot 
)
inline

C-style output function: Broken default.

Definition at line 222 of file boussinesq_elements.h.

◆ output() [3/4]

template<unsigned DIM>
void oomph::BuoyantQCrouzeixRaviartElement< DIM >::output ( std::ostream &  outfile)
inline

Overload the standard output function with the broken default.

Definition at line 167 of file boussinesq_elements.h.

◆ output() [4/4]

template<unsigned DIM>
void oomph::BuoyantQCrouzeixRaviartElement< DIM >::output ( std::ostream &  outfile,
const unsigned &  nplot 
)
inline

Output function: Output x, y, u, v, p, theta at Nplot^DIM plot points.

Definition at line 175 of file boussinesq_elements.h.

◆ output_fct() [1/2]

template<unsigned DIM>
void oomph::BuoyantQCrouzeixRaviartElement< DIM >::output_fct ( std::ostream &  outfile,
const unsigned &  Nplot,
const double &  time,
FiniteElement::UnsteadyExactSolutionFctPt  exact_soln_pt 
)
inline

Output function for a time-dependent exact solution: Broken default.

Definition at line 238 of file boussinesq_elements.h.

◆ output_fct() [2/2]

template<unsigned DIM>
void oomph::BuoyantQCrouzeixRaviartElement< DIM >::output_fct ( std::ostream &  outfile,
const unsigned &  Nplot,
FiniteElement::SteadyExactSolutionFctPt  exact_soln_pt 
)
inline

Output function for an exact solution: Broken default.

Definition at line 228 of file boussinesq_elements.h.

◆ ra()

template<unsigned DIM>
const double& oomph::BuoyantQCrouzeixRaviartElement< DIM >::ra ( ) const
inline

◆ ra_pt()

template<unsigned DIM>
double*& oomph::BuoyantQCrouzeixRaviartElement< DIM >::ra_pt ( )
inline

Access function for the pointer to the Rayleigh number.

Definition at line 108 of file boussinesq_elements.h.

References oomph::BuoyantQCrouzeixRaviartElement< DIM >::Ra_pt.

◆ required_nvalue()

template<unsigned DIM>
unsigned oomph::BuoyantQCrouzeixRaviartElement< DIM >::required_nvalue ( const unsigned &  n) const
inline

The required number of values stored at the nodes is the sum of the required values of the two single-physics elements. Note that this step is generic for any multi-physics element of this type.

Definition at line 95 of file boussinesq_elements.h.

◆ scalar_name_paraview()

template<unsigned DIM>
std::string oomph::BuoyantQCrouzeixRaviartElement< DIM >::scalar_name_paraview ( const unsigned &  i) const
inline

Name of the i-th scalar field. Default implementation returns V1 for the first one, V2 for the second etc. Can (should!) be overloaded with more meaningful names.

Definition at line 160 of file boussinesq_elements.h.

◆ scalar_value_paraview()

template<unsigned DIM>
void oomph::BuoyantQCrouzeixRaviartElement< DIM >::scalar_value_paraview ( std::ofstream &  file_out,
const unsigned &  i,
const unsigned &  nplot 
) const
inline

Write values of the i-th scalar field at the plot points. Broken virtual. Needs to be implemented for each new specific element type. Temporary dummy.

Definition at line 146 of file boussinesq_elements.h.

◆ u_index_adv_diff()

template<unsigned DIM>
unsigned oomph::BuoyantQCrouzeixRaviartElement< DIM >::u_index_adv_diff ( ) const
inline

Overload the index at which the temperature variable is stored. We choose to store it after the fluid velocities.

Definition at line 248 of file boussinesq_elements.h.

Referenced by oomph::BuoyantQCrouzeixRaviartElement< DIM >::fill_in_off_diagonal_jacobian_blocks_analytic(), and oomph::BuoyantQCrouzeixRaviartElement< DIM >::fill_in_off_diagonal_jacobian_blocks_by_fd().

◆ unfix_pressure()

template<unsigned DIM>
void oomph::BuoyantQCrouzeixRaviartElement< DIM >::unfix_pressure ( const unsigned &  p_dof)
inline

Unpin p_dof-th pressure dof.

Definition at line 87 of file boussinesq_elements.h.

Member Data Documentation

◆ Default_Physical_Constant_Value

template<unsigned DIM>
double oomph::BuoyantQCrouzeixRaviartElement< DIM >::Default_Physical_Constant_Value
staticprivate

The static default value of the Rayleigh number.

Definition at line 74 of file boussinesq_elements.h.

Referenced by oomph::BuoyantQCrouzeixRaviartElement< DIM >::BuoyantQCrouzeixRaviartElement().

◆ Ra_pt

template<unsigned DIM>
double* oomph::BuoyantQCrouzeixRaviartElement< DIM >::Ra_pt
private

The documentation for this class was generated from the following file: