gambit is hosted by Hepforge, IPPP Durham
GAMBIT  v1.5.0-2191-ga4742ac
a Global And Modular Bsm Inference Tool
topness.h File Reference
#include <iostream>
#include <string>
#include <cstdlib>
#include <cmath>
#include "gambit/Utils/threadsafe_rng.hpp"
Include dependency graph for topness.h:
This graph shows which files directly or indirectly include this file:

Go to the source code of this file.

Classes

struct  my_func
 
class  my_simplex
 
class  my_Nelder_Mead
 

Functions

double my_lp (double u[4], double v[4])
 
void my_add (double a[4], double b[4], double c[4])
 
double my_dot (double p1[], double p2[], int d)
 
double my_dot (double p1[], double p2[], double a[], int d)
 
double my_enorm (double p[], int d)
 
double topnesscompute (double pb1[4], double pl[4], double MET[4], double sigmat, double sigmaW)
 

Variables

const double mt =172.
 
const double mW =80.4
 
const int DMAX =4
 

Function Documentation

◆ my_add()

void my_add ( double  a[4],
double  b[4],
double  c[4] 
)

Definition at line 51 of file topness.h.

Referenced by my_func::operator()().

52 {
53  for (int i=0; i<4; i++)
54  {
55  c[i]=a[i]+b[i];
56  }
57 }
START_MODEL b
Definition: demo.hpp:270
Here is the caller graph for this function:

◆ my_dot() [1/2]

double my_dot ( double  p1[],
double  p2[],
int  d 
)

Definition at line 121 of file topness.h.

122 {
123 // compute euclidean dot product between p1 and p2
124  double value=0.;
125  for (int i=0; i<d; i++)
126  {
127  // std::cout << "p1: " << p1[i] <<", p2:" << p2[i] << std::endl;
128  value+=p1[i]*p2[i];
129  }
130  //std::cout << "dp value is " << value << std::endl;
131  return value;
132 }

◆ my_dot() [2/2]

double my_dot ( double  p1[],
double  p2[],
double  a[],
int  d 
)

Definition at line 133 of file topness.h.

134 {
135  // compute euclidean dot product between p1 and p2 with weights a
136 
137  double value=0.;
138  for (int i=0; i<d; i++)
139  {
140  value+=a[i]*p1[i]*p2[i];
141  }
142  return value;
143 }

◆ my_enorm()

double my_enorm ( double  p[],
int  d 
)

Definition at line 144 of file topness.h.

145 {
146  // compute euclidean norm of p
147  double value=0.;
148  double norm=0.;
149  for (int i=0; i<d; i++)
150  {
151  value+=p[i]*p[i];
152  }
153  if (value>0.)
154  {
155  norm=sqrt(value);
156  }
157  else
158  {
159  norm=0.;
160  }
161  return norm;
162 }

◆ my_lp()

double my_lp ( double  u[4],
double  v[4] 
)

Definition at line 41 of file topness.h.

References g.

Referenced by my_func::operator()().

42 { // Lorentz product . assume x, y Lorentz vectors ordered as (px,py,pz,E) etc. Initialize with
43  // metric g==(-1,-1,-1,1)
44 
45  double g[]={-1.,-1.,-1.,1.};
46  double s=0;
47  for (int j=0;j<4;j++)
48  {s+= u[j]*g[j]*v[j]; }
49  return s;
50 }
Here is the caller graph for this function:

◆ topnesscompute()

double topnesscompute ( double  pb1[4],
double  pl[4],
double  MET[4],
double  sigmat,
double  sigmaW 
)

Definition at line 744 of file topness.h.

References alpha, beta, my_Nelder_Mead::d, Gambit::Random::draw(), my_Nelder_Mead::find_global_min(), my_Nelder_Mead::xfinal, and my_Nelder_Mead::yfinal.

745 {
746 
747  double alpha=1.0; // parameter for reflection
748  double beta=0.5; // parameter for contraction
749  double gamma=2.0; // parameter for expansion
750 
751  const int d=4; // number of parameters to scan over, not the space-time dimension!
752  const int DIMMAX=d*(d+1); // dimension of simplex
753  double xin[d+1][d]; // starting point
754  double xstart[DIMMAX]; // algorithm stores d+1 points of simplex in a single array
755 
756 
757  double eps=0.000002; // tolerance
758  double Deltastep=20.; // initial spacing of points, in GeV
759  int Ntry=100000; // maximum number of Nelder-Mead cycles to perform for a given initial seed
760  int Nattempts=15; // number of initial starts
761 
762  double edir[]={1,0,0,0,0,1,0,0,0,0,1,0,0,0,0,1};
763  double ybest=1000000000.; // a big number
764  double ybest1=100000000.;
765  double xbest1[]={1000.,1000.,1000.,1000.};
766 // double ybest2=100000000.; // a big number
767 // double xbest2[]={1000.,1000.,1000.,1000.};
768  int i,j,k;
769  // int tid;
770  bool converge1;
771 // converge2;
772 
773  // initialize topness function ***Modified by Yang Zhang 23.1.2018
774  my_func topstat1(pb1,pl,MET,sigmat,sigmaW,1.0);
775 // my_func topstat2(pb2,pb1,pl,MET,sigmat,sigmaW,1.0); // other combination
776 
777  // initialize topness computation
778  my_Nelder_Mead my_check1(d,alpha,beta,gamma, Ntry, eps, &topstat1);
779 // my_Nelder_Mead my_check2(d,alpha,beta,gamma, Ntry, eps, &topstat2);
780  double yl;
781  // begin loop over Nattempts
782  for (k=0; k<Nattempts; k++)
783  {
784  // initialize starting point
785  for (i=0;i<d+1; i++)
786  {
787  for (j=0; j<d; j++)
788  {
789  if (i==0)
790  {
791  xin[0][j]=8000.0*(Gambit::Random::draw()-0.5);
792  }
793  else
794  {
795  xin[i][j]=xin[0][j]+Deltastep*edir[d*(i-1)+j];
796  }
797  }
798  std::copy(xin[i],xin[i]+d, xstart+d*i); // copy initial data into xstart
799  }
800  // now first combination
801  converge1=my_check1.find_global_min(xstart);
802  if (converge1==true)
803  {
804 
805  yl=my_check1.yfinal;
806  if (yl < ybest1)
807  {
808  ybest1=yl;
809  std::copy(my_check1.xfinal,my_check1.xfinal+d,xbest1);
810  }
811  }
812  else
813  {
814  std::cout << " Minimum not found...exiting " << std::endl;
815  }
816 // // now do second combination
817 // converge2=my_check2.find_global_min(xstart);
818 // if (converge2==true)
819 // {
820 // yl=my_check2.yfinal;
821 // if (yl < ybest2)
822 // {
823 // ybest2=yl;
824 // std::copy(my_check2.xfinal,my_check2.xfinal+d,xbest2);
825 // }
826 // }
827 // else
828 // {
829 // std::cout << " Minimum not found...exiting " << std::endl;
830 // }
831 
832  }
833 
834 // if (ybest1 < ybest2)
835 // {
836  ybest=ybest1;
837 // std::copy(xbest1,xbest1+d,xbest);
838 // }
839 // else
840 // {
841 // ybest=ybest2;
842 // std::copy(xbest2,xbest2+d,xbest);
843 // }
844 
845  return ybest;
846 
847 }
START_MODEL beta
Definition: Axions.hpp:36
START_MODEL alpha
static double draw()
Draw a single uniform random deviate from the interval (0,1) using the chosen RNG engine...
Here is the call graph for this function:

Variable Documentation

◆ DMAX

const int DMAX =4

Definition at line 120 of file topness.h.

◆ mt

◆ mW