BrethertonProblem< ELEMENT > Class Template Reference

Bretherton problem. More...

+ Inheritance diagram for BrethertonProblem< ELEMENT >:

Public Member Functions

 BrethertonProblem ()
 Constructor: More...
 
void actions_before_newton_convergence_check ()
 
void actions_before_newton_solve ()
 Update before solve: empty. More...
 
void actions_after_newton_solve ()
 
void fix_pressure (const unsigned &e, const unsigned &l, const double &pvalue)
 Fix pressure value l in element e to value p_value. More...
 
void activate_inflow_dependency ()
 
void parameter_study (const unsigned &nsteps)
 Run a parameter study; perform specified number of steps. More...
 
void doc_solution (DocInfo &doc_info)
 Doc the solution. More...
 
- 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 ()
 

Private Attributes

ELEMENT * Control_element_pt
 Pointer to control element. More...
 
ofstream Trace_file
 Trace file. More...
 
BrethertonSpineMesh< ELEMENT, SpineLineFluidInterfaceElement< ELEMENT > > * Bulk_mesh_pt
 Pointer to bulk mesh. More...
 

Additional Inherited Members

- 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 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 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_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 ()
 
- 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<class ELEMENT>
class BrethertonProblem< ELEMENT >

Bretherton problem.

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

Constructor & Destructor Documentation

◆ BrethertonProblem()

template<class ELEMENT >
BrethertonProblem< ELEMENT >::BrethertonProblem

Constructor:

Problem constructor.

552 {
553 
554  // Number of elements in the deposited film region
555  unsigned nx1=24;
556 
557  // Number of elements in the bottom part of the transition region
558  unsigned nx2=6;
559 
560  // Number of elements in channel region
561  unsigned nx3=12;
562 
563  // Number of elements in the vertical part of the transition region
564  // (=half the number of elements in the liquid filled region ahead
565  // of the finger tip)
566  unsigned nhalf=4;
567 
568  // Number of elements through thickness of deposited film
569  unsigned nh=3;
570 
571  // Thickness of deposited film
572  double h=0.1;
573 
574  // Create wall geom objects
575  GeomObject* lower_wall_pt=new StraightLine(-1.0);
576  GeomObject* upper_wall_pt=new StraightLine( 1.0);
577 
578  // Start coordinate on wall
579  double xi0=-4.0;
580 
581  // End of transition region on wall
582  double xi1=1.5;
583 
584  // End of liquid filled region (inflow) on wall
585  double xi2=5.0;
586 
587  //Now create the mesh
588  Bulk_mesh_pt = new BrethertonSpineMesh<ELEMENT,
590  (nx1,nx2,nx3,nh,nhalf,h,lower_wall_pt,upper_wall_pt,xi0,0.0,xi1,xi2);
591 
592  // Make bulk mesh the global mesh
594 
595  // Store the control element
597 
598 
599  // Set pointers to quantities that determine the inflow profile
600 
601  // Film thickness at outflow on lower wall:
604 
605  // Film thickness at outflow on upper wall:
606  unsigned last_spine=Bulk_mesh_pt->nfree_surface_spines()-1;
608  Bulk_mesh_pt->spine_pt(last_spine)->spine_height_pt()->value_pt(0);
609 
610  // Current y-position on lower wall at inflow
611  unsigned ibound=1;
613  Bulk_mesh_pt->boundary_node_pt(ibound,0)->x_pt(0,1);
614 
615  // Current y-position on upper wall at inflow
616  unsigned nnod=Bulk_mesh_pt->nboundary_node(ibound);
618  Bulk_mesh_pt->boundary_node_pt(ibound,nnod-1)->x_pt(0,1);
619 
620  // Activate dependency of inflow on spine heights at outflow
622 
623  // Set the boundary conditions for this problem: All nodes are
624  // free by default -- just pin the ones that have Dirichlet conditions
625  // here
626 
627  // No slip on boundaries 0 1 and 2
628  for(unsigned ibound=0;ibound<=2;ibound++)
629  {
630  unsigned num_nod=mesh_pt()->nboundary_node(ibound);
631  for (unsigned inod=0;inod<num_nod;inod++)
632  {
633  mesh_pt()->boundary_node_pt(ibound,inod)->pin(0);
634  mesh_pt()->boundary_node_pt(ibound,inod)->pin(1);
635  }
636  }
637 
638  // Uniform, parallel outflow on boundaries 3 and 5
639  for(unsigned ibound=3;ibound<=5;ibound+=2)
640  {
641  unsigned num_nod=mesh_pt()->nboundary_node(ibound);
642  for (unsigned inod=0;inod<num_nod;inod++)
643  {
644  mesh_pt()->boundary_node_pt(ibound,inod)->pin(0);
645  mesh_pt()->boundary_node_pt(ibound,inod)->pin(1);
646  }
647  }
648 
649  // Pin central spine
650  unsigned central_spine=(Bulk_mesh_pt->nfree_surface_spines()-1)/2;
651  Bulk_mesh_pt->spine_pt(central_spine)->spine_height_pt()->pin(0);
652 
653 
654  // No slip in moving frame on boundaries 0 and 2
655  for (unsigned ibound=0;ibound<=2;ibound+=2)
656  {
657  unsigned num_nod=mesh_pt()->nboundary_node(ibound);
658  for (unsigned inod=0;inod<num_nod;inod++)
659  {
660  // Parallel flow
661  mesh_pt()->boundary_node_pt(ibound,inod)->set_value(0,-1.0);
662  mesh_pt()->boundary_node_pt(ibound,inod)->set_value(1, 0.0);
663  }
664  }
665 
666 
667  // Parallel, uniform outflow on boundaries 3 and 5
668  for (unsigned ibound=3;ibound<=5;ibound+=2)
669  {
670  unsigned num_nod=mesh_pt()->nboundary_node(ibound);
671  for (unsigned inod=0;inod<num_nod;inod++)
672  {
673  // Parallel inflow
674  mesh_pt()->boundary_node_pt(ibound,inod)->set_value(0,-1.0);
675  mesh_pt()->boundary_node_pt(ibound,inod)->set_value(1, 0.0);
676  }
677  }
678 
679 
680 
681  //Complete the problem setup to make the elements fully functional
682 
683  //Loop over the elements in the layer
684  unsigned n_bulk=Bulk_mesh_pt->nbulk();
685  for(unsigned i=0;i<n_bulk;i++)
686  {
687  //Cast to a fluid element
688  ELEMENT *el_pt = dynamic_cast<ELEMENT*>(Bulk_mesh_pt->
689  bulk_element_pt(i));
690  //Set the Reynolds number, etc
691  el_pt->re_pt() = &Global_Physical_Variables::Re;
692  el_pt->re_st_pt() = &Global_Physical_Variables::ReSt;
693  el_pt->re_invfr_pt() = &Global_Physical_Variables::ReInvFr;
694  el_pt->g_pt() = &Global_Physical_Variables::G;
695  }
696 
697  //Loop over 1D interface elements and set capillary number
698  unsigned interface_element_pt_range = Bulk_mesh_pt->ninterface_element();
699  for(unsigned i=0;i<interface_element_pt_range;i++)
700  {
701  //Cast to a interface element
705 
706  //Set the Capillary number
708  }
709 
710  //Update the nodal positions, so that the domain is correct
711  //(and therefore the repeated node test passes)
713 
714  //Do equation numbering
715  cout << "Number of unknowns: " << assign_eqn_numbers() << std::endl;
716 
717 }
int i
Definition: BiCGSTAB_step_by_step.cpp:9
void activate_inflow_dependency()
Definition: bretherton.cc:726
ELEMENT * Control_element_pt
Pointer to control element.
Definition: bretherton.cc:535
BrethertonSpineMesh< ELEMENT, SpineLineFluidInterfaceElement< ELEMENT > > * Bulk_mesh_pt
Pointer to bulk mesh.
Definition: bretherton.cc:542
Definition: bretherton_spine_mesh.template.h:46
unsigned long nbulk() const
Number of elements in bulk.
Definition: bretherton_spine_mesh.template.h:99
FiniteElement *& interface_element_pt(const unsigned long &i)
Access functions for pointers to interface elements.
Definition: bretherton_spine_mesh.template.h:81
ELEMENT * control_element_pt()
Definition: bretherton_spine_mesh.template.h:187
unsigned nfree_surface_spines()
Definition: bretherton_spine_mesh.template.h:107
unsigned long ninterface_element() const
Number of elements on interface.
Definition: bretherton_spine_mesh.template.h:87
void pin(const unsigned &i)
Pin the i-th stored variable.
Definition: nodes.h:385
double * value_pt(const unsigned &i) const
Definition: nodes.h:324
void set_value(const unsigned &i, const double &value_)
Definition: nodes.h:271
double *& ca_pt()
Pointer to the Capillary number.
Definition: interface_elements.h:492
Definition: geom_objects.h:101
unsigned long nboundary_node(const unsigned &ibound) const
Return number of nodes on a particular boundary.
Definition: mesh.h:833
Node *& boundary_node_pt(const unsigned &b, const unsigned &n)
Return pointer to node n on boundary b.
Definition: mesh.h:493
double * x_pt(const unsigned &t, const unsigned &i)
Definition: nodes.h:1181
Mesh *& mesh_pt()
Return a pointer to the global mesh.
Definition: problem.h:1280
unsigned long assign_eqn_numbers(const bool &assign_local_eqn_numbers=true)
Definition: problem.cc:1989
Definition: specific_node_update_interface_elements.h:556
Spine *& spine_pt(const unsigned long &i)
Return the i-th spine in the mesh.
Definition: spines.h:623
void node_update(const bool &update_all_solid_nodes=false)
Definition: spines.cc:84
Data *& spine_height_pt()
Access function to Data object that stores the spine height.
Definition: spines.h:156
Definition: geom_objects.h:452
double ReSt
Womersley number.
Definition: rayleigh_instability.cc:56
double * H_lo_pt
Pointer to film thickness at outflow on the lower wall.
Definition: bretherton.cc:63
double * Y_up_pt
Pointer to y-position at inflow on the upper wall.
Definition: bretherton.cc:72
double * H_up_pt
Pointer to film thickness at outflow on the upper wall.
Definition: bretherton.cc:66
double * Y_lo_pt
Pointer to y-position at inflow on the lower wall.
Definition: bretherton.cc:69
double Ca
Capillary number.
Definition: fibre.cc:61
double ReInvFr
Product of Reynolds number and inverse of Froude number.
Definition: fibre.cc:58
double Re
Reynolds number.
Definition: fibre.cc:55
Vector< double > G(3)
Direction of gravity.
Definition: spherical_shell_convection.cc:62

References Global_Physical_Variables::Ca, oomph::FluidInterfaceElement::ca_pt(), Global_Physical_Variables::G, Global_Physical_Variables::H_lo_pt, Global_Physical_Variables::H_up_pt, i, Global_Physical_Variables::Re, Global_Physical_Variables::ReInvFr, Global_Physical_Variables::ReSt, Global_Physical_Variables::Y_lo_pt, and Global_Physical_Variables::Y_up_pt.

Member Function Documentation

◆ actions_after_newton_solve()

template<class ELEMENT >
void BrethertonProblem< ELEMENT >::actions_after_newton_solve ( )
inlinevirtual

Update after solve can remain empty, because the update is performed automatically after every Newton step.

Reimplemented from oomph::Problem.

509 {}

◆ actions_before_newton_convergence_check()

template<class ELEMENT >
void BrethertonProblem< ELEMENT >::actions_before_newton_convergence_check ( )
inlinevirtual

Spine heights/lengths are unknowns in the problem so their values get corrected during each Newton step. However, changing their value does not automatically change the nodal positions, so we need to update all of them

Reimplemented from oomph::Problem.

478  {
479  // Update
481 
482  // Apply inflow on boundary 1
483  unsigned ibound=1;
484  unsigned num_nod=mesh_pt()->nboundary_node(ibound);
485 
486  // Coordinate and velocity
487  Vector<double> x(2);
488  Vector<double> veloc(2);
489 
490  // Loop over all nodes
491  for (unsigned inod=0;inod<num_nod;inod++)
492  {
493  // Get nodal position
494  x[0]=mesh_pt()->boundary_node_pt(ibound,inod)->x(0);
495  x[1]=mesh_pt()->boundary_node_pt(ibound,inod)->x(1);
496 
497  // Get inflow profile
499  mesh_pt()->boundary_node_pt(ibound,inod)->set_value(0,veloc[0]);
500  mesh_pt()->boundary_node_pt(ibound,inod)->set_value(1,veloc[1]);
501  }
502  }
double & x(const unsigned &i)
Return the i-th nodal coordinate.
Definition: nodes.h:1060
void inflow(const Vector< double > &x, Vector< double > &veloc)
Definition: bretherton.cc:76
list x
Definition: plotDoE.py:28

References Global_Physical_Variables::inflow(), and plotDoE::x.

◆ actions_before_newton_solve()

template<class ELEMENT >
void BrethertonProblem< ELEMENT >::actions_before_newton_solve ( )
inlinevirtual

Update before solve: empty.

Reimplemented from oomph::Problem.

505 {}

◆ activate_inflow_dependency()

template<class ELEMENT >
void BrethertonProblem< ELEMENT >::activate_inflow_dependency

Activate the dependency of the inflow velocity on the spine heights at the outflow

Loop over elements on inflow boundary (1)

727 {
728 
729  // Spine heights that affect the inflow
730  Vector<Data*> outflow_spines(2);
731 
732  // First spine
733  outflow_spines[0]=Bulk_mesh_pt->spine_pt(0)->spine_height_pt();
734 
735  // Last proper spine
736  unsigned last_spine=Bulk_mesh_pt->nfree_surface_spines()-1;
737  outflow_spines[1]=Bulk_mesh_pt->spine_pt(last_spine)->spine_height_pt();;
738 
739 
740 
742  unsigned ibound=1;
743  unsigned nel=Bulk_mesh_pt->nboundary_element(ibound);
744  for (unsigned e=0;e<nel;e++)
745  {
746  // Get pointer to element
747  ELEMENT* el_pt=dynamic_cast<ELEMENT*>(Bulk_mesh_pt->
748  boundary_element_pt(ibound,e));
749  // Activate depency on inflow
750  el_pt->activate_inflow_dependency_on_external_data(
751  outflow_spines,ibound,&Global_Physical_Variables::inflow);
752  }
753 
754 }
Array< double, 1, 3 > e(1./3., 0.5, 2.)
unsigned nboundary_element(const unsigned &b) const
Return number of finite elements that are adjacent to boundary b.
Definition: mesh.h:878

References e(), and Global_Physical_Variables::inflow().

◆ doc_solution()

template<class ELEMENT >
void BrethertonProblem< ELEMENT >::doc_solution ( DocInfo doc_info)

Doc the solution.

763 {
764 
765  ofstream some_file;
766  char filename[100];
767 
768  // Number of plot points
769  unsigned npts=5;
770 
771  // Control coordinate: At bubble tip
772  Vector<double> s(2);
773  s[0]=1.0;
774  s[1]=1.0;
775 
776  // Last proper spine
777  unsigned last_spine=Bulk_mesh_pt->nfree_surface_spines()-1;
778 
779  // Doc
781  Trace_file << " " << Bulk_mesh_pt->spine_pt(0)->height();
782  Trace_file << " " << Bulk_mesh_pt->spine_pt(last_spine)->height();
783  Trace_file << " " << 1.3375*pow(Global_Physical_Variables::Ca,2.0/3.0);
784  Trace_file << " " << -Control_element_pt->interpolated_p_nst(s)*
786  Trace_file << " " << 1.0+3.8*pow(Global_Physical_Variables::Ca,2.0/3.0);
787  Trace_file << " " << Control_element_pt->interpolated_u_nst(s,0);
788  Trace_file << " " << Control_element_pt->interpolated_u_nst(s,1);
789  Trace_file << std::endl;
790 
791 
792  // Output solution
793  sprintf(filename,"%s/soln%i.dat",doc_info.directory().c_str(),
794  doc_info.number());
795  some_file.open(filename);
796  Bulk_mesh_pt->output(some_file,npts);
797  some_file.close();
798 
799 
800  // Output boundaries
801  sprintf(filename,"%s/boundaries%i.dat",doc_info.directory().c_str(),
802  doc_info.number());
803  some_file.open(filename);
804  Bulk_mesh_pt->output_boundaries(some_file);
805  some_file.close();
806 
807 }
ofstream Trace_file
Trace file.
Definition: bretherton.cc:538
std::string directory() const
Output directory.
Definition: oomph_utilities.h:524
unsigned & number()
Number used (e.g.) for labeling output files.
Definition: oomph_utilities.h:554
void output(std::ostream &outfile)
Output for all elements.
Definition: mesh.cc:2027
void output_boundaries(std::ostream &outfile)
Output the nodes on the boundaries (into separate tecplot zones)
Definition: mesh.cc:1064
double & height()
Access function to spine height.
Definition: spines.h:150
RealScalar s
Definition: level1_cplx_impl.h:130
EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC bfloat16 pow(const bfloat16 &a, const bfloat16 &b)
Definition: BFloat16.h:625
string filename
Definition: MergeRestartFiles.py:39

References Global_Physical_Variables::Ca, oomph::DocInfo::directory(), MergeRestartFiles::filename, oomph::DocInfo::number(), Eigen::bfloat16_impl::pow(), s, and oomph::Problem_Parameter::Trace_file.

◆ fix_pressure()

template<class ELEMENT >
void BrethertonProblem< ELEMENT >::fix_pressure ( const unsigned e,
const unsigned l,
const double pvalue 
)
inline

Fix pressure value l in element e to value p_value.

514  {
515  //Fix the pressure at that element
516  dynamic_cast<ELEMENT *>(Bulk_mesh_pt->element_pt(e))->
517  fix_pressure(l,pvalue);
518  }
void fix_pressure(const unsigned &e, const unsigned &l, const double &pvalue)
Fix pressure value l in element e to value p_value.
Definition: bretherton.cc:512
GeneralisedElement *& element_pt(const unsigned long &e)
Return pointer to element e.
Definition: mesh.h:448

References e().

◆ parameter_study()

template<class ELEMENT >
void BrethertonProblem< ELEMENT >::parameter_study ( const unsigned nsteps)

Run a parameter study; perform specified number of steps.

Parameter study.

817 {
818 
819  // Increase maximum residual
820  Problem::Max_residuals=100.0;
821 
822  // Set output directory
823  DocInfo doc_info;
824  doc_info.set_directory("RESLT");
825  doc_info.number()=0;
826 
827 
828  // Open trace file
829  char filename[100];
830  sprintf(filename,"%s/trace.dat",doc_info.directory().c_str());
831  Trace_file.open(filename);
832 
833  Trace_file << "VARIABLES=\"Ca\",";
834  Trace_file << "\"h<sub>bottom</sub>\",\"h<sub>too</sub>\",";
835  Trace_file << "\"h<sub>Bretherton</sub>\",\"p<sub>tip</sub>\",";
836  Trace_file << "\"p<sub>tip (Bretherton)</sub>\",\"u<sub>stag</sub>\",";
837  Trace_file << "\"v<sub>stag</sub>\"";
838  Trace_file << "\"<greek>a</greek><sub>bottom</sub>\",";
839  Trace_file << "\"<greek>a</greek><sub>top</sub>\"";
840  Trace_file << std::endl;
841 
842  // Initial scaling factor for Ca reduction
843  double factor=2.0;
844 
845  //Doc initial solution
846  doc_solution(doc_info);
847 
848 //Loop over the steps
849  for (unsigned step=1;step<=nsteps;step++)
850  {
851  cout << std::endl << "STEP " << step << std::endl;
852 
853  // Newton method tends to converge is initial stays bounded below
854  // a certain tolerance: This is cheaper to check than restarting
855  // the Newton method after divergence (we'd also need to back up
856  // the previous solution)
857  double maxres=100.0;
858  while (true)
859  {
860  cout << "Checking max. res for Ca = "
861  << Global_Physical_Variables::Ca << ". Max. residual: ";
862 
863  // Check the maximum residual
864  DoubleVector residuals;
867  get_residuals(residuals);
868  double max_res=residuals.max();
869  cout << max_res;
870 
871  // Check what to do
872  if (max_res>maxres)
873  {
874  cout << ". Too big!" << std::endl;
875  // Reset the capillary number
877  // Reduce the factor
878  factor=1.0+(factor-1.0)/1.5;
879  cout << "New reduction factor: " << factor << std::endl;
880  // Reduce the capillary number
882  // Try again
883  continue;
884  }
885  // Looks promising: Let's proceed to solve
886  else
887  {
888  cout << ". OK" << std::endl << std::endl;
889  // Break out of the Ca adjustment loop
890  break;
891  }
892  }
893 
894  //Solve
895  cout << "Solving for capillary number: "
896  << Global_Physical_Variables::Ca << std::endl;
897  newton_solve();
898 
899  // Doc solution
900  doc_info.number()++;
901  doc_solution(doc_info);
902 
903  // Reduce the capillary number
905  }
906 
907 }
void actions_before_newton_convergence_check()
Definition: bretherton.cc:477
void actions_before_newton_solve()
Update before solve: empty.
Definition: bretherton.cc:505
void doc_solution(DocInfo &doc_info)
Doc the solution.
Definition: bretherton.cc:762
Definition: oomph_utilities.h:499
void set_directory(const std::string &directory)
Definition: oomph_utilities.cc:298
Definition: double_vector.h:58
double max() const
returns the maximum coefficient
Definition: double_vector.cc:604
virtual void get_residuals(DoubleVector &residuals)
Get the total residuals Vector for the problem.
Definition: problem.cc:3714
void newton_solve()
Use Newton method to solve the problem.
Definition: problem.cc:8783

References Global_Physical_Variables::Ca, oomph::DocInfo::directory(), MergeRestartFiles::filename, GlobalFct::get_residuals(), oomph::DoubleVector::max(), oomph::DocInfo::number(), oomph::DocInfo::set_directory(), and oomph::Problem_Parameter::Trace_file.

Member Data Documentation

◆ Bulk_mesh_pt

template<class ELEMENT >
BrethertonSpineMesh<ELEMENT,SpineLineFluidInterfaceElement<ELEMENT> >* BrethertonProblem< ELEMENT >::Bulk_mesh_pt
private

Pointer to bulk mesh.

◆ Control_element_pt

template<class ELEMENT >
ELEMENT* BrethertonProblem< ELEMENT >::Control_element_pt
private

Pointer to control element.

◆ Trace_file

template<class ELEMENT >
ofstream BrethertonProblem< ELEMENT >::Trace_file
private

Trace file.


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