BubbleProfiler
0.3.0
by Peter Athron, Csaba Balazs, Michael Bardsley, Andrew Fowlie, Dylan Harries & Graham White
|
Computes an initial guess for the bubble profile using a kink ansatz. More...
#include <kink_profile_guesser.hpp>
Public Member Functions | |
virtual | ~Kink_profile_guesser ()=default |
virtual Field_profiles | get_profile_guess (const Potential &, const Eigen::VectorXd &, int, double, double, double, double) override |
Calculate an initial guess for the bubble profile. More... | |
void | set_trial_distance (double dist) |
Sets the trial distance used in estimating the potential parameters. More... | |
void | set_max_domain_start (double r) |
Sets the maximum allowed distance if guessing the domain start. More... | |
void | set_min_domain_end (double r) |
Sets the minimum required distance if guessing the domain end. More... | |
![]() | |
virtual | ~Profile_guesser ()=default |
Private Member Functions | |
void | compute_vacuum_distance (const Potential &, const Eigen::VectorXd &) |
void | calculate_potential_parameters (const Potential &, const Eigen::VectorXd &) |
void | fit_ansatz_parameters (const Potential &, const Eigen::VectorXd &, int) |
double | evaluate_ansatz_at (double) const |
double | evaluate_ansatz_deriv_at (double) const |
double | guess_domain_start () const |
double | guess_domain_end () const |
Field_profiles | calculate_field_profiles (int, double, double, double, double) |
Private Attributes | |
double | trial_dist {0.5} |
double | max_domain_start {1.e-8} |
double | min_domain_end {4.} |
double | alpha {0.5} |
double | aE {1.} |
double | phi0 {0.} |
double | delta {0.} |
double | lw {0.} |
double | dist_true_vacuum {0.} |
int | num_fields {0} |
Eigen::MatrixXd | cob_matrix {} |
logging::Basic_logger | logger {} |
double | alpha_threshold = 0.514 |
Minimum value of alpha before throwing a Thin_wall_error. More... | |
Static Private Attributes | |
static const std::array< double, 46 > | alpha_grid_3d |
static const std::array< double, 50 > | alpha_grid_thin_3d |
static const std::array< double, 46 > | phi0_grid_3d |
static const std::array< double, 46 > | delta_grid_3d |
static const std::array< double, 50 > | delta_grid_thin_3d |
static const std::array< double, 46 > | lw_grid_3d |
static const std::array< double, 44 > | alpha_grid_4d |
static const std::array< double, 50 > | alpha_grid_thin_4d |
static const std::array< double, 44 > | phi0_grid_4d |
static const std::array< double, 44 > | delta_grid_4d |
static const std::array< double, 50 > | delta_grid_thin_4d |
static const std::array< double, 44 > | lw_grid_4d |
Computes an initial guess for the bubble profile using a kink ansatz.
Definition at line 34 of file kink_profile_guesser.hpp.
|
virtualdefault |
|
private |
Definition at line 386 of file kink_profile_guesser.cpp.
References cob_matrix, dist_true_vacuum, evaluate_ansatz_at(), guess_domain_end(), guess_domain_start(), num_fields, and BubbleProfiler::Field_profiles::set_number_of_dimensions().
Referenced by get_profile_guess().
|
private |
Definition at line 324 of file kink_profile_guesser.cpp.
References BubbleProfiler::Abs(), aE, alpha, alpha_threshold, BubbleProfiler::calculate_rotation_to_target(), BubbleProfiler::Potential::clone(), cob_matrix, dist_true_vacuum, BubbleProfiler::logging::Basic_logger::log_message(), logger, num_fields, BubbleProfiler::logging::Trace, and trial_dist.
Referenced by fit_ansatz_parameters().
|
private |
Definition at line 218 of file kink_profile_guesser.cpp.
References BubbleProfiler::Abs(), dist_true_vacuum, BubbleProfiler::Potential::get_number_of_fields(), and num_fields.
Referenced by get_profile_guess().
|
private |
Definition at line 437 of file kink_profile_guesser.cpp.
References aE, delta, dist_true_vacuum, BubbleProfiler::Exp(), lw, phi0, BubbleProfiler::Sech(), BubbleProfiler::Sqrt(), and BubbleProfiler::Tanh().
Referenced by calculate_field_profiles().
|
private |
Definition at line 455 of file kink_profile_guesser.cpp.
References aE, delta, dist_true_vacuum, BubbleProfiler::Exp(), lw, phi0, BubbleProfiler::Sech(), and BubbleProfiler::Sqrt().
Referenced by guess_domain_start().
|
private |
Definition at line 238 of file kink_profile_guesser.cpp.
References alpha, alpha_grid_3d, alpha_grid_4d, alpha_grid_thin_3d, alpha_grid_thin_4d, calculate_potential_parameters(), delta, delta_grid_3d, delta_grid_4d, delta_grid_thin_3d, delta_grid_thin_4d, dist_true_vacuum, BubbleProfiler::logging::Basic_logger::log_message(), logger, lw, lw_grid_3d, lw_grid_4d, phi0, phi0_grid_3d, phi0_grid_4d, BubbleProfiler::quadratic_lagrange_interpolation_at_point(), and BubbleProfiler::logging::Trace.
Referenced by get_profile_guess().
|
overridevirtual |
Calculate an initial guess for the bubble profile.
potential | the potential for which the profiles are to be computed |
true_vacuum | the location of the true vacuum in field space |
n_dimensions | the number of space-time dimensions |
domain_start | the minimum value of the radial coordinate at which the profiles are evaluated. If negative, a value is guessed based on the computed ansatz properties. |
domain_end | the maximum value of the radial coordinate at which the profiles are evaluated. If negative, a value is guessed based on the computed ansatz properties. |
initial_step_size | the initial step size to be used in constructing the discretized solution |
interpolation_points_fraction | fraction of grid points to use for interpolation |
The initial guess is based on a kink ansatz corresponding to the fitted solution of an approximate one-dimensional problem. This is obtained by first rotating the given field basis into one in which the vector from the false to the true vacuum lies along the first axis in field space. The function then fits a quartic potential to the given potential along the line between the two minima. The returned guess is obtained from the fitted solution to this one-dimensional problem, rotated back into the initial field basis.
If the start of the integration domain is to be guessed, the guessed value is determined by finding the radial distance at which the derivative of the one-dimensional fitted profile is equal to 1e-5.
If the end of the integration domain is to be guessed, the guessed value is determined by finding the radial distance at which the value of the one-dimensional fitted profile is less than 1e-5.
Implements BubbleProfiler::Profile_guesser.
Definition at line 311 of file kink_profile_guesser.cpp.
References calculate_field_profiles(), compute_vacuum_distance(), and fit_ansatz_parameters().
|
private |
Definition at line 514 of file kink_profile_guesser.cpp.
References aE, BubbleProfiler::ArcTanh(), delta, dist_true_vacuum, BubbleProfiler::logging::Basic_logger::log_message(), logger, lw, min_domain_end, phi0, BubbleProfiler::Sqrt(), and BubbleProfiler::logging::Trace.
Referenced by calculate_field_profiles().
|
private |
Definition at line 471 of file kink_profile_guesser.cpp.
References BubbleProfiler::Abs(), evaluate_ansatz_deriv_at(), BubbleProfiler::logging::Basic_logger::log_message(), logger, max_domain_start, and BubbleProfiler::logging::Trace.
Referenced by calculate_field_profiles().
void BubbleProfiler::Kink_profile_guesser::set_max_domain_start | ( | double | r | ) |
Sets the maximum allowed distance if guessing the domain start.
If the start of the integration domain is to be guessed instead of being explicitly provided, this function sets the maximum distance away from the origin that the returned profiles are allowed to start. If the guessed starting distance would be larger than the value of r
, the provided value of r
is used instead.
r | the maximum allowed distance from the origin to start from |
Definition at line 202 of file kink_profile_guesser.cpp.
References max_domain_start.
Referenced by set_trial_distance().
void BubbleProfiler::Kink_profile_guesser::set_min_domain_end | ( | double | r | ) |
Sets the minimum required distance if guessing the domain end.
If the end of the integration domain is to be guessed instead of being explicitly provided, this function sets the minimum distance away from the origin that the returned profiles are required to reach. If the guessed domain end would otherwise be less than the value of r
, the provided value of r
is used instead.
r | the minimum distance that the profiles must extend to |
Definition at line 210 of file kink_profile_guesser.cpp.
References min_domain_end.
Referenced by set_trial_distance().
|
inline |
Sets the trial distance used in estimating the potential parameters.
The parameters defining the quartic potential from which the ansatz is determined are computed by evaluating the potential at a point on the line between the local and global minima. Initially, the trial point is taken to be halfway along this line. The distance of the point along this line can be changed using this method.
dist | the distance between the minima to use |
Definition at line 73 of file kink_profile_guesser.hpp.
References set_max_domain_start(), set_min_domain_end(), and trial_dist.
|
private |
Absolute value of E in the fitted one-dimensional quartic potential
Definition at line 137 of file kink_profile_guesser.hpp.
Referenced by calculate_potential_parameters(), evaluate_ansatz_at(), evaluate_ansatz_deriv_at(), and guess_domain_end().
|
private |
Value of alpha in the fitted one-dimensional quartic potential
Definition at line 135 of file kink_profile_guesser.hpp.
Referenced by calculate_potential_parameters(), and fit_ansatz_parameters().
|
staticprivate |
Values of alpha used for interpolation in 3 spacetime dimensions
Definition at line 103 of file kink_profile_guesser.hpp.
Referenced by fit_ansatz_parameters().
|
staticprivate |
Values of alpha used for interpolation in 4 spacetime dimensions
Definition at line 116 of file kink_profile_guesser.hpp.
Referenced by fit_ansatz_parameters().
|
staticprivate |
Values of alpha used for thin-wall interpolation in 3 spacetime dimensions
Definition at line 105 of file kink_profile_guesser.hpp.
Referenced by fit_ansatz_parameters().
|
staticprivate |
Values of alpha used for thin-wall interpolation in 4 spacetime dimensions
Definition at line 118 of file kink_profile_guesser.hpp.
Referenced by fit_ansatz_parameters().
|
private |
Minimum value of alpha before throwing a Thin_wall_error.
Definition at line 167 of file kink_profile_guesser.hpp.
Referenced by calculate_potential_parameters().
|
private |
Change of (field space) basis matrix, takes vectors in the original basis to those in the ansatz basis.
Definition at line 153 of file kink_profile_guesser.hpp.
Referenced by calculate_field_profiles(), and calculate_potential_parameters().
|
private |
Computed value of delta in the ansatz solution
Definition at line 141 of file kink_profile_guesser.hpp.
Referenced by evaluate_ansatz_at(), evaluate_ansatz_deriv_at(), fit_ansatz_parameters(), and guess_domain_end().
|
staticprivate |
Value of delta for each value of alpha in 3 spacetime dimensions
Definition at line 109 of file kink_profile_guesser.hpp.
Referenced by fit_ansatz_parameters().
|
staticprivate |
Value of delta for each value of alpha in 4 spacetime dimensions
Definition at line 122 of file kink_profile_guesser.hpp.
Referenced by fit_ansatz_parameters().
|
staticprivate |
Value of delta for each value of alpha in 3 spacetime dimensions, thin walled case
Definition at line 111 of file kink_profile_guesser.hpp.
Referenced by fit_ansatz_parameters().
|
staticprivate |
Value of delta for each value of alpha in 4 spacetime dimensions, thin walled case
Definition at line 124 of file kink_profile_guesser.hpp.
Referenced by fit_ansatz_parameters().
|
private |
Distance to the global minimum from the origin in field space
Definition at line 145 of file kink_profile_guesser.hpp.
Referenced by calculate_field_profiles(), calculate_potential_parameters(), compute_vacuum_distance(), evaluate_ansatz_at(), evaluate_ansatz_deriv_at(), fit_ansatz_parameters(), and guess_domain_end().
|
private |
Definition at line 154 of file kink_profile_guesser.hpp.
Referenced by calculate_potential_parameters(), fit_ansatz_parameters(), guess_domain_end(), and guess_domain_start().
|
private |
Computed value of Lw in the ansatz solution
Definition at line 143 of file kink_profile_guesser.hpp.
Referenced by evaluate_ansatz_at(), evaluate_ansatz_deriv_at(), fit_ansatz_parameters(), and guess_domain_end().
|
staticprivate |
Value of Lw for each value of alpha in 3 spacetime dimensions
Definition at line 113 of file kink_profile_guesser.hpp.
Referenced by fit_ansatz_parameters().
|
staticprivate |
Value of Lw for each value of alpha in 4 spacetime dimensions
Definition at line 126 of file kink_profile_guesser.hpp.
Referenced by fit_ansatz_parameters().
|
private |
Maximum allowed starting distance from origin
Definition at line 131 of file kink_profile_guesser.hpp.
Referenced by guess_domain_start(), and set_max_domain_start().
|
private |
Minimum required endpoint for guessed profile
Definition at line 133 of file kink_profile_guesser.hpp.
Referenced by guess_domain_end(), and set_min_domain_end().
|
private |
Number of fields in the potential
Definition at line 147 of file kink_profile_guesser.hpp.
Referenced by calculate_field_profiles(), calculate_potential_parameters(), and compute_vacuum_distance().
|
private |
Computed value of phi0 in the ansatz solution
Definition at line 139 of file kink_profile_guesser.hpp.
Referenced by evaluate_ansatz_at(), evaluate_ansatz_deriv_at(), fit_ansatz_parameters(), and guess_domain_end().
|
staticprivate |
Value of phi0 for each value of alpha in 3 spacetime dimensions
Definition at line 107 of file kink_profile_guesser.hpp.
Referenced by fit_ansatz_parameters().
|
staticprivate |
Value of phi0 for each value of alpha in 4 spacetime dimensions
Definition at line 120 of file kink_profile_guesser.hpp.
Referenced by fit_ansatz_parameters().
|
private |
Distance of trial point between minima
Definition at line 129 of file kink_profile_guesser.hpp.
Referenced by calculate_potential_parameters(), and set_trial_distance().