oomph::BiharmonicFluidProblem< DIM > Class Template Reference

#include <biharmonic_problem.h>

+ Inheritance diagram for oomph::BiharmonicFluidProblem< DIM >:

Public Types

typedef void(* FluidBCFctPt) (const double &s, Vector< double > &u)
 
- Public Types inherited from oomph::Problem
typedef void(* SpatialErrorEstimatorFctPt) (Mesh *&mesh_pt, Vector< double > &elemental_error)
 Function pointer for spatial error estimator. More...
 
typedef void(* SpatialErrorEstimatorWithDocFctPt) (Mesh *&mesh_pt, Vector< double > &elemental_error, DocInfo &doc_info)
 Function pointer for spatial error estimator with doc. More...
 

Public Member Functions

 BiharmonicFluidProblem ()
 constructor More...
 
void actions_before_newton_solve ()
 actions before solve, performs self test More...
 
void actions_after_newton_solve ()
 action after solve More...
 
void doc_solution (DocInfo &doc_info, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt=0)
 
- Public Member Functions inherited from oomph::Problem
virtual void debug_hook_fct (const unsigned &i)
 
void set_analytic_dparameter (double *const &parameter_pt)
 
void unset_analytic_dparameter (double *const &parameter_pt)
 
bool is_dparameter_calculated_analytically (double *const &parameter_pt)
 
void set_analytic_hessian_products ()
 
void unset_analytic_hessian_products ()
 
bool are_hessian_products_calculated_analytically ()
 
void set_pinned_values_to_zero ()
 
bool distributed () const
 
virtual void actions_before_adapt ()
 
virtual void actions_after_adapt ()
 Actions that are to be performed after a mesh adaptation. More...
 
OomphCommunicatorcommunicator_pt ()
 access function to the oomph-lib communicator More...
 
const OomphCommunicatorcommunicator_pt () const
 access function to the oomph-lib communicator, const version More...
 
 Problem ()
 
 Problem (const Problem &dummy)=delete
 Broken copy constructor. More...
 
void operator= (const Problem &)=delete
 Broken assignment operator. More...
 
virtual ~Problem ()
 Virtual destructor to clean up memory. More...
 
Mesh *& mesh_pt ()
 Return a pointer to the global mesh. More...
 
Mesh *const & mesh_pt () const
 Return a pointer to the global mesh (const version) More...
 
Mesh *& mesh_pt (const unsigned &imesh)
 
Mesh *const & mesh_pt (const unsigned &imesh) const
 Return a pointer to the i-th submesh (const version) More...
 
unsigned nsub_mesh () const
 Return number of submeshes. More...
 
unsigned add_sub_mesh (Mesh *const &mesh_pt)
 
void flush_sub_meshes ()
 
void build_global_mesh ()
 
void rebuild_global_mesh ()
 
LinearSolver *& linear_solver_pt ()
 Return a pointer to the linear solver object. More...
 
LinearSolver *const & linear_solver_pt () const
 Return a pointer to the linear solver object (const version) More...
 
LinearSolver *& mass_matrix_solver_for_explicit_timestepper_pt ()
 
LinearSolvermass_matrix_solver_for_explicit_timestepper_pt () const
 
EigenSolver *& eigen_solver_pt ()
 Return a pointer to the eigen solver object. More...
 
EigenSolver *const & eigen_solver_pt () const
 Return a pointer to the eigen solver object (const version) More...
 
Time *& time_pt ()
 Return a pointer to the global time object. More...
 
Timetime_pt () const
 Return a pointer to the global time object (const version). More...
 
doubletime ()
 Return the current value of continuous time. More...
 
double time () const
 Return the current value of continuous time (const version) More...
 
TimeStepper *& time_stepper_pt ()
 
const TimeSteppertime_stepper_pt () const
 
TimeStepper *& time_stepper_pt (const unsigned &i)
 Return a pointer to the i-th timestepper. More...
 
ExplicitTimeStepper *& explicit_time_stepper_pt ()
 Return a pointer to the explicit timestepper. More...
 
unsigned long set_timestepper_for_all_data (TimeStepper *const &time_stepper_pt, const bool &preserve_existing_data=false)
 
virtual void shift_time_values ()
 Shift all values along to prepare for next timestep. More...
 
AssemblyHandler *& assembly_handler_pt ()
 Return a pointer to the assembly handler object. More...
 
AssemblyHandler *const & assembly_handler_pt () const
 Return a pointer to the assembly handler object (const version) More...
 
doubleminimum_dt ()
 Access function to min timestep in adaptive timestepping. More...
 
doublemaximum_dt ()
 Access function to max timestep in adaptive timestepping. More...
 
unsignedmax_newton_iterations ()
 Access function to max Newton iterations before giving up. More...
 
void problem_is_nonlinear (const bool &prob_lin)
 Access function to Problem_is_nonlinear. More...
 
doublemax_residuals ()
 
booltime_adaptive_newton_crash_on_solve_fail ()
 Access function for Time_adaptive_newton_crash_on_solve_fail. More...
 
doublenewton_solver_tolerance ()
 
void add_time_stepper_pt (TimeStepper *const &time_stepper_pt)
 
void set_explicit_time_stepper_pt (ExplicitTimeStepper *const &explicit_time_stepper_pt)
 
void initialise_dt (const double &dt)
 
void initialise_dt (const Vector< double > &dt)
 
Data *& global_data_pt (const unsigned &i)
 Return a pointer to the the i-th global data object. More...
 
void add_global_data (Data *const &global_data_pt)
 
void flush_global_data ()
 
LinearAlgebraDistribution *const & dof_distribution_pt () const
 Return the pointer to the dof distribution (read-only) More...
 
unsigned long ndof () const
 Return the number of dofs. More...
 
unsigned ntime_stepper () const
 Return the number of time steppers. More...
 
unsigned nglobal_data () const
 Return the number of global data values. More...
 
unsigned self_test ()
 Self-test: Check meshes and global data. Return 0 for OK. More...
 
void enable_store_local_dof_pt_in_elements ()
 
void disable_store_local_dof_pt_in_elements ()
 
unsigned long assign_eqn_numbers (const bool &assign_local_eqn_numbers=true)
 
void describe_dofs (std::ostream &out= *(oomph_info.stream_pt())) const
 
void enable_discontinuous_formulation ()
 
void disable_discontinuous_formulation ()
 
void get_dofs (DoubleVector &dofs) const
 
void get_dofs (const unsigned &t, DoubleVector &dofs) const
 Return vector of the t'th history value of all dofs. More...
 
void set_dofs (const DoubleVector &dofs)
 Set the values of the dofs. More...
 
void set_dofs (const unsigned &t, DoubleVector &dofs)
 Set the history values of the dofs. More...
 
void set_dofs (const unsigned &t, Vector< double * > &dof_pt)
 
void add_to_dofs (const double &lambda, const DoubleVector &increment_dofs)
 Add lambda x incremenet_dofs[l] to the l-th dof. More...
 
doubleglobal_dof_pt (const unsigned &i)
 
doubledof (const unsigned &i)
 i-th dof in the problem More...
 
double dof (const unsigned &i) const
 i-th dof in the problem (const version) More...
 
double *& dof_pt (const unsigned &i)
 Pointer to i-th dof in the problem. More...
 
doubledof_pt (const unsigned &i) const
 Pointer to i-th dof in the problem (const version) More...
 
virtual void get_inverse_mass_matrix_times_residuals (DoubleVector &Mres)
 
virtual void get_dvaluesdt (DoubleVector &f)
 
virtual void get_residuals (DoubleVector &residuals)
 Get the total residuals Vector for the problem. More...
 
virtual void get_jacobian (DoubleVector &residuals, DenseDoubleMatrix &jacobian)
 
virtual void get_jacobian (DoubleVector &residuals, CRDoubleMatrix &jacobian)
 
virtual void get_jacobian (DoubleVector &residuals, CCDoubleMatrix &jacobian)
 
virtual void get_jacobian (DoubleVector &residuals, SumOfMatrices &jacobian)
 
void get_fd_jacobian (DoubleVector &residuals, DenseMatrix< double > &jacobian)
 Get the full Jacobian by finite differencing. More...
 
void get_derivative_wrt_global_parameter (double *const &parameter_pt, DoubleVector &result)
 
void get_hessian_vector_products (DoubleVectorWithHaloEntries const &Y, Vector< DoubleVectorWithHaloEntries > const &C, Vector< DoubleVectorWithHaloEntries > &product)
 
void solve_eigenproblem (const unsigned &n_eval, Vector< std::complex< double >> &eigenvalue, Vector< DoubleVector > &eigenvector, const bool &steady=true)
 Solve the eigenproblem. More...
 
void solve_eigenproblem (const unsigned &n_eval, Vector< std::complex< double >> &eigenvalue, const bool &steady=true)
 
virtual void get_eigenproblem_matrices (CRDoubleMatrix &mass_matrix, CRDoubleMatrix &main_matrix, const double &shift=0.0)
 
void assign_eigenvector_to_dofs (DoubleVector &eigenvector)
 Assign the eigenvector passed to the function to the dofs. More...
 
void add_eigenvector_to_dofs (const double &epsilon, const DoubleVector &eigenvector)
 
void store_current_dof_values ()
 Store the current values of the degrees of freedom. More...
 
void restore_dof_values ()
 Restore the stored values of the degrees of freedom. More...
 
void enable_jacobian_reuse ()
 
void disable_jacobian_reuse ()
 Disable recycling of Jacobian in Newton iteration. More...
 
bool jacobian_reuse_is_enabled ()
 Is recycling of Jacobian in Newton iteration enabled? More...
 
booluse_predictor_values_as_initial_guess ()
 
void newton_solve ()
 Use Newton method to solve the problem. More...
 
void enable_globally_convergent_newton_method ()
 enable globally convergent Newton method More...
 
void disable_globally_convergent_newton_method ()
 disable globally convergent Newton method More...
 
void newton_solve (unsigned const &max_adapt)
 
void steady_newton_solve (unsigned const &max_adapt=0)
 
void copy (Problem *orig_problem_pt)
 
virtual Problemmake_copy ()
 
virtual void read (std::ifstream &restart_file, bool &unsteady_restart)
 
virtual void read (std::ifstream &restart_file)
 
virtual void dump (std::ofstream &dump_file) const
 
void dump (const std::string &dump_file_name) const
 
void delete_all_external_storage ()
 
virtual void symmetrise_eigenfunction_for_adaptive_pitchfork_tracking ()
 
doublebifurcation_parameter_pt () const
 
void get_bifurcation_eigenfunction (Vector< DoubleVector > &eigenfunction)
 
void activate_fold_tracking (double *const &parameter_pt, const bool &block_solve=true)
 
void activate_bifurcation_tracking (double *const &parameter_pt, const DoubleVector &eigenvector, const bool &block_solve=true)
 
void activate_bifurcation_tracking (double *const &parameter_pt, const DoubleVector &eigenvector, const DoubleVector &normalisation, const bool &block_solve=true)
 
void activate_pitchfork_tracking (double *const &parameter_pt, const DoubleVector &symmetry_vector, const bool &block_solve=true)
 
void activate_hopf_tracking (double *const &parameter_pt, const bool &block_solve=true)
 
void activate_hopf_tracking (double *const &parameter_pt, const double &omega, const DoubleVector &null_real, const DoubleVector &null_imag, const bool &block_solve=true)
 
void deactivate_bifurcation_tracking ()
 
void reset_assembly_handler_to_default ()
 Reset the system to the standard non-augemented state. More...
 
double arc_length_step_solve (double *const &parameter_pt, const double &ds, const unsigned &max_adapt=0)
 
double arc_length_step_solve (Data *const &data_pt, const unsigned &data_index, const double &ds, const unsigned &max_adapt=0)
 
void reset_arc_length_parameters ()
 
intsign_of_jacobian ()
 
void explicit_timestep (const double &dt, const bool &shift_values=true)
 Take an explicit timestep of size dt. More...
 
void unsteady_newton_solve (const double &dt)
 
void unsteady_newton_solve (const double &dt, const bool &shift_values)
 
void unsteady_newton_solve (const double &dt, const unsigned &max_adapt, const bool &first, const bool &shift=true)
 
double doubly_adaptive_unsteady_newton_solve (const double &dt, const double &epsilon, const unsigned &max_adapt, const bool &first, const bool &shift=true)
 
double doubly_adaptive_unsteady_newton_solve (const double &dt, const double &epsilon, const unsigned &max_adapt, const unsigned &suppress_resolve_after_spatial_adapt_flag, const bool &first, const bool &shift=true)
 
double adaptive_unsteady_newton_solve (const double &dt_desired, const double &epsilon)
 
double adaptive_unsteady_newton_solve (const double &dt_desired, const double &epsilon, const bool &shift_values)
 
void assign_initial_values_impulsive ()
 
void assign_initial_values_impulsive (const double &dt)
 
void calculate_predictions ()
 Calculate predictions. More...
 
void enable_mass_matrix_reuse ()
 
void disable_mass_matrix_reuse ()
 
bool mass_matrix_reuse_is_enabled ()
 Return whether the mass matrix is being reused. More...
 
void refine_uniformly (const Vector< unsigned > &nrefine_for_mesh)
 
void refine_uniformly (const Vector< unsigned > &nrefine_for_mesh, DocInfo &doc_info)
 
void refine_uniformly_and_prune (const Vector< unsigned > &nrefine_for_mesh)
 
void refine_uniformly_and_prune (const Vector< unsigned > &nrefine_for_mesh, DocInfo &doc_info)
 
void refine_uniformly (DocInfo &doc_info)
 
void refine_uniformly_and_prune (DocInfo &doc_info)
 
void refine_uniformly ()
 
void refine_uniformly (const unsigned &i_mesh, DocInfo &doc_info)
 Do uniform refinement for submesh i_mesh with documentation. More...
 
void refine_uniformly (const unsigned &i_mesh)
 Do uniform refinement for submesh i_mesh without documentation. More...
 
void p_refine_uniformly (const Vector< unsigned > &nrefine_for_mesh)
 
void p_refine_uniformly (const Vector< unsigned > &nrefine_for_mesh, DocInfo &doc_info)
 
void p_refine_uniformly_and_prune (const Vector< unsigned > &nrefine_for_mesh)
 
void p_refine_uniformly_and_prune (const Vector< unsigned > &nrefine_for_mesh, DocInfo &doc_info)
 
void p_refine_uniformly (DocInfo &doc_info)
 
void p_refine_uniformly_and_prune (DocInfo &doc_info)
 
void p_refine_uniformly ()
 
void p_refine_uniformly (const unsigned &i_mesh, DocInfo &doc_info)
 Do uniform p-refinement for submesh i_mesh with documentation. More...
 
void p_refine_uniformly (const unsigned &i_mesh)
 Do uniform p-refinement for submesh i_mesh without documentation. More...
 
void refine_selected_elements (const Vector< unsigned > &elements_to_be_refined)
 
void refine_selected_elements (const Vector< RefineableElement * > &elements_to_be_refined_pt)
 
void refine_selected_elements (const unsigned &i_mesh, const Vector< unsigned > &elements_to_be_refined)
 
void refine_selected_elements (const unsigned &i_mesh, const Vector< RefineableElement * > &elements_to_be_refined_pt)
 
void refine_selected_elements (const Vector< Vector< unsigned >> &elements_to_be_refined)
 
void refine_selected_elements (const Vector< Vector< RefineableElement * >> &elements_to_be_refined_pt)
 
void p_refine_selected_elements (const Vector< unsigned > &elements_to_be_refined)
 
void p_refine_selected_elements (const Vector< PRefineableElement * > &elements_to_be_refined_pt)
 
void p_refine_selected_elements (const unsigned &i_mesh, const Vector< unsigned > &elements_to_be_refined)
 
void p_refine_selected_elements (const unsigned &i_mesh, const Vector< PRefineableElement * > &elements_to_be_refined_pt)
 
void p_refine_selected_elements (const Vector< Vector< unsigned >> &elements_to_be_refined)
 
void p_refine_selected_elements (const Vector< Vector< PRefineableElement * >> &elements_to_be_refined_pt)
 
unsigned unrefine_uniformly ()
 
unsigned unrefine_uniformly (const unsigned &i_mesh)
 
void p_unrefine_uniformly (DocInfo &doc_info)
 
void p_unrefine_uniformly (const unsigned &i_mesh, DocInfo &doc_info)
 Do uniform p-unrefinement for submesh i_mesh without documentation. More...
 
void adapt (unsigned &n_refined, unsigned &n_unrefined)
 
void adapt ()
 
void p_adapt (unsigned &n_refined, unsigned &n_unrefined)
 
void p_adapt ()
 
void adapt_based_on_error_estimates (unsigned &n_refined, unsigned &n_unrefined, Vector< Vector< double >> &elemental_error)
 
void adapt_based_on_error_estimates (Vector< Vector< double >> &elemental_error)
 
void get_all_error_estimates (Vector< Vector< double >> &elemental_error)
 
void doc_errors (DocInfo &doc_info)
 Get max and min error for all elements in submeshes. More...
 
void doc_errors ()
 Get max and min error for all elements in submeshes. More...
 
void enable_info_in_newton_solve ()
 
void disable_info_in_newton_solve ()
 Disable the output of information when in the newton solver. More...
 
- Public Member Functions inherited from oomph::ExplicitTimeSteppableObject
 ExplicitTimeSteppableObject ()
 Empty constructor. More...
 
 ExplicitTimeSteppableObject (const ExplicitTimeSteppableObject &)=delete
 Broken copy constructor. More...
 
void operator= (const ExplicitTimeSteppableObject &)=delete
 Broken assignment operator. More...
 
virtual ~ExplicitTimeSteppableObject ()
 Empty destructor. More...
 
virtual void actions_before_explicit_stage ()
 
virtual void actions_after_explicit_stage ()
 

Protected Member Functions

void impose_solid_boundary_on_edge (const unsigned &b, const double &psi=0)
 
void impose_traction_free_edge (const unsigned &b)
 
void impose_fluid_flow_on_edge (const unsigned &b, FluidBCFctPt u_imposed_fn)
 
- Protected Member Functions inherited from oomph::Problem
unsigned setup_element_count_per_dof ()
 
virtual void sparse_assemble_row_or_column_compressed (Vector< int * > &column_or_row_index, Vector< int * > &row_or_column_start, Vector< double * > &value, Vector< unsigned > &nnz, Vector< double * > &residual, bool compressed_row_flag)
 
virtual void actions_before_newton_convergence_check ()
 
virtual void actions_before_newton_step ()
 
virtual void actions_after_newton_step ()
 
virtual void actions_before_implicit_timestep ()
 
virtual void actions_after_implicit_timestep ()
 
virtual void actions_after_implicit_timestep_and_error_estimation ()
 
virtual void actions_before_explicit_timestep ()
 Actions that should be performed before each explicit time step. More...
 
virtual void actions_after_explicit_timestep ()
 Actions that should be performed after each explicit time step. More...
 
virtual void actions_before_read_unstructured_meshes ()
 
virtual void actions_after_read_unstructured_meshes ()
 
virtual void actions_after_change_in_global_parameter (double *const &parameter_pt)
 
virtual void actions_after_change_in_bifurcation_parameter ()
 
virtual void actions_after_parameter_increase (double *const &parameter_pt)
 
doubledof_derivative (const unsigned &i)
 
doubledof_current (const unsigned &i)
 
virtual void set_initial_condition ()
 
virtual double global_temporal_error_norm ()
 
unsigned newton_solve_continuation (double *const &parameter_pt)
 
unsigned newton_solve_continuation (double *const &parameter_pt, DoubleVector &z)
 
void calculate_continuation_derivatives (double *const &parameter_pt)
 
void calculate_continuation_derivatives (const DoubleVector &z)
 
void calculate_continuation_derivatives_fd (double *const &parameter_pt)
 
bool does_pointer_correspond_to_problem_data (double *const &parameter_pt)
 
void set_consistent_pinned_values_for_continuation ()
 

Private Attributes

unsigned Npoint_element
 

Additional Inherited Members

- Public Attributes inherited from oomph::Problem
bool Shut_up_in_newton_solve
 
- Static Public Attributes inherited from oomph::Problem
static bool Suppress_warning_about_actions_before_read_unstructured_meshes
 
- Protected Types inherited from oomph::Problem
enum  Assembly_method {
  Perform_assembly_using_vectors_of_pairs , Perform_assembly_using_two_vectors , Perform_assembly_using_maps , Perform_assembly_using_lists ,
  Perform_assembly_using_two_arrays
}
 Enumerated flags to determine which sparse assembly method is used. More...
 
- Protected Attributes inherited from oomph::Problem
Vector< Problem * > Copy_of_problem_pt
 
std::map< double *, boolCalculate_dparameter_analytic
 
bool Calculate_hessian_products_analytic
 
LinearAlgebraDistributionDof_distribution_pt
 
Vector< double * > Dof_pt
 Vector of pointers to dofs. More...
 
DoubleVectorWithHaloEntries Element_count_per_dof
 
double Relaxation_factor
 
double Newton_solver_tolerance
 
unsigned Max_newton_iterations
 Maximum number of Newton iterations. More...
 
unsigned Nnewton_iter_taken
 
Vector< doubleMax_res
 Maximum residuals at start and after each newton iteration. More...
 
double Max_residuals
 
bool Time_adaptive_newton_crash_on_solve_fail
 
bool Jacobian_reuse_is_enabled
 Is re-use of Jacobian in Newton iteration enabled? Default: false. More...
 
bool Jacobian_has_been_computed
 
bool Problem_is_nonlinear
 
bool Pause_at_end_of_sparse_assembly
 
bool Doc_time_in_distribute
 
unsigned Sparse_assembly_method
 
unsigned Sparse_assemble_with_arrays_initial_allocation
 
unsigned Sparse_assemble_with_arrays_allocation_increment
 
Vector< Vector< unsigned > > Sparse_assemble_with_arrays_previous_allocation
 
double Numerical_zero_for_sparse_assembly
 
double FD_step_used_in_get_hessian_vector_products
 
bool Mass_matrix_reuse_is_enabled
 
bool Mass_matrix_has_been_computed
 
bool Discontinuous_element_formulation
 
double Minimum_dt
 Minimum desired dt: if dt falls below this value, exit. More...
 
double Maximum_dt
 Maximum desired dt. More...
 
double DTSF_max_increase
 
double DTSF_min_decrease
 
double Minimum_dt_but_still_proceed
 
bool Scale_arc_length
 Boolean to control whether arc-length should be scaled. More...
 
double Desired_proportion_of_arc_length
 Proportion of the arc-length to taken by the parameter. More...
 
double Theta_squared
 
int Sign_of_jacobian
 Storage for the sign of the global Jacobian. More...
 
double Continuation_direction
 
double Parameter_derivative
 Storage for the derivative of the global parameter wrt arc-length. More...
 
double Parameter_current
 Storage for the present value of the global parameter. More...
 
bool Use_continuation_timestepper
 Boolean to control original or new storage of dof stuff. More...
 
unsigned Dof_derivative_offset
 
unsigned Dof_current_offset
 
Vector< doubleDof_derivative
 Storage for the derivative of the problem variables wrt arc-length. More...
 
Vector< doubleDof_current
 Storage for the present values of the variables. More...
 
double Ds_current
 Storage for the current step value. More...
 
unsigned Desired_newton_iterations_ds
 
double Minimum_ds
 Minimum desired value of arc-length. More...
 
bool Bifurcation_detection
 Boolean to control bifurcation detection via determinant of Jacobian. More...
 
bool Bisect_to_find_bifurcation
 Boolean to control wheter bisection is used to located bifurcation. More...
 
bool First_jacobian_sign_change
 Boolean to indicate whether a sign change has occured in the Jacobian. More...
 
bool Arc_length_step_taken
 Boolean to indicate whether an arc-length step has been taken. More...
 
bool Use_finite_differences_for_continuation_derivatives
 
OomphCommunicatorCommunicator_pt
 The communicator for this problem. More...
 
bool Always_take_one_newton_step
 
double Timestep_reduction_factor_after_nonconvergence
 
bool Keep_temporal_error_below_tolerance
 
- Static Protected Attributes inherited from oomph::Problem
static ContinuationStorageScheme Continuation_time_stepper
 Storage for the single static continuation timestorage object. More...
 

Detailed Description

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

Biharmonic Fluid Problem Class - describes stokes flow in 2D. Developed for the topologically rectangular Hermite Element Mesh. Contains functions allowing the following boundary conditions to be applied (on a given edge):

  • wall : v_n = 0 and v_t = 0 (psi must also be prescribed)
  • traction free : v_t = 0
  • flow : v_n and v_t are prescribed NOTE 1 : psi is the stream function
    • fluid velocity normal to boundary v_n = +/- dpsi/dt
    • fluid velocity tangential to boundary v_t = -/+ dpsi/dn NOTE 2 : when a solid wall boundary condition is applied to ensure that v_n = 0 the the streamfunction psi must also be prescribed (and constant)

Member Typedef Documentation

◆ FluidBCFctPt

template<unsigned DIM>
typedef void(* oomph::BiharmonicFluidProblem< DIM >::FluidBCFctPt) (const double &s, Vector< double > &u)

Definition of a dirichlet boundary condition function pointer. Takes the position along a boundary (s) in the macro element coordinate scheme and returns the fluid velocity normal (dpsi/dt) to the boundary (u[0]) and the fluid velocity tangential (dpsidn) to the boundary (u[1]).

Constructor & Destructor Documentation

◆ BiharmonicFluidProblem()

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

constructor

237  {
238  // initialise the number of non bulk elements
239  Npoint_element = 0;
240  }
unsigned Npoint_element
Definition: biharmonic_problem.h:297

References oomph::BiharmonicFluidProblem< DIM >::Npoint_element.

Member Function Documentation

◆ actions_after_newton_solve()

template<unsigned DIM>
void oomph::BiharmonicFluidProblem< DIM >::actions_after_newton_solve ( )
inlinevirtual

action after solve

Reimplemented from oomph::Problem.

260 {}

◆ actions_before_newton_solve()

template<unsigned DIM>
void oomph::BiharmonicFluidProblem< DIM >::actions_before_newton_solve ( )
inlinevirtual

actions before solve, performs self test

Reimplemented from oomph::Problem.

245  {
246 #ifdef PARANOID
247  if (0 == self_test())
248  {
249  oomph_info << "self test passed" << std::endl;
250  }
251  else
252  {
253  oomph_info << "self test failed" << std::endl;
254  }
255 #endif
256  }
unsigned self_test()
Self-test: Check meshes and global data. Return 0 for OK.
Definition: problem.cc:13276
OomphInfo oomph_info
Definition: oomph_definitions.cc:319

References oomph::oomph_info, and oomph::Problem::self_test().

◆ doc_solution()

template<unsigned DIM>
void oomph::BiharmonicFluidProblem< DIM >::doc_solution ( DocInfo doc_info,
FiniteElement::SteadyExactSolutionFctPt  exact_soln_pt = 0 
)

documents the solution, and if an exact solution is provided, then the error between the numerical and exact solution is presented

812  {
813  // create an output stream
814  std::ofstream some_file;
815  std::ostringstream filename;
816 
817  // Number of plot points: npts x npts
818  unsigned npts = 5;
819 
820  // Output solution
821  filename << doc_info.directory() << "/soln_" << doc_info.label() << ".dat";
822  some_file.open(filename.str().c_str());
823  mesh_pt()->output(some_file, npts);
824  some_file.close();
825 
826  // Output fluid velocity solution
827  filename.str("");
828  filename << doc_info.directory() << "/soln_velocity_" << doc_info.label()
829  << ".dat";
830  some_file.open(filename.str().c_str());
831  unsigned n_element = mesh_pt()->nelement();
832  for (unsigned i = 0; i < n_element - Npoint_element; i++)
833  {
834  BiharmonicElement<2>* biharmonic_element_pt =
835  dynamic_cast<BiharmonicElement<2>*>(mesh_pt()->element_pt(i));
836  biharmonic_element_pt->output_fluid_velocity(some_file, npts);
837  }
838  some_file.close();
839 
840  // if the exact solution is not provided
841  if (exact_soln_pt != 0)
842  {
843  // Output exact solution
844  filename.str("");
845  filename << doc_info.directory() << "/exact_soln_" << doc_info.label()
846  << ".dat";
847  some_file.open(filename.str().c_str());
848  mesh_pt()->output_fct(some_file, npts, exact_soln_pt);
849  some_file.close();
850 
851  // Doc error and return of the square of the L2 error
852  double error, norm;
853  filename.str("");
854  filename << doc_info.directory() << "/error_" << doc_info.label()
855  << ".dat";
856  some_file.open(filename.str().c_str());
857  mesh_pt()->compute_error(some_file, exact_soln_pt, error, norm);
858  some_file.close();
859 
860  // Doc L2 error and norm of solution
861  oomph_info << "\nNorm of error : " << sqrt(error) << std::endl;
862  oomph_info << "Norm of solution: " << sqrt(norm) << std::endl
863  << std::endl;
864  }
865  }
AnnoyingScalar sqrt(const AnnoyingScalar &x)
Definition: AnnoyingScalar.h:134
int i
Definition: BiCGSTAB_step_by_step.cpp:9
void output_fct(std::ostream &outfile, const unsigned &n_plot, FiniteElement::SteadyExactSolutionFctPt)
Output a given Vector function at f(n_plot) points in each element.
Definition: mesh.cc:2199
virtual void compute_error(std::ostream &outfile, FiniteElement::UnsteadyExactSolutionFctPt exact_soln_pt, const double &time, double &error, double &norm)
Definition: mesh.h:1140
GeneralisedElement *& element_pt(const unsigned long &e)
Return pointer to element e.
Definition: mesh.h:448
void output(std::ostream &outfile)
Output for all elements.
Definition: mesh.cc:2027
unsigned long nelement() const
Return number of elements in the mesh.
Definition: mesh.h:590
Mesh *& mesh_pt()
Return a pointer to the global mesh.
Definition: problem.h:1280
string filename
Definition: MergeRestartFiles.py:39
int error
Definition: calibrate.py:297

References oomph::DocInfo::directory(), calibrate::error, MergeRestartFiles::filename, i, oomph::DocInfo::label(), oomph::oomph_info, oomph::BiharmonicEquations< DIM >::output_fluid_velocity(), and sqrt().

◆ impose_fluid_flow_on_edge()

template<unsigned DIM>
void oomph::BiharmonicFluidProblem< DIM >::impose_fluid_flow_on_edge ( const unsigned b,
FluidBCFctPt  u_imposed_fn 
)
protected

Impose a prescribed fluid flow comprising the velocity normal to the boundary (u_imposed_fn[0]) and the velocity tangential to the boundary (u_imposed_fn[1])

660  {
661  // number of nodes on boundary b
662  unsigned n_node = mesh_pt()->nboundary_node(b);
663 
664  // fixed faced index for boundary
665  int face_index = mesh_pt()->face_index_at_boundary(b, 0);
666 
667  // Need to get the s_fixed_index
668  unsigned s_fixed_index = 0;
669  // Also set the edge sign
670  int edge_sign = 0;
671 
672  switch (face_index)
673  {
674  case -1:
675  s_fixed_index = 0;
676  edge_sign = -1;
677  break;
678 
679  case 1:
680  s_fixed_index = 0;
681  edge_sign = 1;
682  break;
683 
684  case -2:
685  s_fixed_index = 1;
686  edge_sign = 1;
687  break;
688 
689  case 2:
690  s_fixed_index = 1;
691  edge_sign = -1;
692  break;
693 
694  default:
695  throw OomphLibError("Face Index not +/-1 or +/-2: Need 2D QElements",
698  }
699 
700 
701  // node position along edge b in macro element boundary representation
702  // [-1,1]
703  Vector<double> s(1);
704 
705  // finite difference step
706  const double h = 10e-8;
707 
708  // loop over nodes on boundary b
709  for (unsigned n = 0; n < n_node; n++)
710  {
711  // get m_t and dm_t/ds_t for this node
712  Vector<double> m_N(2);
714 
715  // vectors for dx_i/ds_n and dx_i/ds_t
716  Vector<double> dxds_n(2);
717  Vector<double> dxds_t(2);
718  dxds_n[0] =
719  mesh_pt()->boundary_node_pt(b, n)->x_gen(1 + s_fixed_index, 0);
720  dxds_n[1] =
721  mesh_pt()->boundary_node_pt(b, n)->x_gen(1 + s_fixed_index, 1);
722  dxds_t[0] =
723  mesh_pt()->boundary_node_pt(b, n)->x_gen(2 - s_fixed_index, 0);
724  dxds_t[1] =
725  mesh_pt()->boundary_node_pt(b, n)->x_gen(2 - s_fixed_index, 1);
726 
727  // vector for d2xi/ds_n ds_t
728  Vector<double> d2xds_nds_t(2);
729  d2xds_nds_t[0] = mesh_pt()->boundary_node_pt(b, n)->x_gen(3, 0);
730  d2xds_nds_t[1] = mesh_pt()->boundary_node_pt(b, n)->x_gen(3, 1);
731 
732  // compute dn/ds_n
733  double dnds_n =
734  ((dxds_n[0] * dxds_t[1]) - (dxds_n[1] * dxds_t[0])) /
735  (sqrt(dxds_t[0] * dxds_t[0] + dxds_t[1] * dxds_t[1]) * edge_sign);
736 
737  // compute dt/ds_n
738  double dtds_n = ((dxds_n[0] * dxds_t[0]) + (dxds_n[1] * dxds_t[1])) /
739  sqrt(dxds_t[0] * dxds_t[0] + dxds_t[1] * dxds_t[1]);
740 
741  // compute dt/ds_t
742  double dtds_t = sqrt(dxds_t[0] * dxds_t[0] + dxds_t[1] * dxds_t[1]);
743 
744  // compute d2t/ds_nds_t
745  double d2tds_nds_t =
746  (dxds_t[0] * d2xds_nds_t[0] + dxds_t[1] * d2xds_nds_t[1]) /
747  sqrt(dxds_t[0] * dxds_t[0] + dxds_t[1] * dxds_t[1]);
748 
749  // get imposed velocities
750  Vector<double> u(2);
751  (*u_imposed_fn)(m_N[0], u);
752  u[0] *= edge_sign;
753  u[1] *= -edge_sign;
754 
755  // compute d/dm_t(dpsidn) by finite difference
756  double ddm_tdudn = 0;
757  double ddm_tdudt = 0;
758  Vector<double> u_L(2);
759  Vector<double> u_R(2);
760  // evaluate dudn_fn about current node
761  if (n == 0)
762  {
763  (*u_imposed_fn)(m_N[0], u_L);
764  (*u_imposed_fn)(m_N[0] + h, u_R);
765  }
766  else if (n == n_node - 1)
767  {
768  (*u_imposed_fn)(m_N[0] - h, u_L);
769  (*u_imposed_fn)(m_N[0], u_R);
770  }
771  else
772  {
773  (*u_imposed_fn)(m_N[0] - 0.5 * h, u_L);
774  (*u_imposed_fn)(m_N[0] + 0.5 * h, u_R);
775  }
776  // compute
777  ddm_tdudn = (u_R[1] - u_L[1]) / h;
778  ddm_tdudt = (u_R[0] - u_L[0]) / h;
779 
780  // compute du/ds_t
781  double duds_t = dtds_t * u[0];
782 
783  // compute du/ds_n
784  double duds_n = dnds_n * u[1] + dtds_n * u[0];
785 
786  // compute d2u/ds_n ds_t
787  double d2uds_nds_t = dnds_n * m_N[1] * ddm_tdudn + d2tds_nds_t * u[0] +
788  dtds_n * m_N[1] * ddm_tdudt;
789 
790  // pin du/ds_n dof and set value
791  mesh_pt()->boundary_node_pt(b, n)->pin(1 + s_fixed_index);
792  mesh_pt()->boundary_node_pt(b, n)->set_value(1 + s_fixed_index, duds_n);
793 
794  // pin du/ds_t dof and set value
795  mesh_pt()->boundary_node_pt(b, n)->pin(2 - s_fixed_index);
796  mesh_pt()->boundary_node_pt(b, n)->set_value(2 - s_fixed_index, duds_t);
797 
798  // pin du/ds_t dof and set value
799  mesh_pt()->boundary_node_pt(b, n)->pin(3);
800  mesh_pt()->boundary_node_pt(b, n)->set_value(3, d2uds_nds_t);
801  }
802  }
const unsigned n
Definition: CG3DPackingUnitTest.cpp:11
Array< double, 1, 3 > e(1./3., 0.5, 2.)
Scalar * b
Definition: benchVecAdd.cpp:17
void pin(const unsigned &i)
Pin the i-th stored variable.
Definition: nodes.h:385
void set_value(const unsigned &i, const double &value_)
Definition: nodes.h:271
unsigned long nboundary_node(const unsigned &ibound) const
Return number of nodes on a particular boundary.
Definition: mesh.h:833
int face_index_at_boundary(const unsigned &b, const unsigned &e) const
Definition: mesh.h:896
Node *& boundary_node_pt(const unsigned &b, const unsigned &n)
Return pointer to node n on boundary b.
Definition: mesh.h:493
double & x_gen(const unsigned &k, const unsigned &i)
Definition: nodes.h:1126
virtual void get_coordinates_on_boundary(const unsigned &b, const unsigned &k, Vector< double > &boundary_zeta)
Definition: nodes.cc:2379
RealScalar s
Definition: level1_cplx_impl.h:130
#define OOMPH_EXCEPTION_LOCATION
Definition: oomph_definitions.h:61
#define OOMPH_CURRENT_FUNCTION
Definition: oomph_definitions.h:86

References b, e(), n, OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION, s, and sqrt().

◆ impose_solid_boundary_on_edge()

template<unsigned DIM>
void oomph::BiharmonicFluidProblem< DIM >::impose_solid_boundary_on_edge ( const unsigned b,
const double psi = 0 
)
protected

Imposes a solid boundary on boundary b - no flow into boundary or along boundary v_n = 0 and v_t = 0. User must presribe the streamfunction psi to ensure dpsi/dt = 0 is imposed at all points on the boundary and not just at the nodes

Imposes a solid boundary - no flow into boundary or along boundary v_n = 0 and v_t = 0. User must presribe the streamfunction psi to ensure dpsi/dt = 0 is imposed at all points on the boundary and not just at the nodes

459  {
460  // number of nodes on boundary b
461  unsigned n_node = mesh_pt()->nboundary_node(b);
462 
463  // loop over nodes on boundary b
464  for (unsigned n = 0; n < n_node; n++)
465  {
466  // loop over DOFs to be pinned du/ds_n, du/ds_t and d2u/ds_nds_t
467  for (unsigned k = 1; k < 4; k++)
468  {
469  // pin and zero DOF k
470  mesh_pt()->boundary_node_pt(b, n)->pin(k);
471  mesh_pt()->boundary_node_pt(b, n)->set_value(k, 0.0);
472  }
473  mesh_pt()->boundary_node_pt(b, n)->pin(0);
474  mesh_pt()->boundary_node_pt(b, n)->set_value(0, psi);
475  }
476  }
char char char int int * k
Definition: level2_impl.h:374

References b, k, and n.

◆ impose_traction_free_edge()

template<unsigned DIM>
void oomph::BiharmonicFluidProblem< DIM >::impose_traction_free_edge ( const unsigned b)
protected

Impose a traction free edge - i.e. v_t = 0 or dpsi/dn = 0. In general dpsi/dn = 0 can only be imposed using equation elements to set the DOFs dpsi/ds_n, however in the special case of dt/ds_n = 0, then dpsi/ds_n = 0 and can be imposed using pinning - this is handled automatically in this function. For a more detailed description of the equations see the description of the class BiharmonicFluidBoundaryElement

Impose a traction free edge - i.e. v_t = 0 or dpsi/dn = 0. In general dpsi/dn = 0 can only be imposed using equation elements to set the DOFs dpsi/ds_n, however in the special case of dt/ds_n = 0, then dpsi/ds_n = 0 and can be imposed using pinning - this is handled automatically in this function. For a more detailed description of the equations see the description of the class : BiharmonicFluidBoundaryElement

490  {
491  // fixed faced index for boundary
492  int face_index = mesh_pt()->face_index_at_boundary(b, 0);
493 
494  // Need to get the s_fixed_index
495  unsigned s_fixed_index = 0;
496  switch (face_index)
497  {
498  case -1:
499  case 1:
500  s_fixed_index = 0;
501  break;
502 
503  case -2:
504  case 2:
505  s_fixed_index = 1;
506  break;
507 
508  default:
509  throw OomphLibError("Face Index not +/-1 or +/-2: Need 2D QElements",
512  }
513 
514  // create a point to a hijacked biharmonic element
515  Hijacked<BiharmonicElement<2>>* hijacked_element_pt;
516 
517  // vectors for dx_i/ds_n and dx_i/ds_t
518  Vector<double> dxds_n(2);
519  Vector<double> dxds_t(2);
520 
521  // number of nodes along edge
522  unsigned n_node = mesh_pt()->nboundary_node(b);
523 
524  // loop over boudnary nodes
525  for (unsigned n = 0; n < n_node; n++)
526  {
527  // get dx_i/ds_t and dx_i/ds_n at node n
528  dxds_n[0] =
529  mesh_pt()->boundary_node_pt(b, n)->x_gen(1 + s_fixed_index, 0);
530  dxds_n[1] =
531  mesh_pt()->boundary_node_pt(b, n)->x_gen(1 + s_fixed_index, 1);
532  dxds_t[0] =
533  mesh_pt()->boundary_node_pt(b, n)->x_gen(2 - s_fixed_index, 0);
534  dxds_t[1] =
535  mesh_pt()->boundary_node_pt(b, n)->x_gen(2 - s_fixed_index, 1);
536 
537  // compute dt/ds_n
538  double dtds_n = ((dxds_n[0] * dxds_t[0]) + (dxds_n[1] * dxds_t[1])) /
539  sqrt(dxds_t[0] * dxds_t[0] + dxds_t[1] * dxds_t[1]);
540 
541  // if dt/ds_n = 0 we can impose the traction free edge at this node by
542  // pinning dpsi/ds_n = 0, otherwise the equation elements
543  // BiharmonicFluidBoundaryElement are used
544  if (dtds_n == 0.0)
545  {
546  // pin dpsi/dn
547  mesh_pt()->boundary_node_pt(b, n)->pin(1 + s_fixed_index);
548  mesh_pt()->boundary_node_pt(b, n)->set_value(1 + s_fixed_index, 0.0);
549  }
550 
551  // otherwise impose equation elements
552  else
553  {
554  // hijack DOFs in element either side of node
555  // boundary 0
556  if (b == 0)
557  {
558  // hijack DOFs in left element
559  if (n > 0)
560  {
561  hijacked_element_pt = dynamic_cast<Hijacked<BiharmonicElement<2>>*>(
562  mesh_pt()->boundary_element_pt(b, n - 1));
563  delete hijacked_element_pt->hijack_nodal_value(1,
564  1 + s_fixed_index);
565  }
566  // hijack DOFs in right element
567  if (n < (n_node - 1))
568  {
569  hijacked_element_pt = dynamic_cast<Hijacked<BiharmonicElement<2>>*>(
571  delete hijacked_element_pt->hijack_nodal_value(0,
572  1 + s_fixed_index);
573  }
574  }
575  // boundary 1
576  else if (b == 1)
577  {
578  // hijack DOFs in left element
579  if (n > 0)
580  {
581  hijacked_element_pt = dynamic_cast<Hijacked<BiharmonicElement<2>>*>(
582  mesh_pt()->boundary_element_pt(b, n - 1));
583  delete hijacked_element_pt->hijack_nodal_value(3,
584  1 + s_fixed_index);
585  }
586  // hijack DOFs in right element
587  if (n < n_node - 1)
588  {
589  hijacked_element_pt = dynamic_cast<Hijacked<BiharmonicElement<2>>*>(
591  delete hijacked_element_pt->hijack_nodal_value(1,
592  1 + s_fixed_index);
593  }
594  }
595 
596  // boundary 2
597  else if (b == 2)
598  {
599  // hijack DOFs in left element
600  if (n > 0)
601  {
602  hijacked_element_pt = dynamic_cast<Hijacked<BiharmonicElement<2>>*>(
603  mesh_pt()->boundary_element_pt(b, n - 1));
604  delete hijacked_element_pt->hijack_nodal_value(3,
605  1 + s_fixed_index);
606  }
607  if (n < n_node - 1)
608  {
609  // hijack DOFs in right element
610  hijacked_element_pt = dynamic_cast<Hijacked<BiharmonicElement<2>>*>(
612  delete hijacked_element_pt->hijack_nodal_value(2,
613  1 + s_fixed_index);
614  }
615  }
616  // boundary 3
617  else if (b == 3)
618  {
619  // hijack DOFs in left element
620  if (n > 0)
621  {
622  hijacked_element_pt = dynamic_cast<Hijacked<BiharmonicElement<2>>*>(
623  mesh_pt()->boundary_element_pt(b, n - 1));
624  delete hijacked_element_pt->hijack_nodal_value(2,
625  1 + s_fixed_index);
626  }
627  if (n < n_node - 1)
628  {
629  // hijack DOFs in right element
630  hijacked_element_pt = dynamic_cast<Hijacked<BiharmonicElement<2>>*>(
632  delete hijacked_element_pt->hijack_nodal_value(0,
633  1 + s_fixed_index);
634  }
635  }
636 
637  // create the boundary point element
638  BiharmonicFluidBoundaryElement* boundary_point_element_pt =
639  new BiharmonicFluidBoundaryElement(mesh_pt()->boundary_node_pt(b, n),
640  s_fixed_index);
641 
642  // add element to mesh
643  mesh_pt()->add_element_pt(boundary_point_element_pt);
644 
645  // increment number of non bulk elements
646  Npoint_element++;
647  }
648  }
649  }
FiniteElement * boundary_element_pt(const unsigned &b, const unsigned &e) const
Return pointer to e-th finite element on boundary b.
Definition: mesh.h:840
void add_element_pt(GeneralisedElement *const &element_pt)
Add a (pointer to) an element to the mesh.
Definition: mesh.h:617

References b, oomph::Hijacked< ELEMENT >::hijack_nodal_value(), n, OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION, and sqrt().

Member Data Documentation

◆ Npoint_element

template<unsigned DIM>
unsigned oomph::BiharmonicFluidProblem< DIM >::Npoint_element
private

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