# Solve system of 6 quadratic equations

Hi all,

I’m attempting to solve a system of quadratic equations using the GSLMultiRootFinder (ROOT::Math::GSLMultiRootFinder) following this example:
ROOT: tutorials/math/exampleMultiRoot.C File Reference

I can easily extend the example to a set of three quadratic equations

``````ROOT::Math::MultiRootFinder* r = new ROOT::Math::MultiRootFinder();

TF3 * f1 = new TF3("f1","[0]*[0] - x*x - y*y");
TF3 * f2 = new TF3("f2","[0]*[0] - x*x - z*z");
TF3 * f3 = new TF3("f3","[0]*[0] - y*y - z*z");

double sigma_01 = 1.23803;  double err_sigma_01 = 0.00668558;
double sigma_02 = 0.846751; double err_sigma_02 = 0.00586008;
double sigma_12 = 1.24057;  double err_sigma_12 = 0.00926783;

f1->SetParameter(0,sigma_01);
f2->SetParameter(0,sigma_02);
f3->SetParameter(0,sigma_12);

// wrap the functions
ROOT::Math::WrappedMultiTF1 g1(*f1,3);
ROOT::Math::WrappedMultiTF1 g2(*f2,3);
ROOT::Math::WrappedMultiTF1 g3(*f3,3);
r->SetPrintLevel(printlevel);

// starting point
double x0[3]={0.5, 1.0, 0.6};
r->Solve(x0);

std::cout << "Testing result: " << f1->Eval(x0[0], x0[1], x0[2]) << std::endl;
std::cout << "Testing result: " << f2->Eval(x0[0], x0[1], x0[2]) << std::endl;
std::cout << "Testing result: " << f3->Eval(x0[0], x0[1], x0[2]) << std::endl;

``````

However, I don’t know how to treat the functions in case of six because there is no such thing as a TF6.
How can I feed the 6-dimensional functions to the root finder?

``````    // This does not work -- how to initiate a multi-dimensional function?
ROOT::Math::IMultiGradFunction* f1 = new ROOT::Math::IMultiGradFunction("f1","[0]*[0] - x[0]*x[0] - (x[1]*x[1] + x[2]*x[2] + x[3]*x[3] + x[4]*x[4] + x[5]*x[5])/5");
ROOT::Math::IMultiGradFunction* f2 = new ROOT::Math::IMultiGradFunction("f2","[0]*[0] - x[1]*x[1] - (x[0]*x[0] + x[2]*x[2] + x[3]*x[3] + x[4]*x[4] + x[5]*x[5])/5");
ROOT::Math::IMultiGradFunction* f3 = new ROOT::Math::IMultiGradFunction("f3","[0]*[0] - x[2]*x[2] - (x[0]*x[0] + x[1]*x[1] + x[3]*x[3] + x[4]*x[4] + x[5]*x[5])/5");
ROOT::Math::IMultiGradFunction* f4 = new ROOT::Math::IMultiGradFunction("f4","[0]*[0] - x[3]*x[3] - (x[0]*x[0] + x[1]*x[1] + x[2]*x[2] + x[4]*x[4] + x[5]*x[5])/5");
ROOT::Math::IMultiGradFunction* f5 = new ROOT::Math::IMultiGradFunction("f5","[0]*[0] - x[4]*x[4] - (x[0]*x[0] + x[1]*x[1] + x[2]*x[2] + x[3]*x[3] + x[5]*x[5])/5");
ROOT::Math::IMultiGradFunction* f6 = new ROOT::Math::IMultiGradFunction("f6","[0]*[0] - x[5]*x[5] - (x[0]*x[0] + x[1]*x[1] + x[2]*x[2] + x[3]*x[3] + x[4]*x[4])/5");

double sigma_0 = 1.5;
double sigma_1 = 1.4;
double sigma_2 = 1.5;
double sigma_3 = 1.5;
double sigma_4 = 1.6;
double sigma_5 = 1.5;
f1->SetParameter(0,sigma_0);
f2->SetParameter(1,sigma_1);
f3->SetParameter(2,sigma_2);
f4->SetParameter(3,sigma_3);
f5->SetParameter(4,sigma_4);
f6->SetParameter(5,sigma_5);

ROOT::Math::WrappedMultiTF1 func1(*f1,6);
ROOT::Math::WrappedMultiTF1 func2(*f2,6);
ROOT::Math::WrappedMultiTF1 func3(*f3,6);
ROOT::Math::WrappedMultiTF1 func4(*f4,6);
ROOT::Math::WrappedMultiTF1 func5(*f5,6);
ROOT::Math::WrappedMultiTF1 func6(*f6,6);

// starting point
double x0[6]={1.5, 1.5, 1.5, 1.5, 1.5, 1.5};
r->Solve(x0);
``````

Hi @jenskroeger,

Cheers,
J.

Thanks @jalopezg for the quick reply.

@moneta, @jonas, could you have a quick look at this?
I’m sure it’s a trivial issue and a small technicality only…

Many thanks,
Jens

Hi @jenskroeger,

indeed you can solve your issue easily. First of all, you can avoid creating the function wrapper yourself by using the templated MultiRootFinder::AddFunction() version, where you also pass the dimension directly. This will create a (WrappedMultiFunction)[root/WrappedFunction.h at master · root-project/root · GitHub] for you, and you can pass any function that can be called with a `const double * x` array of input values and returns a `double`.

To create the functions that you are wrapping, you can’t use the `ROOT::Math::IMultiGradFunction` class, because it’s an abstract class. Instead, you can create lambda functions, if you want to use some modern C++:

``````void solveDim6EquationSystem() {

ROOT::Math::MultiRootFinder* r = new ROOT::Math::MultiRootFinder();

double sigma_0 = 1.5;
double sigma_1 = 1.4;
double sigma_2 = 1.5;
double sigma_3 = 1.5;
double sigma_4 = 1.6;
double sigma_5 = 1.5;

auto f1 = [sigma_0](const double * x)
{ return sigma_0*sigma_0 - x[0]*x[0] - (x[1]*x[1] + x[2]*x[2]
+ x[3]*x[3] + x[4]*x[4] + x[5]*x[5])/5; };
auto f2 = [sigma_1](const double * x)
{ return sigma_1*sigma_1 - x[1]*x[1] - (x[0]*x[0] + x[2]*x[2]
+ x[3]*x[3] + x[4]*x[4] + x[5]*x[5])/5; };
auto f3 = [sigma_2](const double * x)
{ return sigma_2*sigma_2 - x[2]*x[2] - (x[0]*x[0] + x[1]*x[1]
+ x[3]*x[3] + x[4]*x[4] + x[5]*x[5])/5; };
auto f4 = [sigma_3](const double * x)
{ return sigma_3*sigma_3 - x[3]*x[3] - (x[0]*x[0] + x[1]*x[1]
+ x[2]*x[2] + x[4]*x[4] + x[5]*x[5])/5; };
auto f5 = [sigma_4](const double * x)
{ return sigma_4*sigma_4 - x[4]*x[4] - (x[0]*x[0] + x[1]*x[1]
+ x[2]*x[2] + x[3]*x[3] + x[5]*x[5])/5; };
auto f6 = [sigma_5](const double * x)
{ return sigma_5*sigma_5 - x[5]*x[5] - (x[0]*x[0] + x[1]*x[1]
+ x[2]*x[2] + x[3]*x[3] + x[4]*x[4])/5; };

// starting point
double x0[6]={1.5, 1.5, 1.5, 1.5, 1.5, 1.5};

int printlevel = 3;
r->SetPrintLevel(printlevel);

r->Solve(x0);
}
``````

Is this minimal solution fine for you? If you have further questions and need advice, please feel free to ask!

Cheers,
Jonas

Great, thanks a lot! This was exactly what I was looking for