DropInChannelProblem< ELEMENT > Class Template Reference

Problem class to simulate viscous drop propagating along 2D channel. More...

+ Inheritance diagram for DropInChannelProblem< ELEMENT >:

Public Member Functions

 DropInChannelProblem ()
 Constructor. More...
 
 ~DropInChannelProblem ()
 Destructor. More...
 
void actions_before_adapt ()
 Actions before adapt: Wipe the mesh of free surface elements. More...
 
void actions_after_adapt ()
 Actions after adapt: Rebuild the mesh of free surface elements. More...
 
void actions_after_newton_solve ()
 Update the after solve (empty) More...
 
void actions_before_newton_solve ()
 Update the problem specs before solve. More...
 
void complete_problem_setup ()
 Set boundary conditions and complete the build of all elements. More...
 
void doc_solution (const std::string &comment="")
 Doc the solution. More...
 
void compute_error_estimate (double &max_err, double &min_err)
 Compute the error estimates and assign to elements for plotting. More...
 
void remove_volume_constraint ()
 Change the boundary conditions to remove the volume constraint. 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
 
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 Types

enum  {
  Inflow_boundary_id =0 , Upper_wall_boundary_id =1 , Outflow_boundary_id =2 , Bottom_wall_boundary_id =3 ,
  First_drop_boundary_id =4 , Second_drop_boundary_id =5
}
 Enumeration of channel boundaries. More...
 

Private Member Functions

void create_free_surface_elements ()
 Create free surface elements. More...
 
void delete_free_surface_elements ()
 Delete free surface elements. More...
 
void create_volume_constraint_elements ()
 Create elements that impose volume constraint on the drop. More...
 
void delete_volume_constraint_elements ()
 Delete volume constraint elements. More...
 

Private Attributes

MeshFree_surface_mesh_pt
 Pointers to mesh of free surface elements. More...
 
MeshVolume_constraint_mesh_pt
 Pointer to mesh containing elements that impose volume constraint. More...
 
RefineableSolidTriangleMesh< ELEMENT > * Fluid_mesh_pt
 Pointer to Fluid_mesh. More...
 
Vector< TriangleMeshPolygon * > Drop_polygon_pt
 Vector storing pointer to the drop polygons. More...
 
TriangleMeshPolygonOuter_boundary_polyline_pt
 Triangle mesh polygon for outer boundary. More...
 
DataDrop_pressure_data_pt
 Pointer to a global drop pressure datum. More...
 
VolumeConstraintElementVol_constraint_el_pt
 Pointer to element that imposes volume constraint for drop. More...
 
double Initial_value_for_drop_pressure
 Backed up drop pressure between adaptations. More...
 
ELEMENT * Hijacked_element_pt
 Pointer to hijacked element. More...
 
bool Use_volume_constraint
 Bool to indicate if volume constraint is applied (only for steady run) 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_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 ()
 
- 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 DropInChannelProblem< ELEMENT >

Problem class to simulate viscous drop propagating along 2D channel.

Member Enumeration Documentation

◆ anonymous enum

template<class ELEMENT >
anonymous enum
private

Enumeration of channel boundaries.

Enumerator
Inflow_boundary_id 
Upper_wall_boundary_id 
Outflow_boundary_id 
Bottom_wall_boundary_id 
First_drop_boundary_id 
Second_drop_boundary_id 
645  {
652  };
@ Inflow_boundary_id
Definition: adaptive_drop_in_channel.cc:646
@ Bottom_wall_boundary_id
Definition: adaptive_drop_in_channel.cc:649
@ First_drop_boundary_id
Definition: adaptive_drop_in_channel.cc:650
@ Outflow_boundary_id
Definition: adaptive_drop_in_channel.cc:648
@ Upper_wall_boundary_id
Definition: adaptive_drop_in_channel.cc:647
@ Second_drop_boundary_id
Definition: adaptive_drop_in_channel.cc:651

Constructor & Destructor Documentation

◆ DropInChannelProblem()

template<class ELEMENT >
DropInChannelProblem< ELEMENT >::DropInChannelProblem

Constructor.

663 {
664  // Output directory
666 
667  // Allocate the timestepper -- this constructs the Problem's
668  // time object with a sufficient amount of storage to store the
669  // previous timsteps.
670  this->add_time_stepper_pt(new BDF<2>);
671 
672 
673 
674  // Build the boundary segments for outer boundary, consisting of
675  //--------------------------------------------------------------
676  // four separate polylines
677  //------------------------
678  Vector<TriangleMeshCurveSection*> boundary_polyline_pt(4);
679 
680  // Each polyline only has two vertices -- provide storage for their
681  // coordinates
682  Vector<Vector<double> > vertex_coord(2);
683  for(unsigned i=0;i<2;i++)
684  {
685  vertex_coord[i].resize(2);
686  }
687 
688  // First polyline: Inflow
689  vertex_coord[0][0]=0.0;
690  vertex_coord[0][1]=0.0;
691  vertex_coord[1][0]=0.0;
692  vertex_coord[1][1]=1.0;
693 
694  // Build the 1st boundary polyline
695  boundary_polyline_pt[0] = new TriangleMeshPolyLine(vertex_coord,
697 
698  // Second boundary polyline: Upper wall
699  vertex_coord[0][0]=0.0;
700  vertex_coord[0][1]=1.0;
701  vertex_coord[1][0]=Problem_Parameter::Length;
702  vertex_coord[1][1]=1.0;
703 
704  // Build the 2nd boundary polyline
705  boundary_polyline_pt[1] = new TriangleMeshPolyLine(vertex_coord,
707 
708  // Third boundary polyline: Outflow
709  vertex_coord[0][0]=Problem_Parameter::Length;
710  vertex_coord[0][1]=1.0;
711  vertex_coord[1][0]=Problem_Parameter::Length;
712  vertex_coord[1][1]=0.0;
713 
714  // Build the 3rd boundary polyline
715  boundary_polyline_pt[2] = new TriangleMeshPolyLine(vertex_coord,
717 
718  // Fourth boundary polyline: Bottom wall
719  vertex_coord[0][0]=Problem_Parameter::Length;
720  vertex_coord[0][1]=0.0;
721  vertex_coord[1][0]=0.0;
722  vertex_coord[1][1]=0.0;
723 
724  // Build the 4th boundary polyline
725  boundary_polyline_pt[3] = new TriangleMeshPolyLine(vertex_coord,
727 
728  // Create the triangle mesh polygon for outer boundary
729  Outer_boundary_polyline_pt = new TriangleMeshPolygon(boundary_polyline_pt);
730 
731 
732  // Now define initial shape of drop(s) with polygon
733  //---------------------------------------------------
734 
735  // We have one drop
736  Drop_polygon_pt.resize(1);
737 
738  // Place it smack in the middle of the channel
739  double x_center = 0.5*Problem_Parameter::Length;
740  double y_center = 0.5;
741  Ellipse * drop_pt = new Ellipse(Problem_Parameter::Radius,
743 
744  // Intrinsic coordinate along GeomObject defining the drop
745  Vector<double> zeta(1);
746 
747  // Position vector to GeomObject defining the drop
748  Vector<double> coord(2);
749 
750  // Number of points defining drop
751  unsigned npoints = 16;
752  double unit_zeta = MathematicalConstants::Pi/double(npoints-1);
753 
754  // This drop is bounded by two distinct boundaries, each
755  // represented by its own polyline
756  Vector<TriangleMeshCurveSection*> drop_polyline_pt(2);
757 
758  // Vertex coordinates
759  Vector<Vector<double> > drop_vertex(npoints);
760 
761  // Create points on boundary
762  for(unsigned ipoint=0; ipoint<npoints;ipoint++)
763  {
764  // Resize the vector
765  drop_vertex[ipoint].resize(2);
766 
767  // Get the coordinates
768  zeta[0]=unit_zeta*double(ipoint);
769  drop_pt->position(zeta,coord);
770 
771  // Shift
772  drop_vertex[ipoint][0]=coord[0]+x_center;
773  drop_vertex[ipoint][1]=coord[1]+y_center;
774  }
775 
776  // Build the 1st drop polyline
777  drop_polyline_pt[0] = new TriangleMeshPolyLine(drop_vertex,
779 
780  // Second boundary of drop
781  for(unsigned ipoint=0; ipoint<npoints;ipoint++)
782  {
783  // Resize the vector
784  drop_vertex[ipoint].resize(2);
785 
786  // Get the coordinates
787  zeta[0]=(unit_zeta*double(ipoint))+MathematicalConstants::Pi;
788  drop_pt->position(zeta,coord);
789 
790  // Shift
791  drop_vertex[ipoint][0]=coord[0]+x_center;
792  drop_vertex[ipoint][1]=coord[1]+y_center;
793  }
794 
795  // Build the 2nd drop polyline
796  drop_polyline_pt[1] = new TriangleMeshPolyLine(drop_vertex,
798 
799  // Create closed polygon from two polylines
801  drop_polyline_pt);
802 
803  // Now build the mesh, based on the boundaries specified by
804  //---------------------------------------------------------
805  // polygons just created
806  //----------------------
807 
808  // Convert to "closed curve" objects
810  unsigned nb=Drop_polygon_pt.size();
811  Vector<TriangleMeshClosedCurve*> drop_closed_curve_pt(nb);
812  for (unsigned i=0;i<nb;i++)
813  {
814  drop_closed_curve_pt[i]=Drop_polygon_pt[i];
815  }
816 
817  // Target area for initial mesh
818  double uniform_element_area=0.2;
819 
820  // Define coordinates of a point inside the drop
821  Vector<double> drop_center(2);
822  drop_center[0]=x_center;
823  drop_center[1]=y_center;
824 
825  // Use the TriangleMeshParameters object for gathering all
826  // the necessary arguments for the TriangleMesh object
827  TriangleMeshParameters triangle_mesh_parameters(
828  outer_closed_curve_pt);
829 
830  // Define the holes on the boundary
831  triangle_mesh_parameters.internal_closed_curve_pt() =
832  drop_closed_curve_pt;
833 
834  // Define the maximum element area
835  triangle_mesh_parameters.element_area() =
836  uniform_element_area;
837 
838  // Define the region
839  triangle_mesh_parameters.add_region_coordinates(1, drop_center);
840 
841  // Create the mesh
842  Fluid_mesh_pt =
843  new RefineableSolidTriangleMesh<ELEMENT>(
844  triangle_mesh_parameters, this->time_stepper_pt());
845 
846  // Set error estimator for bulk mesh
848  Fluid_mesh_pt->spatial_error_estimator_pt()=error_estimator_pt;
849 
850  // Set targets for spatial adaptivity
851  Fluid_mesh_pt->max_permitted_error()=0.005;
852  Fluid_mesh_pt->min_permitted_error()=0.001;
853  Fluid_mesh_pt->max_element_size()=0.2;
854  Fluid_mesh_pt->min_element_size()=0.001;
855 
856  // Use coarser mesh during validation
858  {
859  Fluid_mesh_pt->min_element_size()=0.01;
860  }
861 
862  // Output boundary and mesh initial mesh for information
863  this->Fluid_mesh_pt->output_boundaries("boundaries.dat");
864  this->Fluid_mesh_pt->output("mesh.dat");
865 
866  // Provide initial value for drop pressure
869 
870  // Set boundary condition and complete the build of all elements
872 
873  // Construct the mesh of elements that impose the volume constraint
877 
878  // Construct the mesh of free surface elements
881 
882  // Combine meshes
883  //---------------
884 
885  // Add volume constraint sub mesh
887 
888  // Add Fluid_mesh_pt sub meshes
890 
891  // Add Free_surface sub meshes
892  this->add_sub_mesh(this->Free_surface_mesh_pt);
893 
894  // Build global mesh
895  this->build_global_mesh();
896 
897  // Setup equation numbering scheme
898  cout <<"Number of equations: " << this->assign_eqn_numbers() << std::endl;
899 
900 } // end_of_constructor
int i
Definition: BiCGSTAB_step_by_step.cpp:9
double Initial_value_for_drop_pressure
Backed up drop pressure between adaptations.
Definition: adaptive_drop_in_channel.cc:635
RefineableSolidTriangleMesh< ELEMENT > * Fluid_mesh_pt
Pointer to Fluid_mesh.
Definition: adaptive_drop_in_channel.cc:620
void complete_problem_setup()
Set boundary conditions and complete the build of all elements.
Definition: adaptive_drop_in_channel.cc:1025
void create_volume_constraint_elements()
Create elements that impose volume constraint on the drop.
Definition: adaptive_drop_in_channel.cc:956
Mesh * Free_surface_mesh_pt
Pointers to mesh of free surface elements.
Definition: adaptive_drop_in_channel.cc:614
void create_free_surface_elements()
Create free surface elements.
Definition: adaptive_drop_in_channel.cc:908
bool Use_volume_constraint
Bool to indicate if volume constraint is applied (only for steady run)
Definition: adaptive_drop_in_channel.cc:641
Vector< TriangleMeshPolygon * > Drop_polygon_pt
Vector storing pointer to the drop polygons.
Definition: adaptive_drop_in_channel.cc:623
Mesh * Volume_constraint_mesh_pt
Pointer to mesh containing elements that impose volume constraint.
Definition: adaptive_drop_in_channel.cc:617
TriangleMeshPolygon * Outer_boundary_polyline_pt
Triangle mesh polygon for outer boundary.
Definition: adaptive_drop_in_channel.cc:626
void set_directory(const std::string &directory)
Definition: oomph_utilities.cc:298
Definition: geom_objects.h:644
void position(const Vector< double > &zeta, Vector< double > &r) const
Position Vector at Lagrangian coordinate zeta.
Definition: geom_objects.h:745
Definition: mesh.h:67
void add_time_stepper_pt(TimeStepper *const &time_stepper_pt)
Definition: problem.cc:1545
unsigned add_sub_mesh(Mesh *const &mesh_pt)
Definition: problem.h:1330
void build_global_mesh()
Definition: problem.cc:1493
unsigned long assign_eqn_numbers(const bool &assign_local_eqn_numbers=true)
Definition: problem.cc:1989
TimeStepper *& time_stepper_pt()
Definition: problem.h:1524
Base class defining a closed curve for the Triangle mesh generation.
Definition: unstructured_two_d_mesh_geometry_base.h:1339
Definition: triangle_mesh.template.h:94
Class defining a polyline for use in Triangle Mesh generation.
Definition: unstructured_two_d_mesh_geometry_base.h:868
Class defining a closed polygon for the Triangle mesh generation.
Definition: unstructured_two_d_mesh_geometry_base.h:1451
Definition: oomph-lib/src/generic/Vector.h:58
Definition: error_estimator.h:266
int nb
Definition: level2_impl.h:286
double Pi
Definition: two_d_biharmonic.cc:235
EIGEN_STRONG_INLINE const Eigen::CwiseBinaryOp< Eigen::internal::scalar_zeta_op< typename DerivedX::Scalar >, const DerivedX, const DerivedQ > zeta(const Eigen::ArrayBase< DerivedX > &x, const Eigen::ArrayBase< DerivedQ > &q)
Definition: SpecialFunctionsArrayAPI.h:152
Z2ErrorEstimator * error_estimator_pt
Definition: MortaringCantileverCompareToNonMortaring.cpp:190
DocInfo Doc_info
Doc info.
Definition: refineable_two_layer_interface.cc:291
double Length
Length of the channel.
Definition: refineable_two_layer_interface.cc:323
double Radius
Initial radius of bubble.
Definition: adaptive_bubble_in_channel.cc:384
double Ca
Capillary number.
Definition: refineable_two_layer_interface.cc:314
bool command_line_flag_has_been_set(const std::string &flag)
Definition: oomph_utilities.cc:501

References oomph::TriangleMeshParameters::add_region_coordinates(), Problem_Parameter::Ca, oomph::CommandLineArgs::command_line_flag_has_been_set(), Problem_Parameter::Doc_info, oomph::TriangleMeshParameters::element_area(), MeshRefinement::error_estimator_pt, i, oomph::TriangleMeshParameters::internal_closed_curve_pt(), Problem_Parameter::Length, nb, BiharmonicTestFunctions2::Pi, oomph::Ellipse::position(), Problem_Parameter::Radius, oomph::DocInfo::set_directory(), and Eigen::zeta().

◆ ~DropInChannelProblem()

template<class ELEMENT >
DropInChannelProblem< ELEMENT >::~DropInChannelProblem ( )
inline

Destructor.

418  {
419  // Fluid timestepper
420  delete this->time_stepper_pt(0);
421 
422  // Kill data associated with outer boundary
424  for (unsigned j=0;j<n;j++)
425  {
427  }
429 
430  //Kill data associated with drops
431  unsigned n_drop = Drop_polygon_pt.size();
432  for(unsigned idrop=0;idrop<n_drop;idrop++)
433  {
434  unsigned n=Drop_polygon_pt[idrop]->npolyline();
435  for (unsigned j=0;j<n;j++)
436  {
437  delete Drop_polygon_pt[idrop]->polyline_pt(j);
438  }
439  delete Drop_polygon_pt[idrop];
440  }
441 
442  // Flush element of free surface elements
444  delete Free_surface_mesh_pt;
447 
448  // Delete error estimator
449  delete Fluid_mesh_pt->spatial_error_estimator_pt();
450 
451  // Delete fluid mesh
452  delete Fluid_mesh_pt;
453 
454  // Delete the global pressure drop data
455  delete Drop_pressure_data_pt;
456 
457  // Kill const eqn
459 
460  }
const unsigned n
Definition: CG3DPackingUnitTest.cpp:11
void delete_free_surface_elements()
Delete free surface elements.
Definition: adaptive_drop_in_channel.cc:561
Data * Drop_pressure_data_pt
Pointer to a global drop pressure datum.
Definition: adaptive_drop_in_channel.cc:629
void delete_volume_constraint_elements()
Delete volume constraint elements.
Definition: adaptive_drop_in_channel.cc:585
unsigned npolyline() const
Number of constituent polylines.
Definition: unstructured_two_d_mesh_geometry_base.h:1482
TriangleMeshPolyLine * polyline_pt(const unsigned &i) const
Pointer to i-th constituent polyline.
Definition: unstructured_two_d_mesh_geometry_base.h:1488
ConstitutiveLaw * Constitutive_law_pt
Constitutive law used to determine the mesh deformation.
Definition: jeffery_orbit.cc:82
std::ptrdiff_t j
Definition: tut_arithmetic_redux_minmax.cpp:2

References Problem_Parameter::Constitutive_law_pt, j, and n.

Member Function Documentation

◆ actions_after_adapt()

template<class ELEMENT >
void DropInChannelProblem< ELEMENT >::actions_after_adapt ( )
inlinevirtual

Actions after adapt: Rebuild the mesh of free surface elements.

Reimplemented from oomph::Problem.

478  {
479  // Create the elements that impose the displacement constraint
482 
483 
484  // Rebuild the Problem's global mesh from its various sub-meshes
485  this->rebuild_global_mesh();
486 
487  // Setup the problem again -- remember that fluid mesh has been
488  // completely rebuilt and its element's don't have any
489  // pointers to Re etc. yet
491 
492  }// end of actions_after_adapt
void rebuild_global_mesh()
Definition: problem.cc:1533

◆ actions_after_newton_solve()

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

Update the after solve (empty)

Reimplemented from oomph::Problem.

496 {}

◆ actions_before_adapt()

template<class ELEMENT >
void DropInChannelProblem< ELEMENT >::actions_before_adapt ( )
inlinevirtual

Actions before adapt: Wipe the mesh of free surface elements.

Reimplemented from oomph::Problem.

465  {
466  // Kill the elements and wipe surface mesh
469 
470  // Rebuild the Problem's global mesh from its various sub-meshes
471  this->rebuild_global_mesh();
472 
473  }// end of actions_before_adapt

◆ actions_before_newton_solve()

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

Update the problem specs before solve.

Reimplemented from oomph::Problem.

500  {
501  //Reset the Lagrangian coordinates of the nodes to be the current
502  //Eulerian coordinates -- this makes the current configuration
503  //stress free
504  Fluid_mesh_pt->set_lagrangian_nodal_coordinates();
505  }

◆ complete_problem_setup()

template<class ELEMENT >
void DropInChannelProblem< ELEMENT >::complete_problem_setup

Set boundary conditions and complete the build of all elements.

1026  {
1027  // Map to record if a given boundary is on a drop or not
1028  map<unsigned,bool> is_on_drop_bound;
1029 
1030  // Loop over the drops
1031  unsigned ndrop=Drop_polygon_pt.size();
1032  for(unsigned idrop=0;idrop<ndrop;idrop++)
1033  {
1034  // Get the vector all boundary IDs associated with the polylines that
1035  // make up the closed polygon
1036  Vector<unsigned> drop_bound_id=this->Drop_polygon_pt[idrop]->
1037  polygon_boundary_id();
1038 
1039  // Get the number of boundary
1040  unsigned nbound=drop_bound_id.size();
1041 
1042  // Fill in the map
1043  for(unsigned ibound=0;ibound<nbound;ibound++)
1044  {
1045  // This boundary...
1046  unsigned bound_id=drop_bound_id[ibound];
1047 
1048  // ...is on the drop
1049  is_on_drop_bound[bound_id]=true;
1050  }
1051  }
1052 
1053  // Re-set the boundary conditions for fluid problem: All nodes are
1054  // free by default -- just pin the ones that have Dirichlet conditions
1055  // here.
1056  unsigned nbound=Fluid_mesh_pt->nboundary();
1057  for(unsigned ibound=0;ibound<nbound;ibound++)
1058  {
1059  unsigned num_nod=Fluid_mesh_pt->nboundary_node(ibound);
1060  for (unsigned inod=0;inod<num_nod;inod++)
1061  {
1062  // Get node
1063  Node* nod_pt=Fluid_mesh_pt->boundary_node_pt(ibound,inod);
1064 
1065  //Pin both velocities on inflow (0) and side boundaries (1 and 3)
1066  if((ibound==0) || (ibound==1) || (ibound==3))
1067  {
1068  nod_pt->pin(0);
1069  nod_pt->pin(1);
1070  }
1071 
1072  //If it's the outflow pin only the vertical velocity
1073  if(ibound==2) {nod_pt->pin(1);}
1074 
1075  // Pin pseudo-solid positions apart from drop boundary which
1076  // we allow to move
1077  SolidNode* solid_node_pt = dynamic_cast<SolidNode*>(nod_pt);
1078  if(is_on_drop_bound[ibound])
1079  {
1080  solid_node_pt->unpin_position(0);
1081  solid_node_pt->unpin_position(1);
1082  }
1083  else
1084  {
1085  solid_node_pt->pin_position(0);
1086  solid_node_pt->pin_position(1);
1087  }
1088  }
1089  } // end loop over boundaries
1090 
1091  // Complete the build of all elements so they are fully functional
1092  // Remember that adaptation for triangle meshes involves a complete
1093  // regneration of the mesh (rather than splitting as in tree-based
1094  // meshes where such parameters can be passed down from the father
1095  // element!)
1096  unsigned n_element = Fluid_mesh_pt->nelement();
1097  for(unsigned e=0;e<n_element;e++)
1098  {
1099  // Upcast from GeneralisedElement to the present element
1100  ELEMENT* el_pt = dynamic_cast<ELEMENT*>(Fluid_mesh_pt->element_pt(e));
1101 
1102  // Set the Reynolds number
1103  el_pt->re_pt() = &Problem_Parameter::Re;
1104 
1105  // Set the Womersley number (same as Re since St=1)
1106  el_pt->re_st_pt() = &Problem_Parameter::Re;
1107 
1108  // Set the constitutive law for pseudo-elastic mesh deformation
1109  el_pt->constitutive_law_pt()=Problem_Parameter::Constitutive_law_pt;
1110  }
1111 
1112  //For the elements within the droplet (region 1),
1113  //set the viscosity ratio
1114  n_element = Fluid_mesh_pt->nregion_element(1);
1115  for(unsigned e=0;e<n_element;e++)
1116  {
1117  // Upcast from GeneralisedElement to the present element
1118  ELEMENT* el_pt =
1119  dynamic_cast<ELEMENT*>(Fluid_mesh_pt->region_element_pt(1,e));
1120 
1121  el_pt->viscosity_ratio_pt() = &Problem_Parameter::Viscosity_ratio;
1122  }
1123 
1124 
1125  // Re-apply boundary values on Dirichlet boundary conditions
1126  // (Boundary conditions are ignored when the solution is transferred
1127  // from the old to the new mesh by projection; this leads to a slight
1128  // change in the boundary values (which are, of course, never changed,
1129  // unlike the actual unknowns for which the projected values only
1130  // serve as an initial guess)
1131 
1132  // Set velocity and history values of velocity on walls
1133  nbound=this->Fluid_mesh_pt->nboundary();
1134  for(unsigned ibound=0;ibound<nbound;++ibound)
1135  {
1136  if ((ibound==Upper_wall_boundary_id)||
1137  (ibound==Bottom_wall_boundary_id)||
1138  (ibound==Outflow_boundary_id) ||
1139  (ibound==Inflow_boundary_id))
1140  {
1141  // Loop over nodes on this boundary
1142  unsigned num_nod=this->Fluid_mesh_pt->nboundary_node(ibound);
1143  for (unsigned inod=0;inod<num_nod;inod++)
1144  {
1145  // Get node
1146  Node* nod_pt=this->Fluid_mesh_pt->boundary_node_pt(ibound,inod);
1147 
1148  // Get number of previous (history) values
1149  unsigned n_prev=nod_pt->time_stepper_pt()->nprev_values();
1150 
1151  // Velocity is and was zero at all previous times
1152  for (unsigned t=0;t<=n_prev;t++)
1153  {
1154  if (ibound!=Inflow_boundary_id)
1155  {
1156  // Parallel outflow
1157  if (ibound!=Outflow_boundary_id)
1158  {
1159  nod_pt->set_value(t,0,0.0);
1160  }
1161  nod_pt->set_value(t,1,0.0);
1162  }
1163 
1164  // Nodes have always been there...
1165  nod_pt->x(t,0)=nod_pt->x(0,0);
1166  nod_pt->x(t,1)=nod_pt->x(0,1);
1167  }
1168  }
1169  }
1170  }
1171 
1172  // Re-assign prescribed inflow velocity at inlet
1173  unsigned num_nod=this->Fluid_mesh_pt->nboundary_node(Inflow_boundary_id);
1174  for (unsigned inod=0;inod<num_nod;inod++)
1175  {
1176  // Get node
1177  Node* nod_pt=this->Fluid_mesh_pt->boundary_node_pt(Inflow_boundary_id,
1178  inod);
1179  //Now set the boundary velocity
1180  double y = nod_pt->x(1);
1182  }
1183  }
Array< double, 1, 3 > e(1./3., 0.5, 2.)
void pin(const unsigned &i)
Pin the i-th stored variable.
Definition: nodes.h:385
TimeStepper *& time_stepper_pt()
Return the pointer to the timestepper.
Definition: nodes.h:238
void set_value(const unsigned &i, const double &value_)
Definition: nodes.h:271
Definition: nodes.h:906
double & x(const unsigned &i)
Return the i-th nodal coordinate.
Definition: nodes.h:1060
Definition: nodes.h:1686
void pin_position(const unsigned &i)
Pin the nodal position.
Definition: nodes.h:1816
void unpin_position(const unsigned &i)
Unpin the nodal position.
Definition: nodes.h:1829
virtual unsigned nprev_values() const =0
Number of previous values available: 0 for static, 1 for BDF<1>,...
Scalar * y
Definition: level1_cplx_impl.h:128
double Inflow_veloc_magnitude
Definition: adaptive_bubble_in_channel.cc:392
double Re
Reynolds number.
Definition: jeffery_orbit.cc:59
double Viscosity_ratio
Viscosity ratio of the droplet to the surrounding fluid.
Definition: adaptive_drop_in_channel.cc:375
t
Definition: plotPSD.py:36

References Problem_Parameter::Constitutive_law_pt, e(), Problem_Parameter::Inflow_veloc_magnitude, oomph::TimeStepper::nprev_values(), oomph::Data::pin(), oomph::SolidNode::pin_position(), Problem_Parameter::Re, oomph::Data::set_value(), plotPSD::t, oomph::Data::time_stepper_pt(), oomph::SolidNode::unpin_position(), Problem_Parameter::Viscosity_ratio, oomph::Node::x(), and y.

◆ compute_error_estimate()

template<class ELEMENT >
void DropInChannelProblem< ELEMENT >::compute_error_estimate ( double max_err,
double min_err 
)

Compute the error estimates and assign to elements for plotting.

Compute error estimates and assign to elements for plotting.

1270 {
1271  // Get error estimator
1272  ErrorEstimator* err_est_pt=Fluid_mesh_pt->spatial_error_estimator_pt();
1273 
1274  // Get/output error estimates
1275  unsigned nel=Fluid_mesh_pt->nelement();
1276  Vector<double> elemental_error(nel);
1277 
1278  // We need a dynamic cast, get_element_errors from the Fluid_mesh_pt
1279  // Dynamic cast is used because get_element_errors require a Mesh* ans
1280  // not a SolidMesh*
1281  Mesh* fluid_mesh_pt=dynamic_cast<Mesh*>(Fluid_mesh_pt);
1282  err_est_pt->get_element_errors(fluid_mesh_pt,
1283  elemental_error);
1284 
1285  // Set errors for post-processing and find extrema
1286  max_err=0.0;
1287  min_err=DBL_MAX;
1288  for (unsigned e=0;e<nel;e++)
1289  {
1290  dynamic_cast<MyCrouzeixRaviartElement*>(Fluid_mesh_pt->element_pt(e))->
1291  set_error(elemental_error[e]);
1292 
1293  max_err=std::max(max_err,elemental_error[e]);
1294  min_err=std::min(min_err,elemental_error[e]);
1295  }
1296 
1297 }
Base class for spatial error estimators.
Definition: error_estimator.h:40
void get_element_errors(Mesh *&mesh_pt, Vector< double > &elemental_error)
Definition: error_estimator.h:56
Overload CrouzeixRaviart element to modify output.
Definition: refineable_two_layer_interface.cc:54
#define min(a, b)
Definition: datatypes.h:22
#define max(a, b)
Definition: datatypes.h:23
void set_error(const double &error)
Set error value for post-processing.
Definition: overloaded_element_body.h:432

References e(), oomph::ErrorEstimator::get_element_errors(), max, min, and set_error().

◆ create_free_surface_elements()

template<class ELEMENT >
void DropInChannelProblem< ELEMENT >::create_free_surface_elements
private

Create free surface elements.

Create elements that impose the kinematic and dynamic bcs for the pseudo-solid fluid mesh

909 {
910 
911  //Loop over the free surface boundaries
912  unsigned nb=Fluid_mesh_pt->nboundary();
913  for(unsigned b=First_drop_boundary_id;b<nb;b++)
914  {
915  // Note: region is important
916  // How many bulk fluid elements are adjacent to boundary b in region 0?
917  unsigned n_element = Fluid_mesh_pt->nboundary_element_in_region(b,0);
918 
919  // Loop over the bulk fluid elements adjacent to boundary b?
920  for(unsigned e=0;e<n_element;e++)
921  {
922  // Get pointer to the bulk fluid element that is
923  // adjacent to boundary b in region 0
924  ELEMENT* bulk_elem_pt = dynamic_cast<ELEMENT*>(
925  Fluid_mesh_pt->boundary_element_in_region_pt(b,0,e));
926 
927  //Find the index of the face of element e along boundary b in region 0
928  int face_index = Fluid_mesh_pt->face_index_at_boundary_in_region(b,0,e);
929 
930  // Create new element
933  bulk_elem_pt,face_index);
934 
935  // Add it to the mesh
937 
938  //Add the appropriate boundary number
940 
941  //Specify the capillary number
942  el_pt->ca_pt() = &Problem_Parameter::Ca;
943  }
944  }
945 }
Scalar * b
Definition: benchVecAdd.cpp:17
Specialise the elastic update template class to concrete 1D case.
Definition: specific_node_update_interface_elements.h:1220
void set_boundary_number_in_bulk_mesh(const unsigned &b)
Set function for the boundary number in bulk mesh.
Definition: elements.h:4482
double *& ca_pt()
Pointer to the Capillary number.
Definition: interface_elements.h:492
void add_element_pt(GeneralisedElement *const &element_pt)
Add a (pointer to) an element to the mesh.
Definition: mesh.h:617

References b, Problem_Parameter::Ca, oomph::FluidInterfaceElement::ca_pt(), e(), nb, and oomph::FaceElement::set_boundary_number_in_bulk_mesh().

◆ create_volume_constraint_elements()

template<class ELEMENT >
void DropInChannelProblem< ELEMENT >::create_volume_constraint_elements
private

Create elements that impose volume constraint on the drop.

957 {
958 
959  // Do we need it?
960  if (!Use_volume_constraint) return;
961 
962  // Store pointer to element whose pressure we're trading/hi-jacking:
963  // Element 0 in region 1
964  Hijacked_element_pt= dynamic_cast<ELEMENT*>(
965  Fluid_mesh_pt->region_element_pt(1,0));
966 
967  // Set the global pressure data by hijacking one of the pressure values
968  // from inside the droplet
969  unsigned index_of_traded_pressure=0;
971  hijack_internal_value(0,index_of_traded_pressure);
972 
973  // Build volume constraint element -- pass traded pressure to it
977  index_of_traded_pressure);
978 
979  //Provide a reasonable initial guess for drop pressure (hydrostatics):
980  Drop_pressure_data_pt->set_value(index_of_traded_pressure,
982 
983  // Add volume constraint element to the mesh
985 
986  //Loop over the free surface boundaries
987  unsigned nb=Fluid_mesh_pt->nboundary();
988  for(unsigned b=First_drop_boundary_id;b<nb;b++)
989  {
990  // Note region is important
991  // How many bulk fluid elements are adjacent to boundary b in region 0?
992  unsigned n_element = Fluid_mesh_pt->nboundary_element_in_region(b,0);
993 
994  // Loop over the bulk fluid elements adjacent to boundary b?
995  for(unsigned e=0;e<n_element;e++)
996  {
997  // Get pointer to the bulk fluid element that is
998  // adjacent to boundary b in region 0?
999  ELEMENT* bulk_elem_pt = dynamic_cast<ELEMENT*>(
1000  Fluid_mesh_pt->boundary_element_in_region_pt(b,0,e));
1001 
1002  //Find the index of the face of element e along boundary b in region 0
1003  int face_index = Fluid_mesh_pt->face_index_at_boundary_in_region(b,0,e);
1004 
1005  // Create new element
1008  bulk_elem_pt,face_index);
1009 
1010  //Set the "master" volume constraint element
1012 
1013  // Add it to the mesh
1015  }
1016  }
1017 }
VolumeConstraintElement * Vol_constraint_el_pt
Pointer to element that imposes volume constraint for drop.
Definition: adaptive_drop_in_channel.cc:632
ELEMENT * Hijacked_element_pt
Pointer to hijacked element.
Definition: adaptive_drop_in_channel.cc:638
Definition: constrained_volume_elements.h:372
void set_volume_constraint_element(VolumeConstraintElement *const &vol_constraint_el_pt, const bool &check_nodal_data=true)
Definition: constrained_volume_elements.h:261
Definition: constrained_volume_elements.h:66
double Volume
Definition: adaptive_bubble_in_channel.cc:388

References b, e(), nb, oomph::VolumeConstraintBoundingElement::set_volume_constraint_element(), and Problem_Parameter::Volume.

◆ delete_free_surface_elements()

template<class ELEMENT >
void DropInChannelProblem< ELEMENT >::delete_free_surface_elements ( )
inlineprivate

Delete free surface elements.

562  {
563  // How many surface elements are in the surface mesh
564  unsigned n_element = Free_surface_mesh_pt->nelement();
565 
566  // Loop over the surface elements
567  for(unsigned e=0;e<n_element;e++)
568  {
569  // Kill surface element
571  }
572 
573 
574  // Wipe the mesh
576 
577 
578  } // end of delete_free_surface_elements
void flush_element_and_node_storage()
Definition: mesh.h:407
GeneralisedElement *& element_pt(const unsigned long &e)
Return pointer to element e.
Definition: mesh.h:448
unsigned long nelement() const
Return number of elements in the mesh.
Definition: mesh.h:590

References e().

◆ delete_volume_constraint_elements()

template<class ELEMENT >
void DropInChannelProblem< ELEMENT >::delete_volume_constraint_elements ( )
inlineprivate

Delete volume constraint elements.

586  {
587  if (Volume_constraint_mesh_pt==0) return;
588 
589  // Backup
590  if (Vol_constraint_el_pt!=0)
591  {
593  }
594 
595  // How many surface elements are in the surface mesh
596  unsigned n_element = Volume_constraint_mesh_pt->nelement();
597 
598  // Loop over all
599  unsigned first_el_to_be_killed=0;
600  for(unsigned e=first_el_to_be_killed;e<n_element;e++)
601  {
603  }
604 
605  // We've just killed the volume constraint element
607 
608  // Wipe the mesh
610 
611  } // end of delete_volume_constraint_elements
double p_traded()
Return the traded pressure value.
Definition: constrained_volume_elements.h:139

References e().

◆ doc_solution()

template<class ELEMENT >
void DropInChannelProblem< ELEMENT >::doc_solution ( const std::string &  comment = "")

Doc the solution.

1192 {
1193 
1194  oomph_info << "Docing step: " << Problem_Parameter::Doc_info.number()
1195  << std::endl;
1196 
1197  ofstream some_file;
1198  char filename[100];
1199  sprintf(filename,"%s/soln%i.dat",
1200  Problem_Parameter::Doc_info.directory().c_str(),
1201  Problem_Parameter::Doc_info.number());
1202 
1203  // Number of plot points
1204  unsigned npts;
1205  npts=5;
1206 
1207  // Compute errors and assign to each element for plotting
1208  double max_err;
1209  double min_err;
1210  compute_error_estimate(max_err,min_err);
1211 
1212  // Assemble square of L2 norm
1213  double square_of_l2_norm=0.0;
1214  unsigned nel=Fluid_mesh_pt->nelement();
1215  for (unsigned e=0;e<nel;e++)
1216  {
1218  dynamic_cast<ELEMENT*>(this->Fluid_mesh_pt->element_pt(e))->
1220  }
1221  Problem_Parameter::Norm_file << sqrt(square_of_l2_norm) << std::endl;
1222 
1223 
1224  some_file.open(filename);
1225  some_file << dynamic_cast<ELEMENT*>(this->Fluid_mesh_pt->element_pt(0))
1226  ->variable_identifier();
1227  this->Fluid_mesh_pt->output(some_file,npts);
1228  some_file << "TEXT X = 25, Y = 78, CS=FRAME T = \"Global Step "
1229  << Problem_Parameter::Doc_info.number() << " "
1230  << comment << "\"\n";
1231  some_file.close();
1232 
1233  // Output boundaries
1234  sprintf(filename,"%s/boundaries%i.dat",
1235  Problem_Parameter::Doc_info.directory().c_str(),
1236  Problem_Parameter::Doc_info.number());
1237  some_file.open(filename);
1238  this->Fluid_mesh_pt->output_boundaries(some_file);
1239  some_file.close();
1240 
1241 
1242  // Get max/min area
1243  double max_area;
1244  double min_area;
1245  Fluid_mesh_pt->max_and_min_element_size(max_area, min_area);
1246 
1247  // Write trace file
1249  << this->time_pt()->time() << " "
1250  << Fluid_mesh_pt->nelement() << " "
1251  << max_area << " "
1252  << min_area << " "
1253  << max_err << " "
1254  << min_err << " "
1255  << sqrt(square_of_l2_norm) << " "
1256  << std::endl;
1257 
1258  // Increment the doc_info number
1260 
1261 
1262 }
AnnoyingScalar sqrt(const AnnoyingScalar &x)
Definition: AnnoyingScalar.h:134
void compute_error_estimate(double &max_err, double &min_err)
Compute the error estimates and assign to elements for plotting.
Definition: adaptive_drop_in_channel.cc:1268
unsigned & number()
Number used (e.g.) for labeling output files.
Definition: oomph_utilities.h:554
Time *& time_pt()
Return a pointer to the global time object.
Definition: problem.h:1504
double & time()
Return the current value of the continuous time.
Definition: timesteppers.h:123
string filename
Definition: MergeRestartFiles.py:39
ofstream Norm_file
Definition: refineable_two_layer_interface.cc:341
ofstream Trace_file
Trace file.
Definition: refineable_two_layer_interface.cc:335
OomphInfo oomph_info
Definition: oomph_definitions.cc:319
double square_of_l2_norm()
Get square of L2 norm of velocity.
Definition: overloaded_element_body.h:1031

References Problem_Parameter::Doc_info, e(), MergeRestartFiles::filename, Problem_Parameter::Norm_file, oomph::DocInfo::number(), oomph::oomph_info, sqrt(), square_of_l2_norm(), and Problem_Parameter::Trace_file.

◆ remove_volume_constraint()

template<class ELEMENT >
void DropInChannelProblem< ELEMENT >::remove_volume_constraint ( )
inline

Change the boundary conditions to remove the volume constraint.

520  {
521  // Ignore all volume constraint stuff from here onwards
522  Use_volume_constraint=false;
523 
524  //Unhijack the data in the internal element
525  Hijacked_element_pt->unhijack_all_data();
526 
527  //Delete the volume constraint elements
529 
530  // Internal pressure is gone -- null it out
532 
533  // Kill the mesh too
536 
537  //Remove the sub meshes
538  this->flush_sub_meshes();
539 
540  // Add Fluid_mesh_pt sub meshes
542 
543  // Add Free_surface sub meshes
544  this->add_sub_mesh(this->Free_surface_mesh_pt);
545 
546  //Rebuild the global mesh
547  this->rebuild_global_mesh();
548 
549  //Renumber the equations
550  std::cout << "Removed volume constraint to obtain "
551  << this->assign_eqn_numbers() << " new equation numbers\n";
552  }
void flush_sub_meshes()
Definition: problem.h:1339

Member Data Documentation

◆ Drop_polygon_pt

template<class ELEMENT >
Vector<TriangleMeshPolygon*> DropInChannelProblem< ELEMENT >::Drop_polygon_pt
private

Vector storing pointer to the drop polygons.

◆ Drop_pressure_data_pt

template<class ELEMENT >
Data* DropInChannelProblem< ELEMENT >::Drop_pressure_data_pt
private

Pointer to a global drop pressure datum.

◆ Fluid_mesh_pt

template<class ELEMENT >
RefineableSolidTriangleMesh<ELEMENT>* DropInChannelProblem< ELEMENT >::Fluid_mesh_pt
private

Pointer to Fluid_mesh.

◆ Free_surface_mesh_pt

template<class ELEMENT >
Mesh* DropInChannelProblem< ELEMENT >::Free_surface_mesh_pt
private

Pointers to mesh of free surface elements.

◆ Hijacked_element_pt

template<class ELEMENT >
ELEMENT* DropInChannelProblem< ELEMENT >::Hijacked_element_pt
private

Pointer to hijacked element.

◆ Initial_value_for_drop_pressure

template<class ELEMENT >
double DropInChannelProblem< ELEMENT >::Initial_value_for_drop_pressure
private

Backed up drop pressure between adaptations.

◆ Outer_boundary_polyline_pt

template<class ELEMENT >
TriangleMeshPolygon* DropInChannelProblem< ELEMENT >::Outer_boundary_polyline_pt
private

Triangle mesh polygon for outer boundary.

◆ Use_volume_constraint

template<class ELEMENT >
bool DropInChannelProblem< ELEMENT >::Use_volume_constraint
private

Bool to indicate if volume constraint is applied (only for steady run)

◆ Vol_constraint_el_pt

template<class ELEMENT >
VolumeConstraintElement* DropInChannelProblem< ELEMENT >::Vol_constraint_el_pt
private

Pointer to element that imposes volume constraint for drop.

◆ Volume_constraint_mesh_pt

template<class ELEMENT >
Mesh* DropInChannelProblem< ELEMENT >::Volume_constraint_mesh_pt
private

Pointer to mesh containing elements that impose volume constraint.


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