Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
Public Member Functions | Private Member Functions | Private Attributes | Static Private Attributes | List of all members
Functions::StokesLSingularity Class Reference

#include <deal.II/base/flow_function.h>

Inheritance diagram for Functions::StokesLSingularity:
[legend]

Public Member Functions

 StokesLSingularity ()
 Constructor setting up some data.
 
- Public Member Functions inherited from Functions::FlowFunction< 2 >
 FlowFunction ()
 
virtual ~FlowFunction () override=default
 
void pressure_adjustment (double p)
 
virtual void vector_values (const std::vector< Point< dim >> &points, std::vector< std::vector< double >> &values) const override=0
 
virtual void vector_gradients (const std::vector< Point< dim >> &points, std::vector< std::vector< Tensor< 1, dim >>> &gradients) const override=0
 
virtual void vector_laplacians (const std::vector< Point< dim >> &points, std::vector< std::vector< double >> &values) const=0
 
virtual void vector_laplacian_list (const std::vector< Point< dim >> &points, std::vector< Vector< double >> &values) const override
 
- Public Member Functions inherited from Function< dim >
 Function (const unsigned int n_components=1, const time_type initial_time=0.0)
 
virtual ~Function () override=0
 
Functionoperator= (const Function &f)
 
virtual void value_list (const std::vector< Point< dim >> &points, std::vector< double > &values, const unsigned int component=0) const
 
virtual Tensor< 1, dim, double > gradient (const Point< dim > &p, const unsigned int component=0) const
 
virtual void vector_gradient (const Point< dim > &p, std::vector< Tensor< 1, dim, double >> &gradients) const
 
virtual void gradient_list (const std::vector< Point< dim >> &points, std::vector< Tensor< 1, dim, double >> &gradients, const unsigned int component=0) const
 
virtual void vector_gradients (const std::vector< Point< dim >> &points, std::vector< std::vector< Tensor< 1, dim, double >>> &gradients) const
 
virtual void vector_gradient_list (const std::vector< Point< dim >> &points, std::vector< std::vector< Tensor< 1, dim, double >>> &gradients) const
 
virtual double laplacian (const Point< dim > &p, const unsigned int component=0) const
 
virtual void vector_laplacian (const Point< dim > &p, Vector< double > &values) const
 
virtual void laplacian_list (const std::vector< Point< dim >> &points, std::vector< double > &values, const unsigned int component=0) const
 
virtual SymmetricTensor< 2, dim, double > hessian (const Point< dim > &p, const unsigned int component=0) const
 
virtual void vector_hessian (const Point< dim > &p, std::vector< SymmetricTensor< 2, dim, double >> &values) const
 
virtual void hessian_list (const std::vector< Point< dim >> &points, std::vector< SymmetricTensor< 2, dim, double >> &values, const unsigned int component=0) const
 
virtual void vector_hessian_list (const std::vector< Point< dim >> &points, std::vector< std::vector< SymmetricTensor< 2, dim, double >>> &values) const
 
std::size_t memory_consumption () const
 
- Public Member Functions inherited from FunctionTime< Number >
 FunctionTime (const Number initial_time=Number(0.0))
 
virtual ~FunctionTime ()=default
 
Number get_time () const
 
virtual void set_time (const Number new_time)
 
virtual void advance_time (const Number delta_t)
 
- Public Member Functions inherited from Subscriptor
 Subscriptor ()
 
 Subscriptor (const Subscriptor &)
 
 Subscriptor (Subscriptor &&) noexcept
 
virtual ~Subscriptor ()
 
Subscriptoroperator= (const Subscriptor &)
 
Subscriptoroperator= (Subscriptor &&) noexcept
 
void subscribe (std::atomic< bool > *const validity, const std::string &identifier="") const
 
void unsubscribe (std::atomic< bool > *const validity, const std::string &identifier="") const
 
unsigned int n_subscriptions () const
 
template<typename StreamType >
void list_subscribers (StreamType &stream) const
 
void list_subscribers () const
 
template<class Archive >
void serialize (Archive &ar, const unsigned int version)
 

Private Member Functions

double Psi (double phi) const
 The auxiliary function Psi.
 
double Psi_1 (double phi) const
 The derivative of Psi()
 
double Psi_2 (double phi) const
 The 2nd derivative of Psi()
 
double Psi_3 (double phi) const
 The 3rd derivative of Psi()
 
double Psi_4 (double phi) const
 The 4th derivative of Psi()
 

Private Attributes

const double omega
 The angle of the reentrant corner, set to 3*pi/2.
 
const double coslo
 Cosine of lambda times omega.
 
const double lp
 Auxiliary variable 1+lambda.
 
const double lm
 Auxiliary variable 1-lambda.
 

Static Private Attributes

static const double lambda = 0.54448373678246
 

Additional Inherited Members

- Public Types inherited from Function< dim >
using time_type = typename FunctionTime< typename numbers::NumberTraits< double >::real_type >::time_type
 
- Public Types inherited from FunctionTime< Number >
using time_type = Number
 
- Static Public Member Functions inherited from Subscriptor
static ::ExceptionBaseExcInUse (int arg1, std::string arg2, std::string arg3)
 
static ::ExceptionBaseExcNoSubscriber (std::string arg1, std::string arg2)
 
- Public Attributes inherited from Function< dim >
const unsigned int n_components
 
- Static Public Attributes inherited from Function< dim >
static const unsigned int dimension
 
- Protected Attributes inherited from Functions::FlowFunction< 2 >
double mean_pressure
 

Detailed Description

A singular solution to Stokes' equations on a 2d L-shaped domain.

This function satisfies \(-\triangle \mathbf{u} + \nabla p = 0\) and represents a typical singular solution around a reentrant corner of an L-shaped domain that can be created using GridGenerator::hyper_L(). The velocity vanishes on the two faces of the re-entrant corner and \(\nabla\mathbf{u}\) and \(p\) are singular at the origin while they are smooth in the rest of the domain because they can be written as a product of a smooth function and the term \(r^{\lambda-1}\) where \(r\) is the radius and \(\lambda \approx 0.54448\) is a fixed parameter.

Taken from Houston, Schötzau, Wihler, proceeding ENUMATH 2003.

Author
Guido Kanschat, 2007

Definition at line 247 of file flow_function.h.

Member Data Documentation

◆ lambda

const double Functions::StokesLSingularity::lambda = 0.54448373678246
staticprivate

The exponent of the radius, computed as the solution to \(\sin(\lambda\omega)+\lambda \sin(\omega)=0\)

Definition at line 284 of file flow_function.h.


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