gambit is hosted by Hepforge, IPPP Durham
GAMBIT  v1.5.0-252-gf9a3f78
a Global And Modular Bsm Inference Tool
SpecBit_externaltests.hpp
Go to the documentation of this file.
1 // GAMBIT: Global and Modular BSM Inference Tool
2 // *********************************************
23 
24 #ifndef __SpecBit_tests_hpp__
25 #define __SpecBit_tests_hpp__
26 
29 
30 // Flexible SUSY stuff (should not be needed by the rest of gambit)
31 #include "flexiblesusy/src/ew_input.hpp"
32 #include "flexiblesusy/src/numerics2.hpp"
33 #include "flexiblesusy/src/wrappers.hpp"
34 
35 // Switch test output depending on where this is being compiled
36 #ifdef IN_SPECBIT
37  #define OUTPUT logger()
38  #define TAGerr LogTags::err
39  #define TAGfatal LogTags::fatal
40  #define TAGeom EOM
41 #else
42  #define OUTPUT std::cerr
43  #define TAGerr ""
44  #define TAGfatal ""
45  #define TAGeom ""
46 #endif
47 
48 
49 
50 namespace Gambit
51 {
52 
53  namespace SpecBit
54  {
55 
56  bool print_error(bool pass, std::string get_type, std::string data,
57  double sting_get_out, double data_member,
58  int i = -1, int j = -1)
59  {
60  OUTPUT << " returning fail on test for: " << std::endl;
61  if (i > -1) OUTPUT << "with first index = " << i << std::endl;
62  if (j > -1) OUTPUT << "with second index = " << j << std::endl;
63  OUTPUT << get_type << " with " << data << " string arg." <<std::endl;
64  OUTPUT << "string getter gives = "
65  <<sting_get_out << std::endl;
66  OUTPUT << "data member is = "
67  << data_member << std::endl;
68  OUTPUT << TAGerr << TAGfatal << TAGeom;
69  return pass;
70 
71  }
72 
73  void print_error(std::string get_type, std::string name,
74  double sting_get_out, double data_member,
75  int i = -1, int j = -1)
76  {
77  OUTPUT << " returning fail on test for: " << std::endl;
78  if (i > -1) OUTPUT << "with first index = " << i << std::endl;
79  if (j > -1) OUTPUT << "with second index = " << j << std::endl;
80  OUTPUT << get_type << " with " << name << " string arg."
81  <<std::endl;
82  OUTPUT << "string getter gives = " <<sting_get_out << std::endl;
83  OUTPUT << "data member is = "
84  << data_member << std::endl;
85  OUTPUT << TAGerr << TAGfatal << TAGeom;
86  }
87 
88  bool test_getters(std::string get_type,
89  std::string name, double getter_output,
90  double data_member, int i = -1, int j = -1)
91  {
92  bool pass = flexiblesusy::is_equal(getter_output,data_member);
93  // if(pass == false) return pass;
94  if(pass == false)
95  print_error(get_type, name, getter_output, data_member,i,j);
96  return pass;
97  }
98 
99 
101  // These are not known to Gambit
102 
103  template <class M>
104  bool TestMssmParMass2_0(SubSpectrum * spec, M FSmssm,
105  bool immediate_exit = true)
106  {
107  bool pass = false;
108  //do all in loop
109 
110  // Ben: No initializer lists allowed; changing all these.
111  //std::set<std::pair<std::string,double>> name_value = {
112  // {"BMu", FSmssm.get_BMu()},
113  // {"mHd2", FSmssm.get_mHd2()},
114  // {"mHu2", FSmssm.get_mHu2()}
115  //};
116 
117  std::set<std::pair<std::string,double>> name_value;
118  name_value.insert(std::make_pair( "BMu" , FSmssm.get_BMu() ));
119  name_value.insert(std::make_pair( "mHd2", FSmssm.get_mHd2() ));
120  name_value.insert(std::make_pair( "mHu2", FSmssm.get_mHu2() ));
121 
122  std::set<std::pair<std::string, double>>::iterator iter;
123  for(iter=name_value.begin(); iter != name_value.end(); ++iter)
124  {
125  pass = test_getters("get_mass2_parameter", iter->first,
126  spec->
127  get_mass2_parameter(iter->first),
128  iter->second);
129  if(immediate_exit == true && pass == false) return pass;
130  }
131 
132  return pass;
133  }
134 
135 
136  template <class MI>
137  bool TestMssmParMass2_0(MSSMSpec<MI>& mssm, typename MI::Model& FSmssm,
138  bool immediate_exit = true)
139  {
140  bool pass = false;
141  //do all in loop
142  std::set<std::pair<std::string,double>> name_value;
143  name_value.insert(std::make_pair( "BMu" , FSmssm.get_BMu() ));
144  name_value.insert(std::make_pair( "mHd2", FSmssm.get_mHd2() ));
145  name_value.insert(std::make_pair( "mHu2", FSmssm.get_mHu2() ));
146 
147  std::set<std::pair<std::string, double>>::iterator iter;
148  for(iter=name_value.begin(); iter != name_value.end(); ++iter)
149  {
150  pass = test_getters("get_mass2_parameter", iter->first,
151  mssm.
152  get_mass2_parameter(iter->first),
153  iter->second);
154  if(immediate_exit == true && pass == false) return pass;
155  }
156 
157  return pass;
158  }
159 
160 
161 
162 
163  template <class M>
164  bool TestMssmParMass2_2(SubSpectrum * spec, M FSmssm,
165  bool immediate_exit=true)
166  {
167  bool pass = false;
168 
169  for(int i=1; i<=3; i++){
170  for(int j=1; j<=3; j++){
171  //Would be smarter to take these out of for loop and
172  // use function pointers, but I won't.
173  std::set<std::pair<std::string,double>> name_value;
174  name_value.insert(std::make_pair( "mq2", FSmssm.get_mq2(i-1,j-1) ));
175  name_value.insert(std::make_pair( "mu2", FSmssm.get_mu2(i-1,j-1) ));
176  name_value.insert(std::make_pair( "md2", FSmssm.get_md2(i-1,j-1) ));
177  name_value.insert(std::make_pair( "ml2", FSmssm.get_ml2(i-1,j-1) ));
178  name_value.insert(std::make_pair( "me2", FSmssm.get_me2(i-1,j-1) ));
179 
180  std::set<std::pair<std::string, double>>::iterator iter;
181  for(iter=name_value.begin(); iter != name_value.end(); ++iter)
182  {
183  pass = test_getters("get_mass2_parameter", iter->first,
184  spec->
185  get_mass2_parameter(iter->first,i,j),
186  iter->second, i, j);
187  if(immediate_exit == true && pass == false) return pass;
188  }
189  }
190  }
191  return pass;
192  }
193 
194 
195  template <class MI>
196  bool TestMssmParMass2_2(MSSMSpec<MI>& mssm, typename MI::Model& FSmssm,
197  bool immediate_exit=true)
198  {
199  bool pass = false;
200 
201  for(int i=1; i<=3; i++){
202  for(int j=1; j<=3; j++){
203  //Would be smarter to take these out of for loop and
204  // use function pointers, but I won't.
205  std::set<std::pair<std::string,double>> name_value;
206  name_value.insert(std::make_pair( "mq2", FSmssm.get_mq2(i-1,j-1) ));
207  name_value.insert(std::make_pair( "mu2", FSmssm.get_mu2(i-1,j-1) ));
208  name_value.insert(std::make_pair( "md2", FSmssm.get_md2(i-1,j-1) ));
209  name_value.insert(std::make_pair( "ml2", FSmssm.get_ml2(i-1,j-1) ));
210  name_value.insert(std::make_pair( "me2", FSmssm.get_me2(i-1,j-1) ));
211 
212  std::set<std::pair<std::string, double>>::iterator iter;
213  for(iter=name_value.begin(); iter != name_value.end(); ++iter)
214  {
215  pass = test_getters("get_mass2_parameter", iter->first,
216  mssm.
217  get_mass2_parameter(iter->first,i,j),
218  iter->second, i, j);
219  if(immediate_exit == true && pass == false) return pass;
220  }
221  }
222  }
223  return pass;
224  }
225 
226 
227  template <class M>
228  bool TestMssmParMass1_0(SubSpectrum * spec, M FSmssm,
229  bool immediate_exit=true)
230  {
231  bool pass = false;
232 
233  std::set<std::pair<std::string,double>> name_value;
234  name_value.insert(std::make_pair( "M1", FSmssm.get_MassB() ));
235  name_value.insert(std::make_pair( "M2", FSmssm.get_MassWB() ));
236  name_value.insert(std::make_pair( "M3", FSmssm.get_MassG() ));
237  name_value.insert(std::make_pair( "vu", FSmssm.get_vu() ));
238  name_value.insert(std::make_pair( "vd", FSmssm.get_vd() ));
239 
240  std::set<std::pair<std::string, double>>::iterator iter;
241  for(iter=name_value.begin(); iter != name_value.end(); ++iter)
242  {
243  pass = test_getters("get_mass_parameter", iter->first,
244  spec->
245  get_mass_parameter(iter->first),
246  iter->second);
247  if(immediate_exit == true && pass == false) return pass;
248  }
249 
250  return pass;
251  }
252 
253  template <class MI>
254  bool TestMssmParMass1_0(MSSMSpec<MI>& mssm, typename MI::Model& FSmssm,
255  bool immediate_exit=true)
256  {
257  bool pass = false;
258 
259  std::set<std::pair<std::string,double>> name_value;
260  name_value.insert(std::make_pair( "M1", FSmssm.get_MassB() ));
261  name_value.insert(std::make_pair( "M2", FSmssm.get_MassWB() ));
262  name_value.insert(std::make_pair( "M3", FSmssm.get_MassG() ));
263  name_value.insert(std::make_pair( "vu", FSmssm.get_vu() ));
264  name_value.insert(std::make_pair( "vd", FSmssm.get_vd() ));
265 
266  std::set<std::pair<std::string, double>>::iterator iter;
267  for(iter=name_value.begin(); iter != name_value.end(); ++iter)
268  {
269  pass = test_getters("get_mass_parameter", iter->first,
270  mssm.
271  get_mass_parameter(iter->first),
272  iter->second);
273  if(immediate_exit == true && pass == false) return pass;
274  }
275 
276  return pass;
277  }
278 
279 
280 
281  template <class M>
282  bool TestMssmParMass1_2(SubSpectrum * spec, M FSmssm,
283  bool immediate_exit =true)
284  {
285  bool pass = false;
286  for(int i=1; i<=3; i++){
287  for(int j=1; j<=3; j++){
288  std::set<std::pair<std::string,double>> name_value;
289  name_value.insert(std::make_pair( "TYd", FSmssm.get_TYd(i-1,j-1) ));
290  name_value.insert(std::make_pair( "TYu", FSmssm.get_TYu(i-1,j-1) ));
291  name_value.insert(std::make_pair( "TYe", FSmssm.get_TYe(i-1,j-1) ));
292 
293  std::set<std::pair<std::string, double>>::iterator iter;
294  for(iter=name_value.begin(); iter != name_value.end(); ++iter)
295  {
296  pass = test_getters("get_mass_parameter", iter->first,
297  spec->
298  get_mass_parameter(iter->first,i,j),
299  iter->second, i, j);
300  if(immediate_exit == true && pass == false) return pass;
301  }
302  }
303  }
304  return pass;
305  }
306 
307 
308  template <class MI>
309  bool TestMssmParMass1_2(MSSMSpec<MI>& mssm, typename MI::Model& FSmssm,
310  bool immediate_exit =true)
311  {
312  bool pass = false;
313  for(int i=1; i<=3; i++){
314  for(int j=1; j<=3; j++){
315  std::set<std::pair<std::string,double>> name_value;
316  name_value.insert(std::make_pair( "TYd", FSmssm.get_TYd(i-1,j-1) ));
317  name_value.insert(std::make_pair( "TYu", FSmssm.get_TYu(i-1,j-1) ));
318  name_value.insert(std::make_pair( "TYe", FSmssm.get_TYe(i-1,j-1) ));
319 
320  std::set<std::pair<std::string, double>>::iterator iter;
321  for(iter=name_value.begin(); iter != name_value.end(); ++iter)
322  {
323  pass = test_getters("get_mass_parameter", iter->first,
324  mssm.
325  get_mass_parameter(iter->first,i,j),
326  iter->second, i, j);
327  if(immediate_exit == true && pass == false) return pass;
328  }
329  }
330  }
331  return pass;
332  }
333 
334  template <class M>
335  bool TestMssmParMass0_0(SubSpectrum * spec, M FSmssm,
336  bool immediate_exit =true )
337  {
338  bool pass = false;
339  std::set<std::pair<std::string,double>> name_value;
340  name_value.insert(std::make_pair( "g1", FSmssm.get_g1() ));
341  name_value.insert(std::make_pair( "g2", FSmssm.get_g2() ));
342  name_value.insert(std::make_pair( "g3", FSmssm.get_g3() ));
343 
344  std::set<std::pair<std::string, double>>::iterator iter;
345  for(iter=name_value.begin(); iter != name_value.end(); ++iter)
346  {
347  pass = test_getters("get_dimensionless_parameter", iter->first,
348  spec->
349  get_dimensionless_parameter(iter->first),
350  iter->second);
351  if(immediate_exit == true && pass == false) return pass;
352  }
353 
354  return pass;
355  }
356 
357 
358  template <class MI>
359  bool TestMssmParMass0_0(MSSMSpec<MI>& mssm, typename MI::Model& FSmssm,
360  bool immediate_exit =true )
361  {
362  bool pass = false;
363  std::set<std::pair<std::string,double>> name_value;
364  name_value.insert(std::make_pair( "g1", FSmssm.get_g1() ));
365  name_value.insert(std::make_pair( "g2", FSmssm.get_g2() ));
366  name_value.insert(std::make_pair( "g3", FSmssm.get_g3() ));
367 
368  std::set<std::pair<std::string, double>>::iterator iter;
369  for(iter=name_value.begin(); iter != name_value.end(); ++iter)
370  {
371  pass = test_getters("get_dimensionless_parameter", iter->first,
372  mssm.
373  get_dimensionless_parameter(iter->first),
374  iter->second);
375  if(immediate_exit == true && pass == false) return pass;
376  }
377 
378  return pass;
379  }
380 
381  template <class M>
382  bool TestMssmParMass0_2(SubSpectrum * spec, M FSmssm,
383  bool immediate_exit = true)
384  {
385  bool pass = false;
386  for(int i=1; i<=3; i++){
387  for(int j=1; j<=3; j++){
388 
389  std::set<std::pair<std::string,double>> name_value;
390  name_value.insert(std::make_pair( "Yd", FSmssm.get_Yd(i-1,j-1) ));
391  name_value.insert(std::make_pair( "Yu", FSmssm.get_Yu(i-1,j-1) ));
392  name_value.insert(std::make_pair( "Ye", FSmssm.get_Ye(i-1,j-1) ));
393 
394  std::set<std::pair<std::string, double>>::iterator iter;
395  for(iter=name_value.begin(); iter != name_value.end(); ++iter)
396  {
397  pass = test_getters("get_dimensionless_parameter",
398  iter->first,
399  spec->
400  get_dimensionless_parameter(iter->first,
401  i,j),
402  iter->second, i, j);
403  if(immediate_exit == true && pass == false) return pass;
404  }
405  }
406  }
407  return pass;
408  }
409 
410 
411  template <class MI>
412  bool TestMssmParMass0_2(MSSMSpec<MI>& mssm, typename MI::Model& FSmssm,
413  bool immediate_exit = true)
414  {
415  bool pass = false;
416  for(int i=1; i<=3; i++){
417  for(int j=1; j<=3; j++){
418 
419  std::set<std::pair<std::string,double>> name_value;
420  name_value.insert(std::make_pair( "Yd", FSmssm.get_Yd(i-1,j-1) ));
421  name_value.insert(std::make_pair( "Yu", FSmssm.get_Yu(i-1,j-1) ));
422  name_value.insert(std::make_pair( "Ye", FSmssm.get_Ye(i-1,j-1) ));
423 
424  std::set<std::pair<std::string, double>>::iterator iter;
425  for(iter=name_value.begin(); iter != name_value.end(); ++iter)
426  {
427  pass = test_getters("get_dimensionless_parameter",
428  iter->first,
429  mssm.
430  get_dimensionless_parameter(iter->first,
431  i,j),
432  iter->second, i, j);
433  if(immediate_exit == true && pass == false) return pass;
434  }
435  }
436  }
437  return pass;
438  }
439 
440  template <class M>
441  bool TestMssmPoleGets0(SubSpectrum * spec, M FSmssm,
442  bool immediate_exit = true)
443  {
444  bool pass = false;
445  //do all in loop
446  std::set<std::pair<std::string,double>> name_value;
447  name_value.insert(std::make_pair( "MZ", FSmssm.get_physical().MVZ ));
448  name_value.insert(std::make_pair( "MW", FSmssm.get_physical().MVWm ));
449  name_value.insert(std::make_pair( "MGluino", FSmssm.get_physical().MGlu ));
450  name_value.insert(std::make_pair( "MGluon", FSmssm.get_physical().MVG ));
451  name_value.insert(std::make_pair( "MPhoton", FSmssm.get_physical().MVP ));
452 
453  std::set<std::pair<std::string, double>>::iterator iter;
454  for(iter=name_value.begin(); iter != name_value.end(); ++iter)
455  {
456  pass = test_getters("get_Pole_Mass", iter->first,
457  spec->get_Pole_Mass(iter->first),
458  iter->second);
459  if(immediate_exit == true && pass == false) return pass;
460  }
461  return pass;
462  }
463 
464 
465  template <class MI>
466  bool TestMssmPoleGets0(MSSMSpec<MI>& mssm, typename MI::Model& FSmssm,
467  bool immediate_exit = true)
468  {
469  bool pass = false;
470  //do all in loop
471  std::set<std::pair<std::string,double>> name_value;
472  name_value.insert(std::make_pair( "MZ", FSmssm.get_physical().MVZ ));
473  name_value.insert(std::make_pair( "MW", FSmssm.get_physical().MVWm ));
474  name_value.insert(std::make_pair( "MGluino", FSmssm.get_physical().MGlu ));
475  name_value.insert(std::make_pair( "MGluon", FSmssm.get_physical().MVG ));
476  name_value.insert(std::make_pair( "MPhoton", FSmssm.get_physical().MVP ));
477 
478  std::set<std::pair<std::string, double>>::iterator iter;
479  for(iter=name_value.begin(); iter != name_value.end(); ++iter)
480  {
481  pass = test_getters("get_Pole_Mass", iter->first,
482  mssm.get_Pole_Mass(iter->first),
483  iter->second);
484  if(immediate_exit == true && pass == false) return pass;
485  }
486  return pass;
487  }
488 
489 
490  template <class M>
491  bool TestMssmPoleGets1(SubSpectrum * spec, M FSmssm,
492  bool immediate_exit = true)
493  {
494  bool pass = false;
495 
496  for(int i=1; i<=6; i++){
497  std::set<std::pair<std::string,double>> name_value;
498  name_value.insert(std::make_pair( "MSd", FSmssm.get_physical().MSd(i-1) ));
499  name_value.insert(std::make_pair( "MSu", FSmssm.get_physical().MSu(i-1) ));
500  name_value.insert(std::make_pair( "MSe", FSmssm.get_physical().MSe(i-1) ));
501 
502 
503  std::set<std::pair<std::string, double>>::iterator iter;
504  for(iter=name_value.begin(); iter != name_value.end(); ++iter)
505  {
506  pass = test_getters("get_Pole_Mass", iter->first,
507  spec->get_Pole_Mass(iter->first,i),
508  iter->second, i);
509  if(immediate_exit == true && pass == false) return pass;
510  }
511 
512  }
513  for(int i=1; i<=3; i++){
514 
515  std::set<std::pair<std::string,double>> name_value;
516  name_value.insert(std::make_pair( "MSv", FSmssm.get_physical().MSv(i-1) ));
517  name_value.insert(std::make_pair( "MFd", FSmssm.get_physical().MFd(i-1) ));
518  name_value.insert(std::make_pair( "MFu", FSmssm.get_physical().MFu(i-1) ));
519  name_value.insert(std::make_pair( "MFe", FSmssm.get_physical().MFe(i-1) ));
520 
521  std::set<std::pair<std::string, double>>::iterator iter;
522  for(iter=name_value.begin(); iter != name_value.end(); ++iter)
523  {
524  pass = test_getters("get_Pole_Mass", iter->first,
525  spec->get_Pole_Mass(iter->first,i),
526  iter->second, i);
527  if(immediate_exit == true && pass == false) return pass;
528  }
529 
530  }
531  for(int i=1; i<=2; i++){
532  std::string name = "Mh0";
533  pass = test_getters("get_Pole_Mass", name,
534  spec->get_Pole_Mass(name,i),
535  FSmssm.get_physical().Mhh(i-1), i);
536  if(immediate_exit == true && pass == false) return pass;
537  name = "MCha";
538  pass = test_getters("get_Pole_Mass", name,
539  spec->get_Pole_Mass(name,i),
540  FSmssm.get_physical_slha().MCha(i-1), i);
541  if(immediate_exit == true && pass == false) return pass;
542  }
543  //In the the neutralino and chargino tests I compare against
544  // value in physical_slha struct since the value in
545  // physical may differ by a sign since it stores positive masses
546  // and a complex mixing matrix.
547  for(int i=1; i<=4; i++){
548  std::string name = "MChi";
549  pass = test_getters("get_Pole_Mass", name,
550  spec->get_Pole_Mass(name,i),
551  FSmssm.get_physical_slha().MChi(i-1), i);
552  if(immediate_exit == true && pass == false) return pass;
553  }
554  return pass;
555  }
556 
557 
558  template <class MI>
559  bool TestMssmPoleGets1(MSSMSpec<MI>& mssm, typename MI::Model& FSmssm,
560  bool immediate_exit = true)
561  {
562  bool pass = false;
563 
564  for(int i=1; i<=6; i++){
565  std::set<std::pair<std::string,double>> name_value;
566  name_value.insert(std::make_pair( "MSd", FSmssm.get_physical().MSd(i-1) ));
567  name_value.insert(std::make_pair( "MSu", FSmssm.get_physical().MSu(i-1) ));
568  name_value.insert(std::make_pair( "MSe", FSmssm.get_physical().MSe(i-1) ));
569 
570  std::set<std::pair<std::string, double>>::iterator iter;
571  for(iter=name_value.begin(); iter != name_value.end(); ++iter)
572  {
573  pass = test_getters("get_Pole_Mass", iter->first,
574  mssm.get_Pole_Mass(iter->first,i),
575  iter->second, i);
576  if(immediate_exit == true && pass == false) return pass;
577  }
578 
579  }
580  for(int i=1; i<=3; i++){
581 
582  std::set<std::pair<std::string,double>> name_value;
583  name_value.insert(std::make_pair( "MSv", FSmssm.get_physical().MSv(i-1) ));
584  name_value.insert(std::make_pair( "MFd", FSmssm.get_physical().MFd(i-1) ));
585  name_value.insert(std::make_pair( "MFu", FSmssm.get_physical().MFu(i-1) ));
586  name_value.insert(std::make_pair( "MFe", FSmssm.get_physical().MFe(i-1) ));
587 
588  std::set<std::pair<std::string, double>>::iterator iter;
589  for(iter=name_value.begin(); iter != name_value.end(); ++iter)
590  {
591  pass = test_getters("get_Pole_Mass", iter->first,
592  mssm.get_Pole_Mass(iter->first,i),
593  iter->second, i);
594  if(immediate_exit == true && pass == false) return pass;
595  }
596 
597  }
598  for(int i=1; i<=2; i++){
599  std::string name = "Mh0";
600  pass = test_getters("get_Pole_Mass", name,
601  mssm.get_Pole_Mass(name,i),
602  FSmssm.get_physical().Mhh(i-1), i);
603  if(immediate_exit == true && pass == false) return pass;
604  name = "MCha";
605  pass = test_getters("get_Pole_Mass", name,
606  mssm.get_Pole_Mass(name,i),
607  FSmssm.get_physical_slha().MCha(i-1), i);
608  if(immediate_exit == true && pass == false) return pass;
609  }
610  //In the the neutralino and chargino tests I compare against
611  // value in physical_slha struct since the value in
612  // physical may differ by a sign since it stores positive masses
613  // and a complex mixing matrix.
614  for(int i=1; i<=4; i++){
615  std::string name = "MChi";
616  pass = test_getters("get_Pole_Mass", name,
617  mssm.get_Pole_Mass(name,i),
618  FSmssm.get_physical_slha().MChi(i-1), i);
619  if(immediate_exit == true && pass == false) return pass;
620  }
621  return pass;
622  }
623 
624 
625  template <class M>
626  bool TestMssmPoleMixingGets2(SubSpectrum * spec, M FSmssm,
627  bool immediate_exit = true)
628  {
629  bool pass = false;
630  for(int i=1; i<=6; i++){
631  for(int j=1; j<=6; j++){
632  std::set<std::pair<std::string,double>> name_value;
633  name_value.insert(std::make_pair( "~d", FSmssm.get_physical_slha().ZD(i-1,j-1) ));
634  name_value.insert(std::make_pair( "~u", FSmssm.get_physical_slha().ZU(i-1,j-1) ));
635  name_value.insert(std::make_pair( "~e-", FSmssm.get_physical_slha().ZE(i-1,j-1) ));
636 
637  std::set<std::pair<std::string, double>>::iterator iter;
638  for(iter=name_value.begin(); iter != name_value.end(); ++iter)
639  {
640  pass = test_getters("get_Pole_Mixing", iter->first,
641  spec->
642  get_Pole_Mixing(iter->first,i,j),
643  iter->second, i, j);
644  if(immediate_exit == true && pass == false) return pass;
645  }
646  }
647  }
648 
649 
650  for(int i=1; i<=3; i++){
651  for(int j=1; j<=3; j++){
652  string name = "~nu";
653  pass = test_getters("get_Pole_Mixing", name,
654  spec->get_Pole_Mixing(name,i,j),
655  FSmssm.get_physical_slha().ZV(i-1, j-1), i,j);
656  if(immediate_exit == true && pass == false) return pass;
657 
658  }
659  }
660 
661 
662  for(int i=1; i<=2; i++){
663  for(int j=1; j<=2; j++){
664  std::set<std::pair<std::string,double>> name_value;
665  name_value.insert(std::make_pair( "h0", FSmssm.get_physical_slha().ZH(i-1,j-1) ));
666  name_value.insert(std::make_pair( "A0", FSmssm.get_physical_slha().ZA(i-1,j-1) ));
667  name_value.insert(std::make_pair( "H+", FSmssm.get_physical_slha().ZP(i-1,j-1) ));
668  name_value.insert(std::make_pair( "~chi-", flexiblesusy::Re(FSmssm.get_physical_slha()
669  .UM(i-1,j-1)) ));
670  name_value.insert(std::make_pair( "~chi+", flexiblesusy::Re(FSmssm.get_physical_slha()
671  .UP(i-1,j-1)) ));
672 
673 
674  std::set<std::pair<std::string, double>>::iterator iter;
675  for(iter=name_value.begin(); iter != name_value.end(); ++iter)
676  {
677  pass = test_getters("get_Pole_Mixing", iter->first,
678  spec->
679  get_Pole_Mixing(iter->first,i,j),
680  iter->second, i, j);
681  if(immediate_exit == true && pass == false) return pass;
682  }
683  }
684  }
685 
686  return pass;
687  }
688 
689 
690  template <class MI>
691  bool TestMssmPoleMixingGets2(MSSMSpec<MI>& mssm, typename MI::Model& FSmssm,
692  bool immediate_exit = true)
693  {
694  bool pass = false;
695  for(int i=1; i<=6; i++){
696  for(int j=1; j<=6; j++){
697  std::set<std::pair<std::string,double>> name_value;
698  name_value.insert(std::make_pair( "~d", FSmssm.get_physical_slha().ZD(i-1,j-1) ));
699  name_value.insert(std::make_pair( "~u", FSmssm.get_physical_slha().ZU(i-1,j-1) ));
700  name_value.insert(std::make_pair( "~e-", FSmssm.get_physical_slha().ZE(i-1,j-1) ));
701 
702  std::set<std::pair<std::string, double>>::iterator iter;
703  for(iter=name_value.begin(); iter != name_value.end(); ++iter)
704  {
705  pass = test_getters("get_Pole_Mixing", iter->first,
706  mssm.
707  get_Pole_Mixing(iter->first,i,j),
708  iter->second, i, j);
709  if(immediate_exit == true && pass == false) return pass;
710  }
711  }
712  }
713 
714 
715  for(int i=1; i<=3; i++){
716  for(int j=1; j<=3; j++){
717  string name = "~nu";
718  pass = test_getters("get_Pole_Mixing", name,
719  mssm.get_Pole_Mixing(name,i,j),
720  FSmssm.get_physical_slha().ZV(i-1, j-1), i,j);
721  if(immediate_exit == true && pass == false) return pass;
722 
723  }
724  }
725 
726 
727  for(int i=1; i<=2; i++){
728  for(int j=1; j<=2; j++){
729  std::set<std::pair<std::string,double>> name_value;
730  name_value.insert(std::make_pair( "h0", FSmssm.get_physical_slha().ZH(i-1,j-1) ));
731  name_value.insert(std::make_pair( "A0", FSmssm.get_physical_slha().ZA(i-1,j-1) ));
732  name_value.insert(std::make_pair( "H+", FSmssm.get_physical_slha().ZP(i-1,j-1) ));
733  name_value.insert(std::make_pair( "~chi-", flexiblesusy::Re(FSmssm.get_physical_slha()
734  .UM(i-1,j-1)) ));
735  name_value.insert(std::make_pair( "~chi+", flexiblesusy::Re(FSmssm.get_physical_slha()
736  .UP(i-1,j-1)) ));
737 
738  std::set<std::pair<std::string, double>>::iterator iter;
739  for(iter=name_value.begin(); iter != name_value.end(); ++iter)
740  {
741  pass = test_getters("get_Pole_Mixing", iter->first,
742  mssm.
743  get_Pole_Mixing(iter->first,i,j),
744  iter->second, i, j);
745  if(immediate_exit == true && pass == false) return pass;
746  }
747  }
748  }
749 
750  return pass;
751  }
752 
753 
754  template <class M>
755  bool TestMssmPoleGets(SubSpectrum * spec, M FSmssm)
756  {
757  bool pass = false;
758  pass = TestMssmPoleGets0(spec,FSmssm);
759  if(pass == false) return pass;
760  pass = TestMssmPoleGets1(spec,FSmssm);
761  if(pass == false) return pass;
762  pass = TestMssmPoleMixingGets2(spec,FSmssm);
763  if(pass == false) return pass;
764  return pass;
765  }
766 
767  template <class MI>
768  bool TestMssmPoleGets(MSSMSpec<MI>& mssm, typename MI::Model& FSmssm)
769  {
770  bool pass = false;
771  pass = TestMssmPoleGets0(mssm,FSmssm);
772  if(pass == false) return pass;
773  pass = TestMssmPoleGets1(mssm,FSmssm);
774  if(pass == false) return pass;
775  pass = TestMssmPoleMixingGets2(mssm,FSmssm);
776  if(pass == false) return pass;
777  return pass;
778  }
779 
780  template <class M>
781  bool TestMssmParGets(SubSpectrum * spec, M FSmssm)
782  {
783  bool pass = false;
784  pass = TestMssmParMass2_0(spec,FSmssm);
785  if(pass == false) return pass;
786  pass = TestMssmParMass2_2(spec,FSmssm);
787  if(pass == false) return pass;
788  pass = TestMssmParMass1_0(spec,FSmssm);
789  if(pass == false) return pass;
790  pass = TestMssmParMass1_2(spec,FSmssm);
791  if(pass == false) return pass;
792  pass = TestMssmParMass0_0(spec,FSmssm);
793  if(pass == false) return pass;
794  pass = TestMssmParMass0_2(spec,FSmssm);
795  if(pass == false) return pass;
796 
797  return pass;
798 
799  }
800 
801  template <class MI>
802  bool TestMssmParGets(MSSMSpec<MI>& mssm, typename MI::Model& FSmssm)
803  {
804  bool pass = false;
805  pass = TestMssmParMass2_0(mssm,FSmssm);
806  if(pass == false) return pass;
807  pass = TestMssmParMass2_2(mssm,FSmssm);
808  if(pass == false) return pass;
809  pass = TestMssmParMass1_0(mssm,FSmssm);
810  if(pass == false) return pass;
811  pass = TestMssmParMass1_2(mssm,FSmssm);
812  if(pass == false) return pass;
813  pass = TestMssmParMass0_0(mssm,FSmssm);
814  if(pass == false) return pass;
815  pass = TestMssmParMass0_2(mssm,FSmssm);
816  if(pass == false) return pass;
817 
818  return pass;
819 
820  }
821 
822  template <class Model>
823  void setup(Model& mssm)
824  {
825  Eigen::Matrix<double,3,3> Yu;
826  Eigen::Matrix<double,3,3> Yd;
827  Eigen::Matrix<double,3,3> Ye;
828  double Mu;
829  double g1;
830  double g2;
831  double g3;
832  double vd;
833  double vu;
834  Eigen::Matrix<double,3,3> TYu;
835  Eigen::Matrix<double,3,3> TYd;
836  Eigen::Matrix<double,3,3> TYe;
837  double BMu;
838  Eigen::Matrix<double,3,3> mq2;
839  Eigen::Matrix<double,3,3> ml2;
840  double mHd2;
841  double mHu2;
842  Eigen::Matrix<double,3,3> md2;
843  Eigen::Matrix<double,3,3> mu2;
844  Eigen::Matrix<double,3,3> me2;
845  double MassB;
846  double MassWB;
847  double MassG;
848 
849  // susy parameters
850  Yu << 1.26136e-05, 0, 0,
851  0, 0.00667469, 0,
852  0, 0, 0.857849;
853 
854  Yd << 0.000242026, 0, 0,
855  0, 0.00529911, 0,
856  0, 0, 0.193602;
857 
858  Ye << 2.84161e-05, 0, 0,
859  0, 0.00587557, 0,
860  0, 0, 0.10199;
861 
862  Mu = 627.164;
863  g1 = 0.468171;
864  g2 = 0.642353;
865  g3 = 1.06459;
866  vd = 25.0944;
867  vu = 242.968;
868 
869  // soft parameters
870  TYu << -0.0144387, 0, 0,
871  0, -7.64037, 0,
872  0, 0, -759.305;
873 
874  TYd << -0.336207, 0, 0,
875  0, -7.36109, 0,
876  0, 0, -250.124;
877 
878  TYe << -0.00825134, 0, 0,
879  0, -1.70609, 0,
880  0, 0, -29.4466;
881 
882  BMu = 52140.8;
883 
884  mq2 << 1.03883e+06, 0, 0,
885  0, 1.03881e+06, 0,
886  0, 0, 879135;
887 
888  ml2 << 124856, 0, 0,
889  0, 124853, 0,
890  0, 0, 124142;
891 
892  mHd2 = 92436.9;
893  mHu2 = -380337;
894 
895  md2 << 954454, 0, 0,
896  0, 954439, 0,
897  0, 0, 934727;
898 
899  mu2 << 963422, 0, 0,
900  0, 963400, 0,
901  0, 0, 656621;
902 
903  me2 << 49215.8, 0, 0,
904  0, 49210.9, 0,
905  0, 0, 47759.2;
906 
907  MassB = 210.328;
908  MassWB = 389.189;
909  MassG = 1114.45;
910 
911  // set parameters
912  mssm.set_scale(flexiblesusy::Electroweak_constants::MZ);
913  mssm.set_Yu(Yu);
914  mssm.set_Yd(Yd);
915  mssm.set_Ye(Ye);
916  mssm.set_Mu(Mu);
917  mssm.set_g1(g1);
918  mssm.set_g2(g2);
919  mssm.set_g3(g3);
920  mssm.set_vd(vd);
921  mssm.set_vu(vu);
922  mssm.set_TYu(TYu);
923  mssm.set_TYd(TYd);
924  mssm.set_TYe(TYe);
925  mssm.set_BMu(BMu);
926  mssm.set_mq2(mq2);
927  mssm.set_ml2(ml2);
928  mssm.set_mHd2(mHd2);
929  mssm.set_mHu2(mHu2);
930  mssm.set_md2(md2);
931  mssm.set_mu2(mu2);
932  mssm.set_me2(me2);
933  mssm.set_MassB(MassB);
934  mssm.set_MassWB(MassWB);
935  mssm.set_MassG(MassG);
936  }
937 
938 
939  template <class MI>
940  bool test_exact(MSSMSpec<MI>& mssm, typename MI::Model& FS_model_slha)
941  {
942  bool pass = TestMssmParGets(mssm,FS_model_slha);
943  if(pass == false)
944  {
945  OUTPUT << "TestMssmParGets failing." <<std::endl;
946  return pass;
947  }
948  pass = TestMssmPoleGets(mssm,FS_model_slha);
949  if(pass == false)
950  {
951  OUTPUT << "TestMssmParGets failing." <<std::endl;
952  return pass;
953  }
954 
955 
956  return pass;
957  }
958 
959  template <class M>
960  double test_exact(SubSpectrum * spec, M FS_model_slha)
961  {
962  bool pass = TestMssmParGets(spec,FS_model_slha);
963  if(pass == false)
964  {
965  OUTPUT << "TestMssmParGets failing." <<std::endl;
966  return pass;
967  }
968  pass = TestMssmPoleGets(spec,FS_model_slha);
969  if(pass == false)
970  {
971  OUTPUT << "TestMssmParGets failing." <<std::endl;
972  return pass;
973  }
974 
975  return pass;
976 
977  }
978  //This gives identical results after running up, so don't need messy
979  // Test_close
980  template <class MI>
981  bool running_test(MSSMSpec<MI> & mssm, typename MI::Model & FS_model_slha, double tol)
982  {
983  double highscale = 1e+16;
984  double lowscale = mssm.GetScale();
985  double lowscale2 = FS_model_slha.get_scale();
986  bool pass = flexiblesusy::is_equal(lowscale,lowscale2);
987  if(!pass) {
988  OUTPUT << "test fail: "
989  << "objects not at same scale at start of runtest."
990  << std::endl;
991  return pass;
992  }
993 
994  mssm.RunToScale(highscale);
995  FS_model_slha.run_to(highscale);
996  pass = test_exact(mssm, FS_model_slha);
997  if(!pass) {
998  OUTPUT << "test fail: "
999  << "objects not the same after running to MGUT."
1000  << std::endl;
1001  return pass;
1002  }
1003  mssm.RunToScale(lowscale);
1004  FS_model_slha.run_to(lowscale);
1005  pass = test_exact(mssm, FS_model_slha);
1006  if(!pass) {
1007  OUTPUT << "test fail: "
1008  << "objects not the same after running to lowscale."
1009  << std::endl;
1010  return pass;
1011  }
1012 
1013 
1014  return pass;
1015  }
1016 
1017  template <class Model>
1018  bool running_test(SubSpectrum * spec, Model & FS_model_slha,
1019  double tol)
1020  {
1021  double highscale = 1e+16;
1022  double lowscale = spec->GetScale();
1023  double lowscale2 = FS_model_slha.get_scale();
1024  bool pass = flexiblesusy::is_equal(lowscale,lowscale2);
1025  if(!pass) {
1026  OUTPUT << "test fail: "
1027  << "objects not at same scale at start of runtest."
1028  << std::endl;
1029  return pass;
1030  }
1031 
1032  spec->RunToScale(highscale);
1033  FS_model_slha.run_to(highscale);
1034  pass = test_exact(spec, FS_model_slha);
1035  if(!pass) {
1036  OUTPUT << "test fail: "
1037  << "objects not the same after running to MGUT."
1038  << std::endl;
1039  return pass;
1040  }
1041  spec->RunToScale(lowscale);
1042  FS_model_slha.run_to(lowscale);
1043  pass = test_exact(spec, FS_model_slha);
1044  if(!pass) {
1045  OUTPUT << "test fail: "
1046  << "objects not the same after running to lowscale."
1047  << std::endl;
1048  return pass;
1049  }
1050 
1051  return pass;
1052  }
1053 
1054  // Helper function for tests in Spectrum_test
1055  bool test_within_tol(double a, double b, double tol, std::string label)
1056  {
1057  // Tol is considered as a fraction of a
1058  bool pass = std::abs(a - b) <= std::abs(a*tol); // pass == true
1059  OUTPUT << "TESTING: " << label << std::endl;
1060  if(!pass)
1061  {
1062  OUTPUT << " ******FAIL******" << std::endl
1063  << " Inputs do not match within requested relative tolerance (" << tol << ")"
1064  << std::endl
1065  << " a = " << a << std::endl
1066  << " b = " << b << std::endl
1067  << " |(a - b)/a| = " << std::abs((a-b)/a) << " (greater than tol = " << tol << ")"
1068  << std::endl;
1069  }
1070  else
1071  {
1072  OUTPUT << " Pass: (a="<<a<<", b="<<b<<")" << std::endl
1073  << " |(a - b)/a| = " << std::abs((a-b)/a) << " (less than tol = " << tol << ")"
1074  << std::endl;
1075  }
1076  return pass;
1077  }
1078 
1079  // Test that output of Standard Model wrapper (e.g. QedQcdWrapper) matches
1080  // SMINPUTS sufficiently accurately
1081  // Set flag SLHAonly=1 if SMskeleton and/or MSSMskeleton are being used.
1082  void Spectrum_test(Spectrum matched_spectra, const SubSpectrum* smin, bool SLHAonly=0)
1083  {
1084  // Extract pieces of Spectrum to make it clear what they are supposed to be
1085  SMInputs sminputs = matched_spectra->get_SMInputs();
1086  std::unique_ptr<SubSpectrum> SM = matched_spectra->clone_LE(); // COPIES Spectrum object
1087  // const SubSpectrum* SM = matched_spectra->get_LE(); // Cannot do running on original object.
1088 
1089  double tol = 1e-9; // Demanding matching to 1 part in a billion (pole masses, things that don't change)
1090  double tolg = 1e-4; // Seem to get about this level of precision recovering running couplings from QedQcd object.
1091  double tolm = 1e-3; // " " " " masses " "
1092 
1093  // // SLHA1
1094  // double alphainv; // 1: Inverse electromagnetic coupling at the Z pole in the MSbar scheme (with 5 active flavours)
1095  // double GF; // 2: Fermi constant (in units of GeV^-2)
1096  // double alphaS; // 3: Strong coupling at the Z pole in the MSbar scheme (with 5 active flavours).
1097  // double mZ; // 4: Z pole mass
1098  // double mBmB; // 5: b quark running mass in the MSbar scheme (at mB)
1099  // double mT; // 6: Top quark pole mass
1100  // double mTau; // 7: Tau pole mass
1101 
1102  // // SLHA2
1103  // double mNu3; // 8: Heaviest neutrino pole mass
1104 
1105  // double mE; // 11: Electron pole mass
1106  // double mNu1; // 12: Lightest neutrino pole mass
1107  // double mMu; // 13: Muon pole mass
1108  // double mNu2; // 14: Second lightest neutrino pole mass
1109 
1110  // double mD; // 21: d quark running mass in the MSbar scheme at 2 GeV
1111  // double mU; // 21: u quark running mass in the MSbar scheme at 2 GeV
1112  // double mS; // 21: s quark running mass in the MSbar scheme at 2 GeV
1113  // double mC; // 21: c quark running mass in the MSbar scheme at mC
1114 
1115  // First check pole masses
1116  test_within_tol( sminputs.mZ, SM->get_Pole_Mass("Z"), tol, "Z pole" );
1117  //test_within_tol( sminputs.mW, SM->get_Pole_Mass("W"), tol, "W pole" ); // Whoops, no mW in sminputs.
1118  test_within_tol( sminputs.mT, SM->get_Pole_Mass("t"), tol, "top pole" );
1119  test_within_tol( sminputs.mTau, SM->get_Pole_Mass("tau"), tol, "tau pole" );
1120  test_within_tol( sminputs.mMu, SM->get_Pole_Mass("mu"), tol, "mu pole" );
1121  test_within_tol( sminputs.mE, SM->get_Pole_Mass("e"), tol, "e pole" );
1122  //test_within_tol( sminputs.mNu3, SM->get_Pole_Mass(""), tol );
1123  //test_within_tol( sminputs.mNu2, SM->get_Pole_Mass(""), tol );
1124  //test_within_tol( sminputs.mNu1, SM->get_Pole_Mass(""), tol );
1125 
1126  // Next check running quantities evaluated at Z pole
1127  // Note, numerical errors might creep in depending on how we do the running
1128  // back and forth. Might need to consider some method to "reset" object back
1129  // to original condition (keep a copy of itself inside?)
1130  //SM->RunToScale(sminputs.mZ);
1131  OUTPUT << "Current scale: " << SM->GetScale() << std::endl;
1132  OUTPUT << "Z pole mass : " << SM->get_Pole_Mass("Z") << std::endl;
1133  if(not SLHAonly) test_within_tol( sminputs.alphainv, 1./ SM->get_dimensionless_parameter("alpha"), tol, "1/alpha(mZ)" );
1134  if(not SLHAonly) test_within_tol( sminputs.alphaS, SM->get_dimensionless_parameter("alphaS"), tol, "alphaS(mZ)" );
1135 
1136  // Check running quantities evaluated at 2 GeV
1137  if(not SLHAonly) SM->RunToScale(2);
1138  OUTPUT << "Current scale: " << SM->GetScale() << std::endl;
1139  test_within_tol( sminputs.mU, SM->get_mass_parameter("u"), tolm, "mu(2)" );
1140  test_within_tol( sminputs.mD, SM->get_mass_parameter("d"), tolm, "md(2)" );
1141  test_within_tol( sminputs.mS, SM->get_mass_parameter("s"), tolm, "ms(2)" );
1142 
1143  // Check mC(mC) and mB(mB)
1144  if(not SLHAonly){
1145  SM->RunToScale(sminputs.mCmC);
1146  OUTPUT << "Current scale: " << SM->GetScale() << std::endl;
1147  OUTPUT << "mC (MSbar) : " << SM->get_mass_parameter("c") << std::endl;
1148  test_within_tol( sminputs.mCmC, SM->get_mass_parameter("c"), tolm, "mc(mc)" );
1149  SM->RunToScale(sminputs.mBmB);
1150  OUTPUT << "Current scale: " << SM->GetScale() << std::endl;
1151  OUTPUT << "mB (MSbar) : " << SM->get_mass_parameter("b") << std::endl;
1152  test_within_tol( sminputs.mBmB, SM->get_mass_parameter("b"), tolm, "mb(mb)" );
1153  OUTPUT << EOM;
1154  }
1155 
1156  // Check that pre-extracted SM SubSpectrum* and the one from Spectrum object match
1157  if(not SLHAonly) SM->RunToScale(sminputs.mZ);
1158  if(not SLHAonly) smin->RunToScale(sminputs.mZ);
1159  OUTPUT << "Checking match between SM SubSpectrum* retrieved in different ways..." << std::endl;
1160  test_within_tol(SM->get_Pole_Mass("Z"),
1161  smin->get_Pole_Mass("Z"), tol, "Z pole" );
1162  test_within_tol(SM->get_Pole_Mass("W"),
1163  smin->get_Pole_Mass("W"), tol, "W pole" );
1164  test_within_tol(SM->get_Pole_Mass("t"),
1165  smin->get_Pole_Mass("t"), tol, "top pole" );
1166  test_within_tol(SM->get_Pole_Mass("tau"),
1167  smin->get_Pole_Mass("tau"), tol, "tau pole" );
1168  test_within_tol(SM->get_Pole_Mass("mu"),
1169  smin->get_Pole_Mass("mu"), tol, "mu pole" );
1170  test_within_tol(SM->get_Pole_Mass("e"),
1171  smin->get_Pole_Mass("e"), tol, "e pole" );
1172  if(not SLHAonly) test_within_tol(SM->get_dimensionless_parameter("alpha"),
1173  smin->get_dimensionless_parameter("alpha"), tol, "1/alpha(mZ)" );
1174  if(not SLHAonly) test_within_tol(SM->get_dimensionless_parameter("alphaS"),
1175  smin->get_dimensionless_parameter("alphaS"),tol, "alphaS(mZ)" );
1176  test_within_tol(SM->get_mass_parameter("u"),
1177  smin->get_mass_parameter("u"), tolm, "mu(2)" );
1178  test_within_tol(SM->get_mass_parameter("d"),
1179  smin->get_mass_parameter("d"), tolm, "md(2)" );
1180  test_within_tol(SM->get_mass_parameter("s"),
1181  smin->get_mass_parameter("s"), tolm, "ms(2)" );
1182  if(not SLHAonly) test_within_tol(SM->get_mass_parameter("c"),
1183  smin->get_mass_parameter("c"), tolm, "mc(mc)" );
1184  if(not SLHAonly) test_within_tol(SM->get_mass_parameter("b"),
1185  smin->get_mass_parameter("b"), tolm, "mb(mb)" );
1186 
1187 
1188  // Check light quark mass ratios
1189  OUTPUT << "Checking light quark mass ratios:" << std::endl;
1190 
1191  std::vector<double> scales;
1192  scales.push_back(10);
1193  if(not SLHAonly) scales.push_back(2);
1194  if(not SLHAonly) scales.push_back(1);
1195  if(not SLHAonly) scales.push_back(0.5);
1196  if(not SLHAonly) scales.push_back(0.1);
1197 
1198  for(std::vector<double>::iterator it = scales.begin(); it != scales.end(); ++it)
1199  {
1200  if(not SLHAonly) SM->RunToScale(*it);
1201  double Q = SM->GetScale();
1202  double mu = SM->get_mass_parameter("u");
1203  double md = SM->get_mass_parameter("d");
1204  double ms = SM->get_mass_parameter("s");
1205 
1206  OUTPUT << "---------------------------------" << std::endl;
1207  OUTPUT << "Current scale: " << Q << std::endl;
1208  OUTPUT << "mu("<<Q<<") = " << mu << std::endl;
1209  OUTPUT << "md("<<Q<<") = " << md << std::endl;
1210  OUTPUT << "ms("<<Q<<") = " << ms << std::endl;
1211  OUTPUT << "mu/md = " << mu/md << std::endl;
1212  OUTPUT << "ms/md = " << ms/md << std::endl;
1213  }
1214  OUTPUT << EOM;
1215 
1217  if(not SLHAonly) {
1218 
1219  double Qs[] = {
1220  1.00000000e-02, 1.25892541e-02, 1.58489319e-02,
1221  1.99526231e-02, 2.51188643e-02, 3.16227766e-02,
1222  3.98107171e-02, 5.01187234e-02, 6.30957344e-02,
1223  7.94328235e-02, 1.00000000e-01, 1.25892541e-01,
1224  1.58489319e-01, 1.99526231e-01, 2.51188643e-01,
1225  3.16227766e-01, 3.98107171e-01, 5.01187234e-01,
1226  6.30957344e-01, 7.94328235e-01, 1.00000000e+00,
1227  1.25892541e+00, 1.58489319e+00, 1.99526231e+00,
1228  2.51188643e+00, 3.16227766e+00, 3.98107171e+00,
1229  5.01187234e+00, 6.30957344e+00, 7.94328235e+00,
1230  1.00000000e+01, 1.25892541e+01, 1.58489319e+01,
1231  1.99526231e+01, 2.51188643e+01, 3.16227766e+01,
1232  3.98107171e+01, 5.01187234e+01, 6.30957344e+01,
1233  7.94328235e+01, 1.00000000e+02, 1.25892541e+02,
1234  1.58489319e+02, 1.99526231e+02, 2.51188643e+02,
1235  3.16227766e+02, 3.98107171e+02, 5.01187234e+02,
1236  6.30957344e+02, 7.94328235e+02, 1.00000000e+03,
1237  1.25892541e+03, 1.58489319e+03, 1.99526231e+03,
1238  2.51188643e+03, 3.16227766e+03, 3.98107171e+03,
1239  5.01187234e+03, 6.30957344e+03, 7.94328235e+03,
1240  1.00000000e+04, 1.25892541e+04, 1.58489319e+04,
1241  1.99526231e+04, 2.51188643e+04, 3.16227766e+04,
1242  3.98107171e+04, 5.01187234e+04, 6.30957344e+04,
1243  7.94328235e+04
1244  };
1245 
1246  std::vector<double> Qvec(Qs, Utils::endA(Qs));
1247 
1248  std::ofstream Qout;
1249  Qout.open("SpecBit/light_quark_txt");
1250 
1251  for(std::vector<double>::iterator it = Qvec.begin(); it != Qvec.end(); ++it)
1252  {
1253  // Clone to avoid buildup of errors
1254  std::unique_ptr<SubSpectrum> SMloop = matched_spectra->clone_LE();
1255 
1256  SMloop->RunToScale(*it);
1257  double Q = SMloop->GetScale();
1258  double mu = SMloop->get_mass_parameter("u");
1259  double md = SMloop->get_mass_parameter("d");
1260  double ms = SMloop->get_mass_parameter("s");
1261 
1262  // Write to file
1263  Qout << Q << ", " << md << ", " << mu << ", " << ms << std::endl;
1264  }
1265 
1266  Qout.close();
1267 
1268  } // endif
1269 
1271  // Copy to object to clone the hosted SubSpectrum objects.
1272  // i.e. all copies are deep copies.
1273  Spectrum nonconst_spectra(*matched_spectra);
1274 
1275  // Try to access non-const member functions
1276  OUTPUT << std::endl;
1277  OUTPUT << "Testing non-const access to Spectrum object:" << std::endl;
1278  if(not SLHAonly) nonconst_spectra.RunBothToScale(sminputs.mT); // This is the only non-const function atm.
1279  OUTPUT << "Current SM SubSpectrum* scale: " << nonconst_spectra.get_LE()->GetScale() << std::endl;
1280  OUTPUT << "Current UV SubSpectrum* scale: " << nonconst_spectra.get_UV()->GetScale() << std::endl;
1281  // Make sure nothing happened to the original objects
1282  OUTPUT << "Old SM SubSpectrum* scale: " << matched_spectra->get_LE()->GetScale() << std::endl;
1283  OUTPUT << "Old UV SubSpectrum* scale: " << matched_spectra->get_UV()->GetScale() << std::endl;
1284  // Check some other numbers
1285  OUTPUT << "Current SM SubSpectrum* mu :" << nonconst_spectra.get_LE()->get_mass_parameter("u") << std::endl;
1286  OUTPUT << "Current SM Spectrum* md :" << nonconst_spectra.get_LE()->get_mass_parameter("d") << std::endl;
1287  OUTPUT << "Current SM Spectrum* ms :" << nonconst_spectra.get_LE()->get_mass_parameter("s") << std::endl;
1288  OUTPUT << "Old SM Spectrum* mu :" << matched_spectra->get_LE()->get_mass_parameter("u") << std::endl;
1289  OUTPUT << "Old SM Spectrum* md :" << matched_spectra->get_LE()->get_mass_parameter("d") << std::endl;
1290  OUTPUT << "Old SM Spectrum* ms :" << matched_spectra->get_LE()->get_mass_parameter("s") << std::endl;
1291  OUTPUT << EOM;
1292 
1293  // Check running beyond soft and hard limits (assumes QedQcdWrapper for SM)
1294  // behave = 0 -- If running beyond soft limit requested, halt at soft limit
1295  // (assumes hard limits outside of soft limits; but this is not enforced)
1296  // behave = 1 -- If running beyond soft limit requested, throw warning
1297  // " " hard limit " , throw error
1298  // behave = anything else -- Ignore limits and attempt running to requested scale
1299  if(not SLHAonly) {
1300  OUTPUT << "Testing QedQcdWrapper running limits:" << std::endl;
1301  // behave=0 (default)
1302  // Running halted at soft limit
1303  OUTPUT << "behave=0" << std::endl;
1304  SM->RunToScale(sminputs.mT); // Soft limit (and hard limit)
1305  OUTPUT << "Current scale: " << SM->GetScale() << std::endl;
1306  SM->RunToScale(1.5*sminputs.mT);
1307  OUTPUT << "Current scale: " << SM->GetScale() << std::endl;
1308  OUTPUT << EOM;
1309 
1310  // behave=2
1311  // Should be no errors, just potentially inaccurate running
1312  // EDIT: Whoops, so QedQcd object itself will throw an error if you try
1313  // to run above mT. Remove comments to observe this behaviour.
1314  OUTPUT << "behave=2" << std::endl;
1315  SM->RunToScale(sminputs.mT,2); // Soft limit (and hard limit)
1316  OUTPUT << "Current scale: " << SM->GetScale() << std::endl;
1317  //SM->RunToScale(1.5*sminputs.mT,2);
1318  //OUTPUT << "Current scale: " << SM->GetScale() << std::endl;
1319  OUTPUT << EOM;
1320 
1321  // behave=1
1322  OUTPUT << "behave=1" << std::endl;
1323  SM->RunToScale(sminputs.mT,1); // Soft limit (and hard limit)
1324  OUTPUT << "Current scale: " << SM->GetScale() << std::endl;
1325  SM->RunToScale(1.5*sminputs.mT,1); // Beyond hard limit (error)
1326  OUTPUT << "Current scale: " << SM->GetScale() << std::endl;
1327  OUTPUT << EOM;
1328  }
1329  }
1330 
1331  } // end namespace SpecBit
1332 } // end namespace Gambit
1333 
1334 #endif
bool print_error(bool pass, std::string get_type, std::string data, double sting_get_out, double data_member, int i=-1, int j=-1)
void RunToScale(double scale, const int behave=0)
Run spectrum to a new scale This function is a wrapper for RunToScaleOverride which automatically che...
greatScanData data
Definition: great.cpp:38
MSSM derived version of SubSpectrum class.
bool test_within_tol(double a, double b, double tol, std::string label)
#define TAGeom
bool TestMssmParGets(SubSpectrum *spec, M FSmssm)
bool TestMssmParMass1_0(SubSpectrum *spec, M FSmssm, bool immediate_exit=true)
virtual double GetScale() const
Returns the renormalisation scale of parameters.
Definition: MSSMSpec.hpp:198
std::unique_ptr< SubSpectrum > clone_LE() const
Clone getters Note: If you want to clone the whole Spectrum object, just use copy constructor...
Definition: spectrum.cpp:234
bool TestMssmPoleGets(SubSpectrum *spec, M FSmssm)
virtual double GetScale() const
Returns the renormalisation scale of parameters.
bool running_test(MSSMSpec< MI > &mssm, typename MI::Model &FS_model_slha, double tol)
Flexiblesusy model header includes for SpecBit.
START_MODEL b
Definition: demo.hpp:235
bool TestMssmPoleGets0(SubSpectrum *spec, M FSmssm, bool immediate_exit=true)
bool TestMssmParMass1_2(SubSpectrum *spec, M FSmssm, bool immediate_exit=true)
bool test_exact(MSSMSpec< MI > &mssm, typename MI::Model &FS_model_slha)
std::chrono::milliseconds ms
bool TestMssmParMass0_0(SubSpectrum *spec, M FSmssm, bool immediate_exit=true)
bool TestMssmParMass2_0(SubSpectrum *spec, M FSmssm, bool immediate_exit=true)
Module convenience functions.
#define TAGerr
bool TestMssmPoleGets1(SubSpectrum *spec, M FSmssm, bool immediate_exit=true)
const double mu
Definition: SM_Z.hpp:42
#define TAGfatal
const Logging::endofmessage EOM
Explicit const instance of the end of message struct in Gambit namespace.
Definition: logger.hpp:99
void setup(Model &mssm)
SubSpectrum & get_LE()
Standard getters Return references to internal data members.
Definition: spectrum.cpp:224
Virtual base class for interacting with spectrum generator output.
Definition: subspectrum.hpp:87
START_MODEL M
bool test_getters(std::string get_type, std::string name, double getter_output, double data_member, int i=-1, int j=-1)
bool TestMssmParMass2_2(SubSpectrum *spec, M FSmssm, bool immediate_exit=true)
bool TestMssmPoleMixingGets2(SubSpectrum *spec, M FSmssm, bool immediate_exit=true)
void Spectrum_test(Spectrum matched_spectra, const SubSpectrum *smin, bool SLHAonly=0)
T * endA(T(&arr)[N])
#define OUTPUT
TODO: see if we can use this one:
Definition: Analysis.hpp:33
"Standard Model" (low-energy) plus high-energy model container class
Definition: spectrum.hpp:110
bool TestMssmParMass0_2(SubSpectrum *spec, M FSmssm, bool immediate_exit=true)
SMInputs & get_SMInputs()
Definition: spectrum.cpp:226
void RunBothToScale(double scale)
Linked running Only possible with non-const object.
Definition: spectrum.cpp:153
Container class for Standard Model input information (defined as in SLHA2)
Definition: sminputs.hpp:29