//The purpose of this function is to provide for a parametrization of the VEVs such that they can be generated in more generic ways other than just in a hypercubic box. For example in the 2HDM model we may want to generate h_1=246*cos(beta) and h_2=246*sin(beta), i.e. on a circle, using beta as a parameter.
//Here you must write down such a relation in the form Phi[i]=function(PhiPar[0],PhiPar[1],...). The parameter ranges have already been specified in the mathematica notebook.
// MixInternal[i][j] Used defined mixing matrix for the mixing block. i runs over the mass eigenstates and j runs from 0 to Nmixing-1 (Nmixing is the number or mixing states)
// r random number generator (called as r(mix,max)) in case the user whishes to change here the scan of MixPars[] more flexibly.
///////////////////////////
// Description:
//The purpose of this function is to provide for a parametrization of the mixing matrix when for simpler cases where an analytic parametrization is natural, just as for example in the 2HDM. Since at this point in the scan the VEVs have been generated (and possibly their parameters, if they have been re-parametrized), we assume a general structure where the mixing matrix parametrization can depend on such (constant) quantities (i.e. Phi, the VEVS and PhiPar, the parameters) and a new set of parameters (MixPars) whose ranges have been specified in the Mathematica notebook.
//Contrarily to the VEV re-parametrization, we allow here the use the option to actually change the values of MixPars (i.e. override the values that have been generated in the ranges defined in the mathematica notebook) so that relations between Mixing parameters and VEV parameters are possible (such as beta=alpha-Pi/2 in the decoupling limit of the 2HDM).
//Here you must write down:
// i) (this one is optional) a relation MixPars[i]=function(Phi[0],...,PhiPars[0],...)
MixPars[0]=PhiPar[0]-acos(-1)/2e0;
// ii) Expressions for the mixing matrix elements in the form MixInternal[i][j]=Function(MixPars[0],...,Phi[0],...,PhiPars[0],...), where the index i runs over the scalar eigenstates and the index j runs over the original mixing scalar states according to the order defined in the mathematica file. IMPORTANT NOTE: You are advised to check (in VERBOSE mode) that the ordering is indeed as you think it is. Note for example that the mass eigenstates are only ordered in mass if you enforce this in the Mathematica file! The work around for example in the 2HDM mode is to split the scan in 2 scans: The first, where the SM Higgs is the Heavy (and this state mass[0] fixed to 126) and the other is the ligth one (mass[1] < 126) and then another scan where the SM Higgs is the light (and this state mass[1] fixed at 126) and the other is the Heavy one (mass[0] > 126)
// Mixing[i][j] Full mixing matrix (diagonal directions also present in the lines below the mixing directions) i runs over all mass eigenstates and j runs over all initial (gauge) scalar states according to the Mathematica ordering
// r random number generator (called as r(mix,max)) in case the user whishes to change here the scan of MixPars[] more flexibly.
///////////////////////////
// Description:
//The purpose of this function is to provide a way of imposing extra relations in the scan among the last masses and couplings that are left to be generated at the end, through a parametrization, and/or with all the other quantities that have already been generated in in the earlier stages of the scan step.
//IMPORTANT NOTE: it is crucial in this case, that you check in Verbose mode that the relations you impose below are among independent masses and couplings (in VERBOSE mode there are some lines at the end identifying the last independent quantities that are generated). If not, you may have to re-order your fields and couplings in the input mathematica notebook (the rule is that the prigram tries to leave for last, the last masses and couplings in the numbering adopted in the Mathematica notebook). In any case, below you will have to "declare" which couplings and masses you want to constrain so the code will crash and issue a warning if what you are trying to do is not allowed.
// Since at this point in the scan the VEVs, and Mixings have been generated (and possibly their parameters, if they have been re-parametrized), we assume a general structure where the relations you impose may depend on all such (constant) quantities and a new set of parameters (MLparams) whose ranges have been specified in the Mathematica notebook.
//Again, we allow here the use the option to actually change the values of MLparams (i.e. override the values that have been generated in the ranges defined in the mathematica notebook).
//Here you must write down:
// i) (WARNING: this one is NOT optional) specify here which masses and couplings (in the allowed set of independent couplings) you are going to use in your parametrization
IndependentMs.push_back(1);//Mass[1] here for example will re-parametrised
IndependentMs.push_back(2);//and Mass[2] as well
//and couplings (in this example we comment the next line and use only relations among masses)
////IndependentLs.push_back(7);
// ii) Now write here your relations (in this example, this re-parametrization sets two masses to be equal)
///-1) decide first which one is the heavy mixing Higgs & the alpha
doublemA=Mass[2];//Mass of the pseudo scalar (numbering checked in VERBOSE mode)
doublemHcharged=Mass[3];//Mass of the charged Higgs (numbering checked in VERBOSE mode. Note Mass[4] is the same because it corresponds to the imaginary part).
doublemHlight,mHheavy;// Masses of the light and heavy higgs (calculated belwo)
doublealpha;//This is the angle in the 2x2 mixing matrix (calculated below according to the convention that:
// --
// | phi[2]=cos(alpha).H-sin(alpha).h
// | phi[6]=sin(alpha).H+cos(alpha).h
// --
// Note: Since the user has to decide which state (0 or 1) is the heavy and the light, this may imply swaping columns of the mixing matrix. Thus one has to be careful to check that reflections are not introduced. Below for each choice the reflection factor is worked out for each case. Note that this is convention dependent so it must ALWAYS be thought of by the user.
///-1) decide first which one is the heavy mixing Higgs & the alpha
doublemA=150;//Mass of the pseudo scalar (numbering checked in VERBOSE mode)
doublemHcharged=500;//Mass of the charged Higgs (numbering checked in VERBOSE mode. Note Mass[4] is the same because it corresponds to the imaginary part).
doublemHlight=125,mHheavy=500;// Masses of the light and heavy higgs (calculated belwo)
doublealpha=0.5;//This is the angle in the 2x2 mixing matrix