Skip to content

Commit

Permalink
put two-scale specific solvers into separate files
Browse files Browse the repository at this point in the history
  • Loading branch information
Alexander Voigt authored and Alexander Voigt committed Jul 2, 2015
1 parent 92bae33 commit eab378d
Show file tree
Hide file tree
Showing 11 changed files with 516 additions and 319 deletions.
10 changes: 9 additions & 1 deletion meta/FlexibleSUSY.m
View file
Expand Up @@ -1947,7 +1947,15 @@ corresponding tadpole is real or imaginary (only in models with CP
{{FileNameJoin[{Global`$flexiblesusyTemplateDir, "ewsb_solver.hpp.in"}],
FileNameJoin[{Global`$flexiblesusyOutputDir, FlexibleSUSY`FSModelName <> "_ewsb_solver.hpp"}]},
{FileNameJoin[{Global`$flexiblesusyTemplateDir, "ewsb_solver.cpp.in"}],
FileNameJoin[{Global`$flexiblesusyOutputDir, FlexibleSUSY`FSModelName <> "_ewsb_solver.cpp"}]}
FileNameJoin[{Global`$flexiblesusyOutputDir, FlexibleSUSY`FSModelName <> "_ewsb_solver.cpp"}]},
{FileNameJoin[{Global`$flexiblesusyTemplateDir, "two_scale_ewsb_solver.hpp.in"}],
FileNameJoin[{Global`$flexiblesusyOutputDir, FlexibleSUSY`FSModelName <> "_two_scale_ewsb_solver.hpp"}]},
{FileNameJoin[{Global`$flexiblesusyTemplateDir, "two_scale_ewsb_solver.cpp.in"}],
FileNameJoin[{Global`$flexiblesusyOutputDir, FlexibleSUSY`FSModelName <> "_two_scale_ewsb_solver.cpp"}]},
{FileNameJoin[{Global`$flexiblesusyTemplateDir, "two_scale_soft_mass_ewsb_solver.hpp.in"}],
FileNameJoin[{Global`$flexiblesusyOutputDir, FlexibleSUSY`FSModelName <> "_two_scale_soft_mass_ewsb_solver.hpp"}]},
{FileNameJoin[{Global`$flexiblesusyTemplateDir, "two_scale_soft_mass_ewsb_solver.cpp.in"}],
FileNameJoin[{Global`$flexiblesusyOutputDir, FlexibleSUSY`FSModelName <> "_two_scale_soft_mass_ewsb_solver.cpp"}]}
}];

PrintHeadline["Creating SLHA model"];
Expand Down
286 changes: 31 additions & 255 deletions templates/ewsb_solver.cpp.in
View file
Expand Up @@ -19,299 +19,75 @@
// File generated at @DateAndTime@

#include "@ModelName@_ewsb_solver.hpp"
#include "@ModelName@_mass_eigenstates.hpp"

#include "ewsb_solver.hpp"
#include "fixed_point_iterator.hpp"
#include "functors.hpp"
#include "logger.hpp"
#include "minimizer.hpp"
#include "numerics2.hpp"
#include "root_finder.hpp"
#include "wrappers.hpp"
#include "gsl_utils.hpp"
#include <gsl/gsl_multiroots.h>

#define CLASSNAME @ModelName@_ewsb_solver
#define MODEL get_model()
#define MODELPARAMETER(parameter) MODEL->get_##parameter()
#define INPUT(parameter) model->get_input().parameter
#define LOCALINPUT(parameter) MODEL->get_input().parameter
#define INPUTPARAMETER(parameter) LOCALINPUT(parameter)
#define PHASE(p) MODELPARAMETER(p)

namespace flexiblesusy {

/**
* Solves EWSB at given loop order.
*
* @param loop_order_ temporary loop order
*
* @return return value of solve()
*/
int @ModelName@_ewsb_solver_interface::solve(unsigned loop_order_)
{
const unsigned old_loop_order = loop_order;
loop_order = loop_order_;
solve();
loop_order = old_loop_order;
}

@ModelName@_ewsb_solver_tree_level_soft_masses::@ModelName@_ewsb_solver_tree_level_soft_masses(@ModelName@_mass_eigenstates* model_)
: @ModelName@_ewsb_solver_interface(model_)
{
assert(model_ && "@ModelName@_ewsb_solver_tree_level_soft_masses:"
" model pointer must not be zero");

@saveSoftHiggsMasses@
}

@ModelName@_ewsb_solver_tree_level_soft_masses::~@ModelName@_ewsb_solver_tree_level_soft_masses()
{
@restoreSoftHiggsMasses@
}

int @ModelName@_ewsb_solver_tree_level_soft_masses::solve()
{
int error = 0;

@solveTreeLevelEWSBviaSoftHiggsMasses@

return error;
}



@ModelName@_ewsb_solver::@ModelName@_ewsb_solver(@ModelName@_mass_eigenstates* model_)
: @ModelName@_ewsb_solver_interface(model_)
@ModelName@_ewsb_solver::@ModelName@_ewsb_solver(
@ModelName@_mass_eigenstates* model_
)
: number_of_iterations(100)
, loop_order(2)
, precision(1.0e-5)
, model(model_)
{
}

@ModelName@_ewsb_solver::~@ModelName@_ewsb_solver()
{
}

/**
* This method solves the EWSB conditions iteratively, trying several
* root finding methods until a solution is found.
*/
int @ModelName@_ewsb_solver::solve_iteratively()
void @ModelName@_ewsb_solver::set_loop_order(unsigned n)
{
EWSB_args params = {get_model(), this, get_loop_order()};

EWSB_solver* solvers[] = {
@EWSBSolvers@
};

const std::size_t number_of_solvers = sizeof(solvers)/sizeof(*solvers);
double x_init[number_of_ewsb_equations] = { 0. };
initial_guess(x_init);

#ifdef ENABLE_VERBOSE
std::cout << "Solving EWSB equations ...\n"
"\tInitial guess: x_init =";
for (std::size_t i = 0; i < number_of_ewsb_equations; ++i)
std::cout << ' ' << x_init[i];
std::cout << '\n';
#endif

int status;
for (std::size_t i = 0; i < number_of_solvers; ++i) {
VERBOSE_MSG("\tStarting EWSB iteration using solver " << i);
status = solve_iteratively_with(solvers[i], x_init);
if (status == EWSB_solver::SUCCESS) {
VERBOSE_MSG("\tSolver " << i << " finished successfully!");
break;
}
#ifdef ENABLE_VERBOSE
else {
WARNING("\tSolver " << i << " could not find a solution!"
" (requested precision: " << iteration_precision << ")");
}
#endif
}

if (status == EWSB_solver::SUCCESS) {
MODEL->get_problems().unflag_no_ewsb();
} else {
MODEL->get_problems().flag_no_ewsb();
#ifdef ENABLE_VERBOSE
WARNING("\tCould not find a solution to the EWSB equations!"
" (requested precision: " << iteration_precision << ")");
#endif
}

for_each(solvers, solvers + number_of_solvers, Delete_object());

return status;
loop_order = n;
}

/**
* Solves EWSB equations with given EWSB solver
*
* @param solver EWSB solver
* @param x_init initial values
*
* @return status of the EWSB solver
*/
int @ModelName@_ewsb_solver::solve_iteratively_with(
EWSB_solver* solver,
const double x_init[number_of_ewsb_equations]
)
void @ModelName@_ewsb_solver::set_model(@ModelName@_mass_eigenstates* m)
{
const int status = solver->solve(x_init);

@setEWSBSolution@

return status;
model = m;
}

int @ModelName@_ewsb_solver::solve_iteratively(unsigned loop_order_)
void @ModelName@_ewsb_solver::set_number_of_iterations(std::size_t n)
{
// temporarily set `loop_order' to `loop_order_' and do
// iteration
const unsigned old_loop_order = get_loop_order();
set_loop_order(loop_order_);
const int status = solve_iteratively();
set_loop_order(old_loop_order);
return status;
number_of_iterations = n;
}


int @ModelName@_ewsb_solver::solve_tree_level()
{
int error = 0;

@solveEwsbTreeLevel@

return error;
void @ModelName@_ewsb_solver::set_precision(double p) {
precision = p;
}

int @ModelName@_ewsb_solver::solve()
unsigned @ModelName@_ewsb_solver::get_loop_order() const
{
VERBOSE_MSG("\tSolving EWSB at " << get_loop_order() << "-loop order");

if (get_loop_order() == 0)
return solve_tree_level();

return solve_iteratively(get_loop_order());
return loop_order;
}

void @ModelName@_ewsb_solver::initial_guess(double x_init[number_of_ewsb_equations])
@ModelName@_mass_eigenstates* @ModelName@_ewsb_solver::get_model() const
{
@ewsbInitialGuess@
return model;
}

/**
* Calculates EWSB output parameters including loop-corrections.
*
* @param parameters new EWSB output parameters. \a
* parameters is only modified if all new parameters are finite.
*
* @return GSL_SUCCESS if new EWSB output parameters are finite,
* GSL_EDOM otherwise.
*/
int @ModelName@_ewsb_solver::step(double parameters[number_of_ewsb_equations]) const
unsigned @ModelName@_ewsb_solver::get_number_of_iterations() const
{
int error;
double tadpole[number_of_ewsb_equations] = { 0. };

if (get_loop_order() > 0) {
@calculateOneLoopTadpolesFromModel@
if (get_loop_order() > 1) {
@calculateTwoLoopTadpolesFromModel@
}
}

@solveEwsbWithTadpoles@

if (is_finite) {
error = GSL_SUCCESS;
@fillArrayWithEWSBParameters@
} else {
error = GSL_EDOM;
}

return error;
return number_of_iterations;
}

#undef MODELPARAMETER
#define MODELPARAMETER(parameter) model->get_##parameter()

/**
* Calculates EWSB output parameters including loop-corrections.
*
* @param x old EWSB output parameters
* @param params further function parameters
* @param f new EWSB output parameters
*
* @return Returns status of @ModelName@_ewsb_solver::step
*/
int @ModelName@_ewsb_solver::step(const gsl_vector* x, void* params, gsl_vector* f)
double @ModelName@_ewsb_solver::get_precision() const
{
if (!is_finite(x)) {
gsl_vector_set_all(f, std::numeric_limits<double>::max());
return GSL_EDOM;
}

const @ModelName@_ewsb_solver::EWSB_args* args
= static_cast<@ModelName@_ewsb_solver::EWSB_args*>(params);
@ModelName@_mass_eigenstates* model = args->model;
@ModelName@_ewsb_solver* solver = args->solver;
const unsigned loop_order = args->loop_order;

@getEWSBParametersFromGSLVector@
@setEWSBParametersFromLocalCopies@

if (loop_order > 0)
model->calculate_DRbar_masses();

double parameters[number_of_ewsb_equations] =
@ewsbParametersInitializationList@;

const int status = solver->step(parameters);

for (std::size_t i = 0; i < number_of_ewsb_equations; ++i)
gsl_vector_set(f, i, parameters[i]);

return status;
return precision;
}

/**
* Method which calculates the tadpoles at loop order specified in the
* pointer to the @ModelName@_ewsb_solver::EWSB_args struct.
* Solves EWSB at given loop order.
*
* @param x GSL vector of EWSB output parameters
* @param params pointer to @ModelName@_ewsb_solver::EWSB_args struct
* @param f GSL vector with tadpoles
* @param loop_order_ temporary loop order
*
* @return GSL_EDOM if x contains Nans, GSL_SUCCESS otherwise.
* @return return value of solve()
*/
int @ModelName@_ewsb_solver::tadpole_equations(const gsl_vector* x, void* params, gsl_vector* f)
int @ModelName@_ewsb_solver::solve(unsigned loop_order_)
{
if (!is_finite(x)) {
gsl_vector_set_all(f, std::numeric_limits<double>::max());
return GSL_EDOM;
}

const @ModelName@_ewsb_solver::EWSB_args* args
= static_cast<@ModelName@_ewsb_solver::EWSB_args*>(params);
@ModelName@_mass_eigenstates* model = args->model;
const unsigned loop_order = args->loop_order;

@setEWSBParametersFromGSLVector@

if (loop_order > 0)
model->calculate_DRbar_masses();

double tadpole[number_of_ewsb_equations] = { 0. };

model->tadpole_equations(tadpole);

for (std::size_t i = 0; i < number_of_ewsb_equations; ++i)
gsl_vector_set(f, i, tadpole[i]);

return is_finite<number_of_ewsb_equations>(tadpole) ? GSL_SUCCESS : GSL_EDOM;
const unsigned old_loop_order = loop_order;
loop_order = loop_order_;
solve();
loop_order = old_loop_order;
}

} // namespace flexiblesusy

0 comments on commit eab378d

Please sign in to comment.