Как соединить оба метода?

double function_optimiz(double alpha) {
    double y[N];
    y[0] = 10e-10, y[1] = 0.0, y[2] = 0.0, y[3] = 0.0;
    double t0 = 0.0, max_t = 5.73e-08, h = 1e-10;
    double t = t0;
    double sum = 0;
 
    for(t = t0; t <= max_t; t += h) {
 
        rosenbrock(y, t0, max_t, h, alpha);
        std::cout << t << "\t" << std::setprecision(5) << y[0] << "\t" << y[1] << "\t" << y[2] << "\t" << y[3] << std::endl;
 
        sum += (j_0(y) - j_1(y))*(j_0(y) - j_1(y));
        cout << "sum: " << sum << "\n";
 
    }
    return sum;
}
void Gold(double a, double b, double e, double alpha) {
    long nIter = 1;
 
    double t = (std::sqrt(5) - 1)*0.5;
    double x1 = b - t * (b - a);
    double x2 = a + t * (b - a);
 
    double h1, h2;
 
    while (e <= fabs(b - a)) {
        std::cout << "\t**\t Iteration: " << (nIter++) << " \t **" << std::endl;
        std::cout << std::setprecision(4) << "[ a ; b ]: [ " << a << " ; " << b << " ]\t" << " ALPHA_1: " << x1 << "\tALPHA_2: " << x2 << std::endl;
        h2 = function_optimiz(x2);
        h1 = function_optimiz(x1);
        if (h2 < h1) a = x1;
        else b = x2;
 
        x1 = b - t * (b - a);
        x2 = a + t * (b - a);
 
    }
 
    std::cout << "\t**\t Iteration: " << (nIter++) << " \t **" << std::endl;
    std::cout << std::setprecision(4) << "[ a ; b ]: [ " << a << " ; " << b << " ]\t" << " ALPHA_1: " << x1 << "\tALPHA_2: " << x2 << std::endl;
 
    std::cout << "\n\t**\t Result \t **";
    std::cout << "\nALPHA_min =  " << (alpha = (a + b) / 2);
 
}

void descent_method(func_ptr f, double* vars, double eps, int max_steps_count)
{
double B = f(vars), A = 0; // значение целевой функции на текущей и пердидущей итерациях
bool was_counted = false; // был ли определен минимум за максимально допустимое количество итераций
int stpes_ellapsed = 0; // количество итераций, затраченных на поиск минимума (статистическая информация)
double delta = 0.0; // полученная погрешность (статистическая информация)
for(int i = 0; i < max_steps_count; i++){
A = B; // предыдущее значение целевой функции
 
for(int var_index = 0; var_index < var_count; var_index++) // проходим по каждому аргументу целевой функции
vars[var_index] = golden_section(f, vars, var_index, eps, -5000, 5000, max_steps_count); // и находим минимум одномерной функции по этой пе
 
B = f(vars); // текущее значение целевой функции
 
delta = fabs(A - B);
 
if(delta <= eps)
{
stpes_ellapsed = i + 1;
was_counted = true;
break;
}
}
 
std::cout << "Результат поиска минимума функции " << "<тут пишем функцию>" << std::endl;
 
if(!was_counted)
std::cout << "За максимально указанное количество шагов ( " << max_steps_count << " ) минимум не был посчитан." << std::endl;
else {
std::cout << "Количество итераций: " << stpes_ellapsed << std::endl;
std::cout << "Погрешность: " << delta << std::endl;
}
 
std::cout << "Точка: X(";
for(int i = 0; i < var_count; i++){
std::cout << vars[i] << ", ";
 
}
 
std::cout << "\b\b" << ")" << std::endl;
std::cout << "Значение функции f(X): " << std::setprecision(10) << f(vars) << std::endl;
}

typedef double (*func_ptr)(double*);
double function(double* variables);

double function_optimiz(double alpha) {
    double y[N];
    // Что за N? Где инициализация
}

#include <iostream>
#include <iomanip>
#include <fstream>
#include <functional>
#include <vector>
#include <boost/multiprecision/cpp_bin_float.hpp>
 
using namespace boost::multiprecision;
using namespace std;
 
typedef number<cpp_bin_float<200>> cpp_dec_float_unlim;
 
#define double cpp_dec_float_unlim
#define N 4
 
// система уравнений
void dydt(double y[N], double f[N], double t, double alpha) {
    f[0] = 8.85*10e6 - y[2] * y[0] * (0.2 - y[0]) / (1. - y[1]) + y[0] / 10e-8;
    f[1] = y[0] * y[2];
    f[2] = y[0] * y[3];
    f[3] = 10e30*(sqrt(y[0]))*exp((-alpha)/ y[0]);
}
 
// Функция для получания кофактора A[p][q] в temp[][]. n текущее 
// измерение A[][] 
void getCofactor(double A[N][N], double temp[N][N], int p, int q, int n)
{
    int i = 0, j = 0;
    // Циклы для каждого элемента матрицы
    for (int row = 0; row < n; row++)
    {
        for (int col = 0; col < n; col++)
        {
            //  Копирование во временную матрицу только тех элементов, 
            // которых нет в данной строке и столбце
            if (row != p && col != q)
            {
                temp[i][j++] = A[row][col];
                // Строка заполнена, поэтому увеличьте индекс строки и 
                // сбросьте индекс col
                if (j == n - 1)
                {
                    j = 0;
                    i++;
                }
            }
        }
    }
}
 
/* Рекурсивная функция для нахождения определителя матрицы.
   n текущее измерение A[][]. */
double determinant(double A[N][N], int n)
{
    double D = 0; // Initialize result 
    //  Base case : if matrix contains single element 
    if (n == 1)
        return A[0][0];
    double temp[N][N]; // To store cofactors 
    int sign = 1;  // To store sign multiplier 
     // Iterate for each element of first row 
    for (int f = 0; f < n; f++)
    {
        // Getting Cofactor of A[0][f] 
        getCofactor(A, temp, 0, f, n);
        D += sign * A[0][f] * determinant(temp, n - 1);
 
        // terms are to be added with alternate sign 
        sign = -sign;
    }
    return D;
}
 
// Function to get adjoint of A[N][N] in adj[N][N]. 
void adjoint(double A[N][N], double adj[N][N])
{   
    if (N == 1)
    {
        adj[0][0] = 1;
        return;
    }
 
    // temp is used to store cofactors of A[][] 
 
    int sign = 1;
    double temp[N][N];
    for (int i = 0; i < N; i++)
    {
        for (int j = 0; j < N; j++)
        {
            // Get cofactor of A[i][j] 
            getCofactor(A, temp, i, j, N);
 
            // sign of adj[j][i] positive if sum of row 
            // and column indexes is even. 
            sign = ((i + j) % 2 == 0) ? 1 : -1;
 
            // Interchanging rows and columns to get the 
            // transpose of the cofactor matrix 
            adj[j][i] = (sign)*(determinant(temp, N - 1));
        }
    }
}
 
// Function to calculate and store inverse, returns false if 
// matrix is singular 
bool inverse(double A[N][N], double inverse[N][N])
{
    // Find determinant of A[][] 
    double det = determinant(A, N);
    if (det == 0)
    {
        cout << "Singular matrix, can't find its inverse";
        return false;
    }
 
 
    // Find adjoint 
    double adj[N][N];
 
    adjoint(A, adj);
 
    // Find Inverse using formula "inverse(A) = adj(A)/det(A)" 
    for (int i = 0; i < N; i++)
        for (int j = 0; j < N; j++)
            inverse[i][j] = adj[i][j] / double(det);
    return true;
}
 
// Функция вывода матрицы
void PrintMatr(double mas[N][N], int n) {
    for (int i = 0; i < n; i++) {
        for (int j = 0; j < n; j++)
            std::cout << mas[i][j] << " ";
        std::cout << endl;
    }
}
 
//Умножение матрицы на вектор
void mul_1(double matr[N][N], double v[N], double result[N]) {
    for (int i = 0; i < N; ++i) {
        result[i] = 0;
        for (int j = 0; j < N; ++j) {
            result[i] += matr[i][j] * v[j];
        }
    }
}
 
//Умножение вектора на число
void mul_2(double v[N], double h, double result[N]) {
    for (int i = 0; i < N; i++) {
        result[i] = v[i] * h;
    }
}
 
//Сложение векторов
void add(double v1[N], double v2[N], double result[N]) {
    for (int i = 0; i < N; i++) {
        result[i] = v1[i] + v2[i];
    }
}
 
//Якобиан системы
void SetJackobian(double y[N], double jac[N][N], double alpha) {
    jac[0][0] = (2. * (y[2] * (y[0] - 0.1))) / (1.0 - y[1]) + 10e8;
    jac[0][1] = ((y[0] - 0.2)*y[0] * y[2]) / pow((1.0 - y[1]), 2.);
    jac[0][2] = ((y[0] - 0.2)*y[0]) / (1.0 - y[1]);
    jac[0][3] = 0.;
    jac[1][0] = y[2];
    jac[1][1] = 0.;
    jac[1][2] = y[0];
    jac[1][3] = 0.;
    jac[2][0] = y[3];
    jac[2][1] = 0.;
    jac[2][2] = 0.;
    jac[2][3] = y[0];
    jac[3][0] = 5.e29*(exp(-alpha / y[0]))*(0.5 + y[0]) / pow(y[0], 3. / 2.);
    jac[3][1] = 0.;
    jac[3][2] = 0.;
    jac[3][3] = 0.;
}
 
//Одностадийный метод Розенброка
void rosenbrock(double y[N], double t0, double max_t, double h, double alpha) {
    double t = t0;
    double ynew[N];
    double f[N];
    double jac[N][N];
    double adj[N][N];  // To adjoint of A[][] 
    double inv[N][N]; // To store inverse of A[][] 
    double a[N][N];
    
    //единичная матрица
    for (int i = 0; i < N; i++) {
        for (int j = 0; j < N; j++) {
            if (i == j) a[i][j] = 1;
            else a[i][j] = 0;
        }
    }
 
    SetJackobian(y, jac, alpha);
 
    //умножение матрицы на число
    double umn[N][N];
    for (int i = 0; i < N; i++) {
        for (int j = 0; j < N; j++) {
            jac[i][j] *= (h / 2);
            umn[i][j] = jac[i][j];
        }
    }
 
    // вычитание матрицы
    double temp[N][N];
    for (int i = 0; i < N; i++) {
        for (int j = 0; j < N; j++) {
            temp[i][j] = a[i][j] - umn[i][j];
        }
    }
 
    adjoint(temp, adj);
 
    inverse(temp, inv);
 
    dydt(y, f, t, alpha);
    mul_1(inv, f, ynew);
    mul_2(ynew, h, ynew);
    add(ynew, y, y);
 
}
 
double j_0(double y[N]) {
    return -(2 * 3 * 10e-2*((8.85*10e6) - (y[0] / 10e-8) - y[2])*(1 - y[1])) / (0.2 - y[0]);
}
 
double j_1(double y[N]) {
    return j_0(y)*0.1;
}
 
double function_optimiz(double alpha) {
    double y[N];
    y[0] = 10e-10, y[1] = 0.0, y[2] = 0.0, y[3] = 0.0;
    double t0 = 0.0, max_t = 5.73e-08, h = 1e-10;
    double t = t0;
    double sum = 0;
 
    for(t = t0; t <= max_t; t += h) {
 
        rosenbrock(y, t0, max_t, h, alpha);
        std::cout << t << "\t" << std::setprecision(5) << y[0] << "\t" << y[1] << "\t" << y[2] << "\t" << y[3] << std::endl;
 
        sum += (j_0(y) - j_1(y))*(j_0(y) - j_1(y));
        cout << "sum: " << sum << "\n";
 
    }
    return sum;
}
 
double Gold(double a, double b, double e, int max_steps_count) {
    long nIter = 1; // переменная для показа номера итерации
    double res = 0.0; // результат, возвращаемый функцией
    double t = (std::sqrt(5) - 1)/2;// магический коэффициент золотого сечения
    double x1 = b - t * (b - a), x2 = a + t * (b - a); // предполагаемая точка экстремума и точка симметричная ей на отрезке [a, b]
 
    int step = 0; // текущий шаг 
    double h1, h2;
 
    while (e <= fabs(b - a)) {  // повторяем цикл до тех пор, пока разность концов отрезков не будет меньше заданной погрешности
        std::cout << "\t**\t Iteration: " << (nIter++) << " \t **" << std::endl;
        std::cout << std::setprecision(4) << "[ a ; b ]: [ " << a << " ; " << b << " ]\t" << " ALPHA_1: " << x1 << "\tALPHA_2: " << x2 << std::endl;
        h2 = function_optimiz(x2);  // считаем значение функции в точке симетричной точке экстремума
        h1 = function_optimiz(x1); // считаем значение функции в точке экстремума
        // сужаем интервал неопределенности
        if (h2 < h1) a = x1;
        else b = x2;
 
        x1 = b - t * (b - a); // пересчитываем предполагаемую точку экстремума
        x2 = a + t * (b - a); // пересчитываем предполагаемую точку симетричную экстремуму
 
        step++; 
        if (step > max_steps_count) break;
    }
 
    std::cout << "\t**\t Iteration: " << (nIter++) << " \t **" << std::endl;
    std::cout << std::setprecision(4) << "[ a ; b ]: [ " << a << " ; " << b << " ]\t" << " ALPHA_1: " << x1 << "\tALPHA_2: " << x2 << std::endl;
 
    std::cout << "\n\t**\t Result \t **";
    res = (a + b) / 2;
    std::cout << "\n res = " << res; 
 
    return res;
}
 
 
void descent_method(double a, double b,  double eps, int max_steps_count)
{
    double alpha = ;
    double B = function_optimiz(alpha), A = 0;
    bool was_counted = false;
    int stpes_ellapsed = 0;
    double delta = 0.0;
 
 
    for (int i = 0; i < max_steps_count; i++) {
        A = B;
        alpha = Gold(a, b, eps, max_steps_count);
 
        B = function_optimiz(alpha);
 
        if (fabs(A - B) <= eps)
        {
            stpes_ellapsed = i + 1;
            was_counted = true;
            break;
        }
    }
}
 
// Driver program 
 
int main()
{
    double a, b, e;
    cout << "Input a: "; cin >> a;   // 8e-2
    cout << "Input b: "; cin >> b;   // 5e-1
    cout << "Input eps: "; cin >> e; // 1e-3
 
    int max_steps_count;
    cout << "Input max steps count: "; cin >> max_steps_count; // 20
 
    cout << "Gold method: " << endl;
    Gold(a, b, e, max_steps_count);
 
    cout << endl;
 
    system("pause");
 
    return 0;
}

double Gold(double a, double b, double e, int max_steps_count) {
    double res = 0.0; // результат, возвращаемый функцией
    double t = (std::sqrt(5) - 1)/2;// магический коэффициент золотого сечения
    double x1 = b - t * (b - a), x2 = a + t * (b - a); // предполагаемая точка экстремума и точка симметричная ей на отрезке [a, b]
 
    int step = 0; // текущий шаг
    double h1, h2;
 
    while (e <= fabs(b - a)) {  // повторяем цикл до тех пор, пока разность концов отрезков не будет меньше заданной погрешности
        h2 = function_optimiz(x2);  // считаем значение функции в точке симетричной точке экстремума
        h1 = function_optimiz(x1); // считаем значение функции в точке экстремума
        // сужаем интервал неопределенности
        if (h2 < h1) a = x1;
        else b = x2;
 
        x1 = b - t * (b - a); // пересчитываем предполагаемую точку экстремума
        x2 = a + t * (b - a); // пересчитываем предполагаемую точку симетричную экстремуму
 
        step++;
        if (step > max_steps_count) break;
    }
 
    res = (a + b) / 2;
    return res;
}

void descent_method(double a, double b,  double eps, int max_steps_count)
{
    double alpha = Gold(a, b, eps, max_steps_count);
    double B = function_optimiz(alpha), A = 0;
    bool was_counted = false;
    int stpes_ellapsed = 0;
 
    for (int i = 0; i < max_steps_count; i++) {
        A = B;
        alpha = Gold(a, b, eps, max_steps_count);
 
        B = function_optimiz(alpha);
 
        if (fabs(A - B) <= eps)
        {
            stpes_ellapsed = i + 1;
            was_counted = true;
            break;
        }
    }
}

double Gold(double &a, double &b, double e, int max_steps_count) {
    double res = 0.0; // результат, возвращаемый функцией
    double t = (std::sqrt(5) - 1)/2;// магический коэффициент золотого сечения
    double x1 = b - t * (b - a), x2 = a + t * (b - a); // предполагаемая точка экстремума и точка симметричная ей на отрезке [a, b]
 
    int step = 0; // текущий шаг
    double h1, h2;
 
    while (e <= fabs(b - a)) {  // повторяем цикл до тех пор, пока разность концов отрезков не будет меньше заданной погрешности
        h2 = function_optimiz(x2);  // считаем значение функции в точке симетричной точке экстремума
        h1 = function_optimiz(x1); // считаем значение функции в точке экстремума
        // сужаем интервал неопределенности
        if (h2 < h1) a = x1;
        else b = x2;
 
        x1 = b - t * (b - a); // пересчитываем предполагаемую точку экстремума
        x2 = a + t * (b - a); // пересчитываем предполагаемую точку симетричную экстремуму
 
        step++;
        if (step > max_steps_count) break;
    }
 
    res = (a + b) / 2;
    return res;
}
 
void descent_method(double &a, double &b,  double eps, int max_steps_count)
{
    double alpha = Gold(a, b, eps, max_steps_count);
    double B = function_optimiz(alpha), A = 0;
    bool was_counted = false;
    int stpes_ellapsed = 0;
 
    for (int i = 0; i < max_steps_count; i++) {
        A = B;
        alpha = Gold(a, b, eps, max_steps_count);
 
        B = function_optimiz(alpha);
 
        if (fabs(A - B) <= eps)
        {
            stpes_ellapsed = i + 1;
            was_counted = true;
            break;
        }
    }
}

// сужаем интервал неопределенности
if (h2 < h1) a = x1;
else b = x2;