29 #ifndef OOMPH_PROBLEM_CLASS_HEADER
30 #define OOMPH_PROBLEM_CLASS_HEADER
35 #include <oomph-lib-config.h>
67 class RefineableElement;
70 class PRefineableElement;
76 class GeneralisedElement;
85 class AssemblyHandler;
91 class InvertedElementError;
159 template<
unsigned NNODE_1D>
254 oomph_info <<
"Called empty hook fct with i=" <<
i << std::endl;
269 std::map<double*, bool>::iterator it =
281 double*
const& parameter_pt)
329 const double& epsilon,
330 const unsigned& max_adapt,
331 const unsigned& suppress_resolve_after_spatial_adapt,
333 const bool& shift =
true);
369 bool compressed_row_flag);
382 bool compressed_row_flag);
395 bool compressed_row_flag);
408 bool compressed_row_flag);
421 bool compressed_row_flag);
438 double*
const& parameter_pt);
446 unsigned& n_unrefined,
447 const unsigned& bifurcation_type,
448 const bool& actually_adapt =
true);
472 const unsigned& el_lo,
473 const unsigned& el_hi,
492 const bool& do_external_halos);
499 bool& actually_removed_some_data);
686 bool compressed_row_flag);
921 unsigned& max_refinement_level_overall);
943 const unsigned& max_refinement_level_overall,
945 flat_packed_refinement_info_for_root);
963 const unsigned& max_refinement_level_overall,
1014 unsigned n = last_el_for_assembly.size();
1015 for (
unsigned i = 0;
i < n;
i++)
1030 unsigned n = last_el_for_assembly.size();
1031 for (
unsigned i = 0;
i < n;
i++)
1033 last_el_for_assembly[
i] -= 1;
1134 #ifdef OOMPH_HAS_MPI
1154 double*
const& parameter_pt)
1220 std::ostringstream warn_message;
1222 <<
"Warning: We're using the default (empty) set_initial_condition().\n"
1223 <<
"If the initial conditions isn't re-assigned after a mesh adaption "
1225 <<
"the initial conditions will represent the interpolation of the \n"
1226 <<
"initial conditions that were assigned on the original coarse "
1229 "Problem::set_initial_condition()",
1230 OOMPH_EXCEPTION_LOCATION);
1253 "The global_temporal_error_norm function will be problem-specific:\n";
1254 error_message +=
"Please write your own in your Problem class";
1257 error_message, OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION);
1375 #ifdef OOMPH_HAS_MPI
1387 std::string error_message =
"You are likely to have reached this error "
1388 "from Problem::load_balance()\n";
1390 "The build_mesh() function is problem-specific and must be supplied:\n";
1392 "Please write your own in your Problem class, ensuring that\n";
1393 error_message +=
"any (sub)meshes are built before calling "
1394 "Problem::build_global_mesh().\n";
1396 error_message, OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION);
1409 bool report_stats =
false;
1416 doc_info, report_stats, input_target_domain_for_local_non_halo_element);
1434 doc_info, report_stats, input_target_domain_for_local_non_halo_element);
1449 doc_info, report_stats, input_target_domain_for_local_non_halo_element);
1461 const bool& report_stats,
1481 const unsigned& max_level_overall);
1539 double time()
const;
1549 OOMPH_CURRENT_FUNCTION,
1550 OOMPH_EXCEPTION_LOCATION);
1562 OOMPH_CURRENT_FUNCTION,
1563 OOMPH_EXCEPTION_LOCATION);
1584 const bool& preserve_existing_data =
false);
1745 const bool& assign_local_eqn_numbers =
true);
1796 #ifdef OOMPH_HAS_MPI
1803 std::ostringstream error_stream;
1805 <<
"Direct access to global dofs in a distributed problem is only\n"
1806 <<
"possible if the function setup_dof_halo_scheme() has been "
1808 <<
"You can access local dofs via Problem::dof(i).\n";
1810 OOMPH_CURRENT_FUNCTION,
1811 OOMPH_EXCEPTION_LOCATION);
1817 const unsigned first_row_local =
1819 const unsigned n_row_local =
1823 if ((
i >= first_row_local) && (
i < first_row_local + n_row_local))
1825 return this->Dof_pt[
i - first_row_local];
1838 return this->Dof_pt[
i];
1849 double dof(
const unsigned&
i)
const
1911 std::ostringstream error_stream;
1913 <<
"You must overload this function in your problem class to specify\n"
1914 <<
"which matrices should be summed to give the final Jacobian.";
1916 error_stream.str(), OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION);
1958 Vector<std::complex<double>>& alpha,
1962 const bool& steady =
true);
1975 Vector<std::complex<double>>& eigenvalue,
1978 const bool& steady =
true);
1990 Vector<std::complex<double>>& eigenvalue,
1991 const bool& make_timesteppers_steady =
true)
2000 make_timesteppers_steady);
2015 Vector<std::complex<double>>& alpha,
2017 const bool& steady =
true)
2023 n_eval, alpha, beta, eigenvector_real, eigenvector_imag, steady);
2030 Vector<std::complex<double>>& eigenvalue,
2033 const bool& steady =
true);
2040 Vector<std::complex<double>>& eigenvalue,
2041 const bool& steady =
true)
2047 n_eval, eigenvalue, eigenvector_real, eigenvector_imag, steady);
2055 const double& shift = 0.0);
2123 const double& half_residual_squared_old,
2126 double& half_residual_squared,
2127 const double& stpmax);
2166 virtual void read(std::ifstream& restart_file,
bool& unsteady_restart);
2171 virtual void read(std::ifstream& restart_file)
2173 bool unsteady_restart;
2174 read(restart_file, unsteady_restart);
2179 virtual void dump(std::ofstream& dump_file)
const;
2185 std::ofstream dump_stream(dump_file_name.c_str());
2187 if (!dump_stream.is_open())
2189 std::string err =
"Couldn't open file " + dump_file_name;
2191 err, OOMPH_EXCEPTION_LOCATION, OOMPH_CURRENT_FUNCTION);
2197 #ifdef OOMPH_HAS_MPI
2212 void synchronise_dofs(
const bool& do_halos,
const bool& do_external_halos);
2232 const bool& report_stats =
false);
2237 const bool& report_stats =
false);
2242 const bool& report_stats =
false);
2255 const bool& report_stats =
false);
2264 const bool& report_stats);
2392 const bool& block_solve =
true);
2398 const bool& block_solve =
true);
2406 const bool& block_solve =
true);
2427 const bool& block_solve =
true);
2438 const bool& block_solve =
true);
2451 const double& omega,
2454 const bool& block_solve =
true);
2475 const unsigned& max_adapt);
2493 const unsigned& max_adapt = 0);
2505 const unsigned& data_index,
2507 const unsigned& max_adapt = 0);
2559 const unsigned& max_adapt,
2561 const bool& shift =
true);
2577 const double& epsilon,
2578 const unsigned& max_adapt,
2580 const bool& shift =
true)
2584 unsigned suppress_resolve_after_spatial_adapt_flag = 0;
2589 suppress_resolve_after_spatial_adapt_flag,
2611 const double& epsilon,
2612 const unsigned& max_adapt,
2613 const unsigned& suppress_resolve_after_spatial_adapt_flag,
2615 const bool& shift =
true)
2622 suppress_resolve_after_spatial_adapt_flag,
2637 const double& epsilon);
2649 const double& epsilon,
2650 const bool& shift_values);
2722 unsigned nmesh = std::max(
unsigned(1),
nsub_mesh());
2735 unsigned nmesh = std::max(
unsigned(1),
nsub_mesh());
2789 "Problem::p_refine_uniformly_and_prune()",
2790 OOMPH_EXCEPTION_LOCATION);
2804 "Problem::p_refine_uniformly_and_prune()",
2805 OOMPH_EXCEPTION_LOCATION);
2815 unsigned nmesh = std::max(
unsigned(1),
nsub_mesh());
2829 "Problem::p_refine_uniformly_and_prune()",
2830 OOMPH_EXCEPTION_LOCATION);
2832 unsigned nmesh = std::max(
unsigned(1),
nsub_mesh());
2880 const unsigned& i_mesh,
2915 const unsigned& i_mesh,
2955 void adapt(
unsigned& n_refined,
unsigned& n_unrefined);
2966 unsigned n_refined, n_unrefined;
2967 adapt(n_refined, n_unrefined);
2977 void p_adapt(
unsigned& n_refined,
unsigned& n_unrefined);
2988 unsigned n_refined, n_unrefined;
2989 p_adapt(n_refined, n_unrefined);
3002 unsigned& n_refined,
3003 unsigned& n_unrefined,
3017 unsigned n_refined, n_unrefined;
3072 OOMPH_CURRENT_FUNCTION,
3073 OOMPH_EXCEPTION_LOCATION),
3083 OOMPH_CURRENT_FUNCTION,
3084 OOMPH_EXCEPTION_LOCATION),
3094 OOMPH_CURRENT_FUNCTION,
3095 OOMPH_EXCEPTION_LOCATION),
A class that is used to define the functions used to assemble the elemental contributions to the resi...
A custom linear solver class that is used to solve a block-factorised version of the Fold bifurcation...
A custom linear solver class that is used to solve a block-factorised version of the PitchFork bifurc...
A custom linear solver class that is used to solve a block-factorised version of the Hopf bifurcation...
A custom linear solver class that is used to solve a block-factorised version of the PitchFork bifurc...
//////////////////////////////////////////////////////////////// ////////////////////////////////////...
A class for compressed row matrices. This is a distributable object.
GeneralisedTimestepper used to store the arclength derivatives and pervious solutions required in con...
A class that represents a collection of data; each Data object may contain many different individual ...
Class of matrices containing doubles, and stored as a DenseMatrix<double>, but with solving functiona...
Information for documentation of results: Directory and file number to enable output in the form RESL...
void disable_doc()
Disable documentation.
A class that stores the halo/haloed entries required when using a DoubleVectorWithHaloEntries....
LinearAlgebraDistribution *& distribution_pt()
Return the pointer to the distirbution used to setup the halo information.
unsigned local_index(const unsigned &global_eqn)
Return the local index associated with the global equation.
===================================================================== An extension of DoubleVector th...
A vector in the mathematical sense, initially developed for linear algebra type applications....
Base class for all EigenProblem solves. This simply defines standard interfaces so that different sol...
Class for objects than can be advanced in time by an Explicit Timestepper. WARNING: For explicit time...
A Base class for explicit timesteppers.
A class that is used to assemble the augmented system that defines a fold (saddle-node) or limit poin...
A class that is used to assemble the augmented system that defines a Hopf bifurcation....
Describes the distribution of a distributable linear algebra type object. Typically this is a contain...
unsigned first_row() const
access function for the first row on this processor. If not distributed then this is just zero.
unsigned nrow() const
access function to the number of global rows.
unsigned nrow_local() const
access function for the num of local rows on this processor. If no MPI then Nrow is returned.
Base class for all linear solvers. This merely defines standard interfaces for linear solvers,...
A class to handle errors in the Newton solver.
bool linear_solver_error()
Access function to the error in the linear solver.
NewtonSolverError(unsigned Passed_iterations, double Passed_maxres)
Constructor that passes number of iterations and residuals.
unsigned Iterations
Max. # of iterations performed when the Newton solver died.
unsigned iterations()
Access function to Max. # of iterations performed when the Newton solver died.
double Maxres
Max. residual when Newton solver died.
bool Linear_solver_error
Error in the linear solver.
NewtonSolverError()
Default constructor, does nothing.
double maxres()
Access function to Max. residual when Newton solver died.
NewtonSolverError(const bool &Passed_linear_failure)
Constructor that passes a failure of the linear solver.
An oomph-lib wrapper to the MPI_Comm communicator object. Just contains an MPI_Comm object (which is ...
std::ostream *& stream_pt()
Access function for the stream pointer.
An OomphLibError object which should be thrown when an run-time error is encountered....
An OomphLibWarning object which should be created as a temporary object to issue a warning....
A class that is used to assemble and solve the augmented system of equations associated with calculat...
A class that is used to assemble the augmented system that defines a pitchfork (symmetry-breaking) bi...
////////////////////////////////////////////////////////////////// //////////////////////////////////...
EigenSolver *& eigen_solver_pt()
Return a pointer to the eigen solver object.
virtual void actions_after_implicit_timestep()
Actions that should be performed after each implicit time step. This is needed when one wants to solv...
bool Always_take_one_newton_step
Boolean to indicate whether a Newton step should be taken even if the initial residuals are below the...
bool Jacobian_reuse_is_enabled
Is re-use of Jacobian in Newton iteration enabled? Default: false.
void solve_eigenproblem(const unsigned &n_eval, Vector< std::complex< double >> &eigenvalue, const bool &make_timesteppers_steady=true)
Solve an eigenproblem as assembled by the Problem's constituent elements but only return the eigenval...
double *& dof_pt(const unsigned &i)
Pointer to i-th dof in the problem.
virtual void actions_after_newton_solve()
Any actions that are to be performed after a complete Newton solve, e.g. post processing....
void parallel_sparse_assemble(const LinearAlgebraDistribution *const &dist_pt, Vector< int * > &column_or_row_index, Vector< int * > &row_or_column_start, Vector< double * > &value, Vector< unsigned > &nnz, Vector< double * > &residuals)
Helper method to assemble CRDoubleMatrices from distributed on multiple processors.
void describe_dofs(std::ostream &out= *(oomph_info.stream_pt())) const
Function to describe the dofs in terms of the global equation number, i.e. what type of value (nodal ...
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)
Protected helper function that is used to assemble the Jacobian matrix in the case when the storage i...
bool Store_local_dof_pt_in_elements
Boolean to indicate whether local dof pointers should be stored in the elements.
void clear_elemental_assembly_time()
Clear storage of most-recent elemental assembly times (used for load balancing). Next load balancing ...
virtual void actions_before_newton_step()
Any actions that are to be performed before each individual Newton step. Most likely to be used for d...
void remove_duplicate_data(Mesh *const &mesh_pt, bool &actually_removed_some_data)
Private helper function to remove repeated data in external haloed elements in specified mesh....
bool Bifurcation_detection
Boolean to control bifurcation detection via determinant of Jacobian.
void p_refine_uniformly(const Vector< unsigned > &nrefine_for_mesh, DocInfo &doc_info)
p-refine p-refineable sub-meshes, each as many times as specified in the vector and rebuild problem; ...
void adapt()
Adapt problem: Perform mesh adaptation for (all) refineable (sub)mesh(es), based on their own error e...
virtual void actions_before_newton_solve()
Any actions that are to be performed before a complete Newton solve (e.g. adjust boundary conditions)...
Vector< unsigned > First_el_for_assembly
First element to be assembled by given processor for non-distributed problem (only kept up to date wh...
unsigned add_sub_mesh(Mesh *const &mesh_pt)
Add a submesh to the problem and return its number, i, by which it can be accessed via mesh_pt(i).
virtual void get_jacobian(DoubleVector &residuals, SumOfMatrices &jacobian)
Dummy virtual function that must be overloaded by the problem to specify which matrices should be sum...
unsigned long assign_eqn_numbers(const bool &assign_local_eqn_numbers=true)
Assign all equation numbers for problem: Deals with global data (= data that isn't attached to any el...
void refine_uniformly_aux(const Vector< unsigned > &nrefine_for_mesh, DocInfo &doc_info, const bool &prune)
Helper function to do compund refinement of (all) refineable (sub)mesh(es) uniformly as many times as...
void build_global_mesh()
Build the global mesh by combining the all the submeshes. Note: The nodes boundary information refers...
unsigned unrefine_uniformly()
Refine (all) refineable (sub)mesh(es) uniformly and rebuild problem. Return 0 for success,...
void assign_initial_values_impulsive()
Initialise data and nodal positions to simulate impulsive start from initial configuration/solution.
double & dof_current(const unsigned &i)
Access function to the current value of the i-th (local) dof at the start of a continuation step.
LinearSolver * mass_matrix_solver_for_explicit_timestepper_pt() const
Return a pointer to the linear solver object used for explicit time stepping (const version).
void unset_default_partition_in_load_balance()
Do not use the default partition in the load balance.
double Theta_squared
Value of the scaling parameter required so that the parameter occupies the desired proportion of the ...
void refine_uniformly(const Vector< unsigned > &nrefine_for_mesh, DocInfo &doc_info)
Refine refineable sub-meshes, each as many times as specified in the vector and rebuild problem; doc ...
EigenSolver *const & eigen_solver_pt() const
Return a pointer to the eigen solver object (const version)
LinearAlgebraDistribution *const & dof_distribution_pt() const
Return the pointer to the dof distribution (read-only)
virtual void get_eigenproblem_matrices(CRDoubleMatrix &mass_matrix, CRDoubleMatrix &main_matrix, const double &shift=0.0)
Get the matrices required by a eigensolver. If the shift parameter is non-zero the second matrix will...
void refine_uniformly(const Vector< unsigned > &nrefine_for_mesh)
Refine refineable sub-meshes, each as many times as specified in the vector and rebuild problem.
Vector< double > Dof_derivative
Storage for the derivative of the problem variables wrt arc-length.
virtual void sparse_assemble_row_or_column_compressed_with_lists(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)
Private helper function that is used to assemble the Jacobian matrix in the case when the storage is ...
virtual void sparse_assemble_row_or_column_compressed_with_two_arrays(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)
Private helper function that is used to assemble the Jacobian matrix in the case when the storage is ...
void synchronise_dofs(const bool &do_halos, const bool &do_external_halos)
Synchronise the degrees of freedom by overwriting the haloed values with their non-halo counterparts ...
virtual void actions_before_distribute()
Actions to be performed before a (mesh) distribution.
void load_balance(const bool &report_stats)
Balance the load of a (possibly non-uniformly refined) problem that has already been distributed,...
Distributed_problem_matrix_distribution Dist_problem_matrix_distribution
The distributed matrix distribution method 1 - Automatic - the Problem distribution is employed,...
DoubleVectorHaloScheme * Halo_scheme_pt
Pointer to the halo scheme for any global vectors that have the Dof_distribution.
virtual void actions_after_change_in_global_parameter(double *const ¶meter_pt)
Actions that are to be performed when the global parameter addressed by parameter_pt has been changed...
double dof(const unsigned &i) const
i-th dof in the problem (const version)
void p_refine_selected_elements(const Vector< unsigned > &elements_to_be_refined)
p-refine (one and only!) mesh by refining the elements identified by their numbers relative to the pr...
void setup_dof_halo_scheme()
Function that is used to setup the halo scheme.
unsigned long ndof() const
Return the number of dofs.
void p_unrefine_uniformly(DocInfo &doc_info)
p-unrefine (all) p-refineable (sub)mesh(es) uniformly and rebuild problem.
virtual void build_mesh()
Function to build the Problem's base mesh; this must be supplied by the user if they wish to use the ...
void(* SpatialErrorEstimatorWithDocFctPt)(Mesh *&mesh_pt, Vector< double > &elemental_error, DocInfo &doc_info)
Function pointer for spatial error estimator with doc.
static bool Suppress_warning_about_actions_before_read_unstructured_meshes
Flag to allow suppression of warning messages re reading in unstructured meshes during restart.
void get_first_and_last_element_for_assembly(Vector< unsigned > &first_el_for_assembly, Vector< unsigned > &last_el_for_assembly) const
Get first and last elements for parallel assembly of non-distributed problem.
void solve_adjoint_eigenproblem(const unsigned &n_eval, Vector< std::complex< double >> &eigenvalue, Vector< DoubleVector > &eigenvector_real, Vector< DoubleVector > &eigenvector_imag, const bool &steady=true)
Solve an adjoint eigenvalue problem using the same procedure as solve_eigenproblem....
void set_default_partition_in_load_balance()
Set the use of the default partition in the load balance.
bool Use_default_partition_in_load_balance
Flag to use "default partition" during load balance. Should only be set to true when run in validatio...
void bifurcation_adapt_helper(unsigned &n_refined, unsigned &n_unrefined, const unsigned &bifurcation_type, const bool &actually_adapt=true)
A function that is used to adapt a bifurcation-tracking problem, which requires separate interpolatio...
double & max_residuals()
Access function to max residuals in Newton iterations before giving up.
Vector< double > Elemental_assembly_time
Storage for assembly times (used for load balancing)
bool Jacobian_has_been_computed
Has a Jacobian been computed (and can therefore be re-used if required)? Default: false.
bool Bypass_increase_in_dof_check_during_pruning
Boolean to bypass check of increase in dofs during pruning.
virtual void get_dvaluesdt(DoubleVector &f)
Get the time derivative of all values (using get_inverse_mass_matrix_times_residuals(....
void enable_store_local_dof_pt_in_elements()
Insist that local dof pointers are set up in each element when equation numbering takes place.
void enable_info_in_newton_solve()
Enable the output of information when in the newton solver (Default)
friend class BlockFoldLinearSolver
Problem(const Problem &dummy)=delete
Broken copy constructor.
void p_refine_uniformly(const unsigned &i_mesh)
Do uniform p-refinement for submesh i_mesh without documentation.
void initialise_dt(const double &dt)
Set all timesteps to the same value, dt, and assign weights for all timesteppers in the problem.
void add_to_dofs(const double &lambda, const DoubleVector &increment_dofs)
Add lambda x incremenet_dofs[l] to the l-th dof.
void p_refine_uniformly(DocInfo &doc_info)
p-refine (all) p-refineable (sub)mesh(es) uniformly and rebuild problem; doc refinement process.
double Timestep_reduction_factor_after_nonconvergence
What it says: If temporally adaptive Newton solver fails to to converge, reduce timestep by this fact...
Mesh *const & mesh_pt(const unsigned &imesh) const
Return a pointer to the i-th submesh (const version)
Vector< double > Dof_current
Storage for the present values of the variables.
virtual void debug_hook_fct(const unsigned &i)
Hook for debugging. Can be overloaded in driver code; argument allows identification of where we're c...
double Desired_proportion_of_arc_length
Proportion of the arc-length to taken by the parameter.
void adapt_based_on_error_estimates(Vector< Vector< double >> &elemental_error)
Adapt problem: Perform mesh adaptation for (all) refineable (sub)mesh(es), based on the error estimat...
AssemblyHandler *& assembly_handler_pt()
Return a pointer to the assembly handler object.
virtual void actions_after_implicit_timestep_and_error_estimation()
Actions that should be performed after each implicit time step. This is needed if your actions_after_...
void refine_uniformly_and_prune(const Vector< unsigned > &nrefine_for_mesh)
Refine refineable sub-meshes, each as many times as specified in the vector and rebuild problem....
void disable_mass_matrix_reuse()
Turn off recyling of the mass matrix in explicit timestepping schemes.
EigenSolver * Default_eigen_solver_pt
Pointer to the default eigensolver.
Mesh *const & mesh_pt() const
Return a pointer to the global mesh (const version)
void p_refine_uniformly(const Vector< unsigned > &nrefine_for_mesh)
p-refine p-refineable sub-meshes, each as many times as specified in the vector and rebuild problem
unsigned Nnewton_iter_taken
Actual number of Newton iterations taken during the most recent iteration.
double * bifurcation_parameter_pt() const
Return pointer to the parameter that is used in the bifurcation detection. If we are not tracking a b...
void copy_haloed_eqn_numbers_helper(const bool &do_halos, const bool &do_external_halos)
A private helper function to copy the haloed equation numbers into the halo equation numbers,...
virtual void sparse_assemble_row_or_column_compressed_with_vectors_of_pairs(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)
Private helper function that is used to assemble the Jacobian matrix in the case when the storage is ...
void send_data_to_be_sent_during_load_balancing(Vector< int > &send_n, Vector< double > &send_data, Vector< int > &send_displacement)
Load balance helper routine: Send data to other processors during load balancing.
bool Time_adaptive_newton_crash_on_solve_fail
Bool to specify what to do if a Newton solve fails within a time adaptive solve. Default (false) is t...
Problem()
Constructor: Allocate space for one time stepper and set all pointers to NULL and set defaults for al...
const TimeStepper * time_stepper_pt() const
Access function for the pointer to the first (presumably only) timestepper.
void enable_problem_distributed()
Enable problem distributed.
bool is_dparameter_calculated_analytically(double *const ¶meter_pt)
Function to determine whether the parameter derivatives are calculated analytically.
void flush_sub_meshes()
Flush the problem's collection of sub-meshes. Must be followed by call to rebuild_global_mesh().
bool & use_predictor_values_as_initial_guess()
unsigned Sparse_assembly_method
Flag to determine which sparse assembly method to use By default we use assembly by vectors of pairs.
void set_analytic_hessian_products()
Function to turn on analytic calculation of the parameter derivatives in continuation and bifurcation...
Time *& time_pt()
Return a pointer to the global time object.
void calculate_predictions()
Calculate predictions.
Vector< unsigned > Last_el_plus_one_for_assembly
Last element (plus one) to be assembled by given processor for non-distributed problem (only kept up ...
AssemblyHandler *const & assembly_handler_pt() const
Return a pointer to the assembly handler object (const version)
double * global_dof_pt(const unsigned &i)
Return a pointer to the dof, indexed by global equation number which may be haloed or stored locally....
virtual void actions_after_read_unstructured_meshes()
Actions that are to be performed before reading in restart data for problems involving unstructured b...
void copy(Problem *orig_problem_pt)
Copy Data values, nodal positions etc from specified problem. Note: This is not a copy constructor....
virtual void get_jacobian(DoubleVector &residuals, DenseDoubleMatrix &jacobian)
Return the fully-assembled Jacobian and residuals for the problem Interface for the case when the Jac...
void set_explicit_time_stepper_pt(ExplicitTimeStepper *const &explicit_time_stepper_pt)
Set the explicit timestepper for the problem. The function will automatically create or resize the Ti...
Distributed_problem_matrix_distribution & distributed_problem_matrix_distribution()
accesss function to the distributed matrix distribution method 1 - Automatic - the Problem distributi...
virtual void actions_after_distribute()
Actions to be performed after a (mesh) distribution.
virtual void actions_before_implicit_timestep()
Actions that should be performed before each implicit time step. This is needed when one wants to sol...
LinearSolver *const & linear_solver_pt() const
Return a pointer to the linear solver object (const version)
void p_refine_uniformly_and_prune(const Vector< unsigned > &nrefine_for_mesh)
p-refine p-refineable sub-meshes, each as many times as specified in the vector and rebuild problem....
LinearSolver * Linear_solver_pt
Pointer to the linear solver for the problem.
virtual void actions_after_explicit_timestep()
Actions that should be performed after each explicit time step.
void disable_jacobian_reuse()
Disable recycling of Jacobian in Newton iteration.
bool Doc_time_in_distribute
Protected boolean flag to provide comprehensive timimings during problem distribution....
unsigned Max_newton_iterations
Maximum number of Newton iterations.
Vector< Vector< unsigned > > Sparse_assemble_with_arrays_previous_allocation
the number of elements in each row of a compressed matrix in the previous matrix assembly.
virtual void actions_after_parameter_increase(double *const ¶meter_pt)
Empty virtual function; provides hook to perform actions after the increase in the arclength paramete...
virtual void read(std::ifstream &restart_file)
Read refinement pattern of all refineable meshes and refine them accordingly, then read all Data and ...
void enable_doc_imbalance_in_parallel_assembly()
Enable doc of load imbalance in parallel assembly of distributed problem.
unsigned Dof_current_offset
Storage for the offset for the current values of dofs from the original dofs (default is 2,...
void calculate_continuation_derivatives(double *const ¶meter_pt)
A function to calculate the derivatives wrt the arc-length. This version of the function actually doe...
void store_current_dof_values()
Store the current values of the degrees of freedom.
bool Problem_has_been_distributed
Has the problem been distributed amongst multiple processors?
bool jacobian_reuse_is_enabled()
Is recycling of Jacobian in Newton iteration enabled?
void synchronise_all_dofs()
Perform all required synchronisation in solvers.
virtual void actions_after_change_in_bifurcation_parameter()
Actions that are to be performed after a change in the parameter that is being varied as part of the ...
bool Empty_actions_before_read_unstructured_meshes_has_been_called
Boolean to indicate that empty actions_before_read_unstructured_meshes() function has been called.
OomphCommunicator * communicator_pt()
access function to the oomph-lib communicator
bool Mass_matrix_has_been_computed
Has the mass matrix been computed (and can therefore be reused) Default: false.
unsigned nglobal_data() const
Return the number of global data values.
virtual void actions_before_adapt()
Actions that are to be performed before a mesh adaptation. These might include removing any additiona...
void newton_solve()
Use Newton method to solve the problem.
void(* SpatialErrorEstimatorFctPt)(Mesh *&mesh_pt, Vector< double > &elemental_error)
Function pointer for spatial error estimator.
bool First_jacobian_sign_change
Boolean to indicate whether a sign change has occured in the Jacobian.
void reset_arc_length_parameters()
Reset the "internal" arc-length continuation parameters, so as to allow continuation in another param...
void refine_distributed_base_mesh(Vector< Vector< Vector< unsigned >>> &to_be_refined_on_each_root, const unsigned &max_level_overall)
Load balance helper routine: refine each new base (sub)mesh based upon the elements to be refined wit...
void calculate_continuation_derivatives_fd_helper(double *const ¶meter_pt)
A function that performs the guts of the continuation derivative calculation in arc-length continuati...
LinearSolver *& linear_solver_pt()
Return a pointer to the linear solver object.
double Continuation_direction
The direction of the change in parameter that will ensure that a branch is followed in one direction ...
unsigned Parallel_sparse_assemble_previous_allocation
The amount of data allocated during the previous parallel sparse assemble. Used to optimise the next ...
void enable_mass_matrix_reuse()
Enable recycling of the mass matrix in explicit timestepping schemes. Useful for timestepping on fixe...
virtual void actions_before_explicit_timestep()
Actions that should be performed before each explicit time step.
double & target_error_safety_factor()
Access function to the safety factor in adaptive timestepping.
void steady_newton_solve(unsigned const &max_adapt=0)
Solve a steady problem using adaptive Newton's method, but in the context of an overall unstady probl...
void add_global_data(Data *const &global_data_pt)
Add Data to the Problem's global data – the Problem will perform equation numbering etc....
bool Scale_arc_length
Boolean to control whether arc-length should be scaled.
double doubly_adaptive_unsteady_newton_solve_helper(const double &dt, const double &epsilon, const unsigned &max_adapt, const unsigned &suppress_resolve_after_spatial_adapt, const bool &first, const bool &shift=true)
Private helper function that actually performs the unsteady "doubly" adaptive Newton solve....
bool Empty_actions_after_read_unstructured_meshes_has_been_called
Boolean to indicate that empty actions_after_read_unstructured_meshes() function has been called.
virtual void get_residuals(DoubleVector &residuals)
Return the fully-assembled residuals Vector for the problem: Virtual so it can be overloaded in for m...
Vector< GeneralisedElement * > Base_mesh_element_pt
Vector to store the correspondence between a root element and its element number within the global me...
double & dof(const unsigned &i)
i-th dof in the problem
void get_fd_jacobian(DoubleVector &residuals, DenseMatrix< double > &jacobian)
Return the fully-assembled Jacobian and residuals, generated by finite differences.
double doubly_adaptive_unsteady_newton_solve(const double &dt, const double &epsilon, const unsigned &max_adapt, const bool &first, const bool &shift=true)
Unsteady "doubly" adaptive Newton solve: Does temporal adaptation first, i.e. we try to do a timestep...
virtual void sparse_assemble_row_or_column_compressed_with_two_vectors(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)
Private helper function that is used to assemble the Jacobian matrix in the case when the storage is ...
TimeStepper *& time_stepper_pt(const unsigned &i)
Return a pointer to the i-th timestepper.
virtual Problem * make_copy()
Make and return a pointer to the copy of the problem. A virtual function that must be filled in by th...
void disable_problem_distributed()
Disable problem distributed.
int & sign_of_jacobian()
Access function for the sign of the global jacobian matrix. This will be set by the linear solver,...
bool & time_adaptive_newton_crash_on_solve_fail()
Access function for Time_adaptive_newton_crash_on_solve_fail.
double Maximum_dt
Maximum desired dt.
unsigned long set_timestepper_for_all_data(TimeStepper *const &time_stepper_pt, const bool &preserve_existing_data=false)
Set all problem data to have the same timestepper (timestepper_pt) Return the new number of dofs in t...
void activate_pitchfork_tracking(double *const ¶meter_pt, const DoubleVector &symmetry_vector, const bool &block_solve=true)
Turn on pitchfork tracking using the augmented system specified in the PitchForkHandler class....
virtual void dump(std::ofstream &dump_file) const
Dump refinement pattern of all refineable meshes and all generic Problem data to file for restart.
bool Use_finite_differences_for_continuation_derivatives
Boolean to specify which scheme to use to calculate the continuation derivatievs.
void p_refine_uniformly_and_prune(DocInfo &doc_info)
p-refine (all) p-refineable (sub)mesh(es) uniformly and rebuild problem; doc refinement process.
bool Arc_length_step_taken
Boolean to indicate whether an arc-length step has been taken.
void set_pinned_values_to_zero()
Set all pinned values to zero. Used to set boundary conditions to be homogeneous in the copy of the p...
bool Default_set_initial_condition_called
Has default set_initial_condition function been called? Default: false.
void get_bifurcation_eigenfunction(Vector< DoubleVector > &eigenfunction)
Return the eigenfunction calculated as part of a bifurcation tracking process. If we are not tracking...
double Relaxation_factor
Relaxation fator for Newton method (only a fractional Newton correction is applied if this is less th...
void add_time_stepper_pt(TimeStepper *const &time_stepper_pt)
Add a timestepper to the problem. The function will automatically create or resize the Time object so...
void refine_selected_elements(const Vector< unsigned > &elements_to_be_refined)
Refine (one and only!) mesh by splitting the elements identified by their numbers relative to the pro...
void prune_halo_elements_and_nodes(DocInfo &doc_info, const bool &report_stats)
(Irreversibly) prune halo(ed) elements and nodes, usually after another round of refinement,...
virtual void get_inverse_mass_matrix_times_residuals(DoubleVector &Mres)
Return the residual vector multiplied by the inverse mass matrix Virtual so that it can be overloaded...
void check_halo_schemes()
Check the halo/haloed node/element schemes.
double Parameter_current
Storage for the present value of the global parameter.
void enable_jacobian_reuse()
Enable recycling of Jacobian in Newton iteration (if the linear solver allows it)....
void assign_eigenvector_to_dofs(DoubleVector &eigenvector)
Assign the eigenvector passed to the function to the dofs in the problem so that it can be output by ...
Distributed_problem_matrix_distribution
enum for distribution of distributed jacobians. 1 - Automatic - the Problem distribution is employed,...
@ Uniform_matrix_distribution
@ Default_matrix_distribution
@ Problem_matrix_distribution
unsigned setup_element_count_per_dof()
Function that populates the Element_counter_per_dof vector with the number of elements that contribut...
void disable_globally_convergent_newton_method()
disable globally convergent Newton method
void problem_is_nonlinear(const bool &prob_lin)
Access function to Problem_is_nonlinear.
LinearSolver *& mass_matrix_solver_for_explicit_timestepper_pt()
Return a pointer to the linear solver object used for explicit time stepping.
void disable_doc_imbalance_in_parallel_assembly()
Disable doc of load imbalance in parallel assembly of distributed problem.
OomphCommunicator * Communicator_pt
The communicator for this problem.
double Newton_solver_tolerance
The Tolerance below which the Newton Method is deemed to have converged.
void get_flat_packed_refinement_pattern_for_load_balancing(const Vector< unsigned > &old_domain_for_base_element, const Vector< unsigned > &new_domain_for_base_element, const unsigned &max_refinement_level_overall, std::map< unsigned, Vector< unsigned >> &flat_packed_refinement_info_for_root)
Get flat-packed refinement pattern for each root element in current mesh (labeled by unique number of...
void set_consistent_pinned_values_for_continuation()
Private helper function that is used to set the appropriate pinned values for continuation.
void activate_bifurcation_tracking(double *const ¶meter_pt, const DoubleVector &eigenvector, const bool &block_solve=true)
Activate generic bifurcation tracking for a single (real) eigenvalue where the initial guess for the ...
void prune_halo_elements_and_nodes(const bool &report_stats=false)
(Irreversibly) prune halo(ed) elements and nodes, usually after another round of refinement,...
bool Discontinuous_element_formulation
Is the problem a discontinuous one, i.e. can the elemental contributions be treated independently....
bool Doc_imbalance_in_parallel_assembly
Boolean to switch on assessment of load imbalance in parallel assembly of distributed problem.
long synchronise_eqn_numbers(const bool &assign_local_eqn_numbers=true)
Classify any non-classified nodes into halo/haloed and synchronise equation numbers....
void create_new_linear_algebra_distribution(LinearAlgebraDistribution *&dist_pt)
Get new linear algebra distribution (you're in charge of deleting it!)
void p_refine_uniformly_aux(const Vector< unsigned > &nrefine_for_mesh, DocInfo &doc_info, const bool &prune)
Helper function to do compund p-refinement of (all) p-refineable (sub)mesh(es) uniformly as many time...
void calculate_continuation_derivatives_helper(const DoubleVector &z)
A function that performs the guts of the continuation derivative calculation in arc length continuati...
Vector< Problem * > Copy_of_problem_pt
Vector of pointers to copies of the problem used in adaptive bifurcation tracking problems (ALH: TEMP...
unsigned & max_newton_iterations()
Access function to max Newton iterations before giving up.
Vector< unsigned > distribute(const Vector< unsigned > &element_partition, DocInfo &doc_info, const bool &report_stats=false)
Distribute the problem and doc, using the specified partition; returns a vector which details the par...
double & newton_solver_tolerance()
Access function to tolererance of the Newton solver, i.e. the maximum value of the residuals that wil...
bool problem_has_been_distributed()
Access to Problem_has_been_distributed flag.
Mesh * Mesh_pt
The mesh pointer.
double Parameter_derivative
Storage for the derivative of the global parameter wrt arc-length.
void add_eigenvector_to_dofs(const double &epsilon, const DoubleVector &eigenvector)
Add the eigenvector passed to the function scaled by the constat epsilon to the dofs in the problem s...
Time * time_pt() const
Return a pointer to the global time object (const version).
Vector< double > elemental_assembly_time()
Return vector of most-recent elemental assembly times (used for load balancing). Zero sized if no Jac...
void flush_global_data()
Flush the Problem's global data – resizes container to zero. Data objects are not deleted!
Vector< Data * > Global_data_pt
Vector of global data: "Nobody" (i.e. none of the elements etc.) is "in charge" of this Data so it wo...
void recompute_load_balanced_assembly()
Helper function to re-assign the first and last elements to be assembled by each processor during par...
Vector< double * > Dof_pt
Vector of pointers to dofs.
void p_adapt()
p-adapt problem: Perform mesh adaptation for (all) refineable (sub)mesh(es), based on their own error...
double DTSF_max_increase
Maximum possible increase of dt between time-steps in adaptive schemes.
bool Must_recompute_load_balance_for_assembly
Boolean indicating that the division of elements over processors during the assembly process must be ...
virtual void shift_time_values()
Shift all values along to prepare for next timestep.
void p_refine_uniformly_and_prune(const Vector< unsigned > &nrefine_for_mesh, DocInfo &doc_info)
p-refine p-refineable sub-meshes, each as many times as specified in the vector and rebuild problem; ...
double Target_error_safety_factor
Safety factor to ensure we are aiming for a target error, TARGET, that is below our tolerance: TARGET...
Mesh *& mesh_pt()
Return a pointer to the global mesh.
void set_default_first_and_last_element_for_assembly()
Set default first and last elements for parallel assembly of non-distributed problem.
bool Keep_temporal_error_below_tolerance
Boolean to decide if a timestep is to be rejected if the error estimate post-solve (computed by globa...
bool Shut_up_in_newton_solve
Boolean to indicate if all output is suppressed in Problem::newton_solve(). Defaults to false.
void set_dofs(const DoubleVector &dofs)
Set the values of the dofs.
void setup_base_mesh_info_after_pruning()
Helper function to re-setup the Base_mesh enumeration (used during load balancing) after pruning.
Vector< Mesh * > Sub_mesh_pt
Vector of pointers to submeshes.
Vector< double > Max_res
Maximum residuals at start and after each newton iteration.
EigenSolver * Eigen_solver_pt
Pointer to the eigen solver for the problem.
bool Bisect_to_find_bifurcation
Boolean to control wheter bisection is used to located bifurcation.
void explicit_timestep(const double &dt, const bool &shift_values=true)
Take an explicit timestep of size dt and optionally shift any stored values of the time history.
void delete_all_external_storage()
Wrapper function to delete external storage for each submesh of the problem.
double DTSF_min_decrease
Minimum allowed decrease of dt between time-steps in adaptive schemes. Lower scaling values will reje...
virtual void symmetrise_eigenfunction_for_adaptive_pitchfork_tracking()
Virtual function that is used to symmetrise the problem so that the current solution exactly satisfie...
void load_balance(DocInfo &doc_info, const bool &report_stats)
Balance the load of a (possibly non-uniformly refined) problem that has already been distributed,...
void solve_eigenproblem(const unsigned &n_eval, Vector< std::complex< double >> &alpha, Vector< double > &beta, const bool &steady=true)
Solve an eigenproblem as assembled by the Problem's constituent elements but only return the eigenval...
double Minimum_dt_but_still_proceed
If Minimum_dt_but_still_proceed positive, then dt will not be reduced below this value during adaptiv...
double adaptive_unsteady_newton_solve(const double &dt_desired, const double &epsilon)
Attempt to advance timestep by dt_desired. If the solution fails the timestep will be halved until co...
unsigned Desired_newton_iterations_ds
The desired number of Newton Steps to reach convergence at each step along the arc.
void refine_uniformly_and_prune(DocInfo &doc_info)
Refine (all) refineable (sub)mesh(es) uniformly and rebuild problem; doc refinement process.
void rebuild_global_mesh()
If one of the submeshes has changed (e.g. by mesh adaptation) we need to update the global mesh....
void disable_store_local_dof_pt_in_elements()
Insist that local dof pointers are NOT set up in each element when equation numbering takes place (th...
ExplicitTimeStepper *& explicit_time_stepper_pt()
Return a pointer to the explicit timestepper.
bool Calculate_hessian_products_analytic
Map used to determine whether the hessian products should be computed using finite differences....
void set_first_and_last_element_for_assembly(Vector< unsigned > &first_el_for_assembly, Vector< unsigned > &last_el_for_assembly)
Manually set first and last elements for parallel assembly of non-distributed problem.
double & maximum_dt()
Access function to max timestep in adaptive timestepping.
void activate_hopf_tracking(double *const ¶meter_pt, const bool &block_solve=true)
Turn on Hopf bifurcation tracking using the augmented system specified in the HopfHandler class....
Vector< double * > Halo_dof_pt
Storage for the halo degrees of freedom (only required) when accessing via the global equation number...
void deactivate_bifurcation_tracking()
Deactivate all bifuraction tracking, by reseting the assembly handler to the default.
void globally_convergent_line_search(const Vector< double > &x_old, const double &half_residual_squared_old, DoubleVector &gradient, DoubleVector &newton_dir, double &half_residual_squared, const double &stpmax)
Line search helper for globally convergent Newton method.
void get_hessian_vector_products(DoubleVectorWithHaloEntries const &Y, Vector< DoubleVectorWithHaloEntries > const &C, Vector< DoubleVectorWithHaloEntries > &product)
Return the product of the global hessian (derivative of Jacobian matrix with respect to all variables...
virtual double global_temporal_error_norm()
Function to calculate a global error norm, used in adaptive timestepping to control the change in tim...
Assembly_method
Enumerated flags to determine which sparse assembly method is used.
@ Perform_assembly_using_two_arrays
@ Perform_assembly_using_maps
@ Perform_assembly_using_two_vectors
@ Perform_assembly_using_vectors_of_pairs
@ Perform_assembly_using_lists
void enable_globally_convergent_newton_method()
enable globally convergent Newton method
Mesh *& mesh_pt(const unsigned &imesh)
Return a pointer to the i-th submesh. If there are no submeshes, the 0-th submesh is the global mesh ...
int Sign_of_jacobian
Storage for the sign of the global Jacobian.
std::map< GeneralisedElement *, unsigned > Base_mesh_element_number_plus_one
Map which stores the correspondence between a root element and its element number (plus one) within t...
double Max_residuals
Maximum desired residual: if the maximum residual exceeds this value, the program will exit.
unsigned nsub_mesh() const
Return number of submeshes.
void unset_analytic_hessian_products()
Function to turn off analytic calculation of the parameter derivatives in continuation and bifurcatio...
double & time()
Return the current value of continuous time.
unsigned self_test()
Self-test: Check meshes and global data. Return 0 for OK.
unsigned ntime_stepper() const
Return the number of time steppers.
void activate_fold_tracking(double *const ¶meter_pt, const bool &block_solve=true)
Turn on fold tracking using the augmented system specified in the FoldHandler class....
double * dof_pt(const unsigned &i) const
Pointer to i-th dof in the problem (const version)
bool Use_globally_convergent_newton_method
Use the globally convergent newton method.
LinearSolver * Default_linear_solver_pt
Pointer to the default linear solver.
double Minimum_ds
Minimum desired value of arc-length.
void get_data_to_be_sent_during_load_balancing(const Vector< unsigned > &element_domain_on_this_proc, Vector< int > &send_n, Vector< double > &send_data, Vector< int > &send_displacement, Vector< unsigned > &old_domain_for_base_element, Vector< unsigned > &new_domain_for_base_element, unsigned &max_refinement_level_overall)
Load balance helper routine: Get data to be sent to other processors during load balancing and other ...
void restore_dof_values()
Restore the stored values of the degrees of freedom.
double Numerical_zero_for_sparse_assembly
A tolerance used to determine whether the entry in a sparse matrix is zero. If it is then storage nee...
double FD_step_used_in_get_hessian_vector_products
virtual void partition_global_mesh(Mesh *&global_mesh_pt, DocInfo &doc_info, Vector< unsigned > &element_domain, const bool &report_stats=false)
Partition the global mesh, return vector specifying the processor number for each element....
unsigned Sparse_assemble_with_arrays_initial_allocation
the number of elements to initially allocate for a matrix row within the sparse_assembly_with_two_arr...
void bifurcation_adapt_doc_errors(const unsigned &bifurcation_type)
A function that is used to document the errors used in the adaptive solution of bifurcation problems.
double Ds_current
Storage for the current step value.
virtual void set_initial_condition()
Set initial condition (incl previous timesteps). We need to establish this interface because I....
LinearSolver * Mass_matrix_solver_for_explicit_timestepper_pt
Pointer to the linear solver used for explicit time steps (this is likely to be different to the line...
void get_all_halo_data(std::map< unsigned, double * > &map_of_halo_data)
Get pointers to all possible halo data indexed by global equation number in a map.
void load_balance()
Balance the load of a (possibly non-uniformly refined) problem that has already been distributed,...
double & dof_derivative(const unsigned &i)
Access function to the derivative of the i-th (local) dof with respect to the arc length,...
double Minimum_dt
Minimum desired dt: if dt falls below this value, exit.
TimeStepper *& time_stepper_pt()
Access function for the pointer to the first (presumably only) timestepper.
double arc_length_step_solve(double *const ¶meter_pt, const double &ds, const unsigned &max_adapt=0)
Solve a steady problem using arc-length continuation, when the parameter that becomes a variable corr...
virtual void actions_after_adapt()
Actions that are to be performed after a mesh adaptation.
void p_refine_uniformly()
p-refine (all) p-refineable (sub)mesh(es) uniformly and rebuild problem
bool Use_continuation_timestepper
Boolean to control original or new storage of dof stuff.
void calculate_continuation_derivatives_fd(double *const ¶meter_pt)
A function to calculate the derivatives with respect to the arc-length required for continuation by f...
LinearAlgebraDistribution * Dof_distribution_pt
The distribution of the DOFs in this problem. This object is created in the Problem constructor and s...
void refine_uniformly()
Refine (all) refineable (sub)mesh(es) uniformly and rebuild problem.
static ContinuationStorageScheme Continuation_time_stepper
Storage for the single static continuation timestorage object.
bool Problem_is_nonlinear
Boolean flag indicating if we're dealing with a linear or nonlinear Problem – if set to false the New...
virtual void sparse_assemble_row_or_column_compressed_with_maps(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)
Private helper function that is used to assemble the Jacobian matrix in the case when the storage is ...
void refine_uniformly_and_prune(const Vector< unsigned > &nrefine_for_mesh, DocInfo &doc_info)
Refine refineable sub-meshes, each as many times as specified in the vector and rebuild problem; doc ...
void send_refinement_info_helper(Vector< unsigned > &old_domain_for_base_element, Vector< unsigned > &new_domain_for_base_element, const unsigned &max_refinement_level_overall, std::map< unsigned, Vector< unsigned >> &refinement_info_for_root_local, Vector< Vector< Vector< unsigned >>> &refinement_info_for_root_elements)
Send refinement information between processors.
double & minimum_dt()
Access function to min timestep in adaptive timestepping.
void operator=(const Problem &)=delete
Broken assignment operator.
AssemblyHandler * Default_assembly_handler_pt
Pointer to the default assembly handler.
void get_my_eqns(AssemblyHandler *const &assembly_handler_pt, const unsigned &el_lo, const unsigned &el_hi, Vector< unsigned > &my_eqns)
Helper method that returns the (unique) global equations to which the elements in the range el_lo to ...
std::map< double *, bool > Calculate_dparameter_analytic
Map used to determine whether the derivatives with respect to a parameter should be finite difference...
virtual void read(std::ifstream &restart_file, bool &unsteady_restart)
Read refinement pattern of all refineable meshes and refine them accordingly, then read all Data and ...
virtual ~Problem()
Virtual destructor to clean up memory.
void get_dofs(DoubleVector &dofs) const
Return the vector of dofs, i.e. a vector containing the current values of all unknowns.
Time * Time_pt
Pointer to global time for the problem.
DoubleVectorWithHaloEntries Element_count_per_dof
Counter that records how many elements contribute to each dof. Used to determine the (discrete) arc-l...
void solve_adjoint_eigenproblem(const unsigned &n_eval, Vector< std::complex< double >> &eigenvalue, const bool &steady=true)
Solve an adjoint eigenvalue problem using the same procedure as solve_eigenproblem but only return th...
void refine_uniformly(DocInfo &doc_info)
Refine (all) refineable (sub)mesh(es) uniformly and rebuild problem; doc refinement process.
Data *& global_data_pt(const unsigned &i)
Return a pointer to the the i-th global data object.
virtual void actions_before_newton_convergence_check()
Any actions that are to be performed before the residual is checked in the Newton method,...
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)
Unsteady "doubly" adaptive Newton solve: Does temporal adaptation first, i.e. we try to do a timestep...
unsigned newton_solve_continuation(double *const ¶meter_pt)
Perform a basic arc-length continuation step using Newton's method. Returns number of Newton steps ta...
double Max_permitted_error_for_halo_check
Threshold for error throwing in Problem::check_halo_schemes()
unsigned Sparse_assemble_with_arrays_allocation_increment
the number of elements to add to a matrix row when the initial allocation is exceeded within the spar...
void reset_assembly_handler_to_default()
Reset the system to the standard non-augemented state.
void set_analytic_dparameter(double *const ¶meter_pt)
Function to turn on analytic calculation of the parameter derivatives in continuation and bifurcation...
virtual void actions_after_newton_step()
Any actions that are to be performed after each individual Newton step. Most likely to be used for di...
bool Pause_at_end_of_sparse_assembly
Protected boolean flag to halt program execution during sparse assemble process to assess peak memory...
void remove_null_pointers_from_external_halo_node_storage()
Consolidate external halo node storage by removing nulled out pointers in external halo and haloed sc...
Vector< double > * Saved_dof_pt
Pointer to vector for backup of dofs.
void unsteady_newton_solve(const double &dt)
Advance time by dt and solve by Newton's method. This version always shifts time values.
AssemblyHandler * Assembly_handler_pt
void disable_info_in_newton_solve()
Disable the output of information when in the newton solver.
unsigned Dof_derivative_offset
Storage for the offset for the continuation derivatives from the original dofs (default is 1,...
void get_all_error_estimates(Vector< Vector< double >> &elemental_error)
Return the error estimates computed by (all) refineable (sub)mesh(es) in the elemental_error structur...
virtual void actions_before_read_unstructured_meshes()
Actions that are to be performed before reading in restart data for problems involving unstructured b...
const OomphCommunicator * communicator_pt() const
access function to the oomph-lib communicator, const version
Vector< TimeStepper * > Time_stepper_pt
The Vector of time steppers (there could be many different ones in multiphysics problems)
void doc_errors()
Get max and min error for all elements in submeshes.
bool Mass_matrix_reuse_is_enabled
Is re-use of the mass matrix in explicit timestepping enabled Default:false.
void get_derivative_wrt_global_parameter(double *const ¶meter_pt, DoubleVector &result)
Get the derivative of the entire residuals vector wrt a global parameter, used in continuation proble...
void enable_discontinuous_formulation()
Indicate that the problem involves discontinuous elements This allows for a more efficiently assembly...
void adapt_based_on_error_estimates(unsigned &n_refined, unsigned &n_unrefined, Vector< Vector< double >> &elemental_error)
Adapt problem: Perform mesh adaptation for (all) refineable (sub)mesh(es), based on the error estimat...
bool distributed() const
If we have MPI return the "problem has been distributed" flag, otherwise it can't be distributed so r...
bool are_hessian_products_calculated_analytically()
Function to determine whether the hessian products are calculated analytically.
bool does_pointer_correspond_to_problem_data(double *const ¶meter_pt)
Return a boolean flag to indicate whether the pointer parameter_pt refers to values stored in a Data ...
ExplicitTimeStepper * Explicit_time_stepper_pt
Pointer to a single explicit timestepper.
void dump(const std::string &dump_file_name) const
Dump refinement pattern of all refineable meshes and all generic Problem data to file for restart.
bool mass_matrix_reuse_is_enabled()
Return whether the mass matrix is being reused.
void disable_discontinuous_formulation()
Disable the use of a discontinuous-element formulation. Note that the methods will all still work eve...
double arc_length_step_solve_helper(double *const ¶meter_pt, const double &ds, const unsigned &max_adapt)
Private helper function that actually contains the guts of the arc-length stepping,...
void refine_uniformly(const unsigned &i_mesh)
Do uniform refinement for submesh i_mesh without documentation.
void solve_eigenproblem(const unsigned &n_eval, Vector< std::complex< double >> &alpha, Vector< double > &beta, Vector< DoubleVector > &eigenvector_real, Vector< DoubleVector > &eigenvector_imag, const bool &steady=true)
Solve an eigenproblem as assembled by the Problem's constituent elements. Calculate (at least) n_eval...
bool Use_predictor_values_as_initial_guess
Use values from the time stepper predictor as an initial guess.
void unset_analytic_dparameter(double *const ¶meter_pt)
Function to turn off analytic calculation of the parameter derivatives in continuation and bifurcatio...
Class for a matrix of the form M = S + G + H + ... where S is the main matrix and G,...
////////////////////////////////////////////////////////////////////// //////////////////////////////...
Class to keep track of discrete/continous time. It is essential to have a single Time object when usi...
A slight extension to the standard template vector class so that we can include "graceful" array rang...
std::string string(const unsigned &i)
Return the i-th string or "" if the relevant string hasn't been defined.
//////////////////////////////////////////////////////////////////// ////////////////////////////////...
OomphInfo oomph_info
Single (global) instantiation of the OomphInfo object – this is used throughout the library as a "rep...