//code fo bisection

//bisection method
#include<iostream>
#include<cmath>
#include<iomanip>
using namespace std;
double f(double x);    //declare the function for the given equation
double f(double x)    //define the function here, ie give the equation
{
    double a=pow(x,3)-(6*pow(x,2))+(11.0*x)-6;    //write the equation whose roots are to be determined
    return a;
}
int main()
{
    cout.precision(4);        //set the precision
    cout.setf(ios::fixed);
    double a,b,c,e,fa,fb,fc;    //declare some needed variables
    a:cout<<"Enter the initial guesses:\na=";    //Enter the value of a(set a label('a:') for later use with goto)
    cin>>a;
    cout<<"\nb=";            //Enter the value of b
    cin>>b;
    cout<<"\nEnter the degree of accuracy desired, eg E=0.000002"<<endl;    //Enter the accuracy
    cin>>e;                //e stands for  accuracy
    if (f(a)*f(b)>0)        //Check if a root exists between a and b
    {                //If f(a)*f(b)>0 then the root does not exist between a and b
        cout<<"No root exist, enter another guess please"<<endl;
        goto a;            //go back to 'a' ie 17 and ask for different values of a and b
    }
    else                //else a root exists between a and b
    {

    cout<<"There is a root, I am calculating now: hold on " << endl;
    while (fabs(a-b)>=e)        /*if the mod of a-b is greater than the accuracy desired keep                         bisecting the interval*/
    {
        c=(a+b)/2.0;        //bisect the interval and find the value of c
        fa=f(a);
        fb=f(b);
        fc=f(c);
        cout<<"a="<<a<<"     "<<"b="<<b<<"     "<<"c="<<c<<"      fc="<<fc<<endl;/*print the                             values of a,b,c and fc  after each iteration*/
        if (fc==0)        //if f(c)=0, that means we have found the root of the equation
        {
            cout<<"The root of the equation is "<<c;    /*print the root of the equation                                         and break out of the loop*/
            break;
        }

        if (fa*fc>0)    //if f(a)xf(c)>0, that means no root exist between a and c
        {
            a=c;    /*hence make a=c, ie make c the starting point of the interval and b the                     end point*/
        }
        else if (fa*fc<0)
        {
            b=c;    /*this means that a root exist between a and c therfore make c the end                     point of the interval*/
        }


    }
    }            //The loop ends when the difference between a and b becomes less than the desired accuracy ie now the value stored in 'c' can be called the approximate root of the equation
    cout<<"The root of the equation is "<<c;    //print the root
    return 0;
}










//##################################################################################################################
// this program for solving linear equations
// C++ program to demostrate working of Guassian Elimination
// method
#include<bits/stdc++.h>
using namespace std;
#include <iomanip>

#define N 3  // Number of unknowns

// function to reduce matrix to r.e.f. Returns a value to
// indicate whether matrix is singular or not
int forwardElim(double mat[N][N+1]);

// function to calculate the values of the unknowns
void backSub(double mat[N][N+1]);

// function to get matrix content
void gaussianElimination(double mat[N][N+1])
{
    /* reduction into r.e.f. */
    int singular_flag = forwardElim(mat);

    /* if matrix is singular */
    if (singular_flag != -1)
    {
        printf("Singular Matrix.\n");

        /* if the RHS of equation corresponding to
        zero row is 0, * system has infinitely
        many solutions, else inconsistent*/
        if (mat[singular_flag][N])
            printf("Inconsistent System.");
        else
            printf("May have infinitely many "
                "solutions.");

        return;
    }

    /* get solution to system and print it using
    backward substitution */
    backSub(mat);
}

// function for elemntary operation of swapping two rows
void swap_row(double mat[N][N+1], int i, int j)
{
    //printf("Swapped rows %d and %d\n", i, j);

    for (int k=0; k<=N; k++)
    {
        double temp = mat[i][k];
        mat[i][k] = mat[j][k];
        mat[j][k] = temp;
    }
}

// function to print matrix content at any stage
void print(double mat[N][N+1])
{
    for (int i=0; i<N; i++, printf("\n"))
        for (int j=0; j<=N; j++)
            printf("%lf ", mat[i][j]);

    printf("\n");
}

// function to reduce matrix to r.e.f.
int forwardElim(double mat[N][N+1])
{
    for (int k=0; k<N; k++)
    {
        // Initialize maximum value and index for pivot
        int i_max = k;
        int v_max = mat[i_max][k];

        /* find greater amplitude for pivot if any */
        for (int i = k+1; i < N; i++)
            if (abs(mat[i][k]) > v_max)
                v_max = mat[i][k], i_max = i;

        /* if a prinicipal diagonal element is zero,
        * it denotes that matrix is singular, and
        * will lead to a division-by-zero later. */
        if (!mat[k][i_max])
            return k; // Matrix is singular

        /* Swap the greatest value row with current row */
        if (i_max != k)
            swap_row(mat, k, i_max);


        for (int i=k+1; i<N; i++)
        {
            /* factor f to set current row kth elemnt to 0,
            * and subsequently remaining kth column to 0 */
            double f = mat[i][k]/mat[k][k];

            /* subtract fth multiple of corresponding kth
            row element*/
            for (int j=k+1; j<=N; j++)
                mat[i][j] -= mat[k][j]*f;

            /* filling lower triangular matrix with zeros*/
            mat[i][k] = 0;
        }

        //print(mat);    //for matrix state
    }
    //print(mat);        //for matrix state
    return -1;
}

// function to calculate the values of the unknowns
void backSub(double mat[N][N+1])
{
    double x[N]; // An array to store solution

    /* Start calculating from last equation up to the
    first */
    for (int i = N-1; i >= 0; i--)
    {
        /* start with the RHS of the equation */
        x[i] = mat[i][N];

        /* Initialize j to i+1 since matrix is upper
        triangular*/
        for (int j=i+1; j<N; j++)
        {
            /* subtract all the lhs values
            * except the coefficient of the variable
            * whose value is being calculated */
            x[i] -= mat[i][j]*x[j];
        }

        /* divide the RHS by the coefficient of the
        unknown being calculated */
        x[i] = x[i]/mat[i][i];
    }

    printf("\nSolution for the system:\n");
    for (int i=0; i<N; i++)
        printf("%lf\n", x[i]);
}

// Driver program
int main()
{

    //don't touch it, it works   check this link: https://pics.me.me/most-of-our-code-the-parts-copied-from-stack-overflow-38269718.png
    /* input matrix */

    cout<<fixed << setprecision(2)<<showpoint ;
    double mat[N][N+1] = {{8.0, 2.0,-4.0, 84.0},
                        {2.0, 3.0, 3.0, 15.0},
                        {5.0, -3,7.0, 14.0}
                        };

    gaussianElimination(mat);

    return 0;
}