Bisection method
#include <iostream>
#include <cmath>
using namespace std;
double f(double x) {
return x*x*x - x - 2;
double bisection(double a, double b, double tolerance) {
if (f(a) * f(b) >= 0) {
cout << "Invalid interval: f(a) and f(b) must have opposite signs!" << endl;
return -1;
double mid;
while ((b - a) >= tolerance) {
mid = (a + b) / 2.0;
if (f(mid) == 0.0)
break;
if (f(a) * f(mid) < 0)
b = mid;
else
a = mid;
}
return mid;
int main() {
double a = 1, b = 2;
double tolerance = 0.0001;
double root = bisection(a, b, tolerance);
if (root != -1) {
cout << "Approximate root = " << root << endl;
cout << "f(root) = " << f(root) << endl;
return 0;
REGULA FALSI METHOD
#include <iostream>
#include <cmath>
using namespace std;
double f(double x) {
return x*x*x - x - 2;
double regulaFalsi(double a, double b, double tolerance, int maxIter) {
double fa = f(a);
double fb = f(b);
if (fa * fb >= 0) {
cout << "Invalid interval: f(a) and f(b) must have opposite signs!" << endl;
return -1;
double c = a, fc;
for (int iter = 1; iter <= maxIter; ++iter) {
if (fb - fa == 0) {
cout << "Division by zero in formula. Stopping." << endl;
return -1;
c = (a * fb - b * fa) / (fb - fa);
fc = f(c);
cout << "Iter " << iter << " | a=" << a << " b=" << b
<< " c=" << c << " f(c)=" << fc << endl;
if (fabs(fc) < tolerance || fabs(b - a) < tolerance) {
cout << "Converged after " << iter << " iterations." << endl;
return c;
if (fa * fc < 0) {
b = c;
fb = fc;
} else {
a = c;
fa = fc;
cout << "Max iterations reached. Returning best c." << endl;
return c;
int main() {
double a = 1.0, b = 2.0;
double tolerance = 1e-6;
int maxIter = 100;
double root = regulaFalsi(a, b, tolerance, maxIter);
if (root != -1) {
cout << "\nApproximate root = " << root << endl;
cout << "f(root) = " << f(root) << endl;
return 0;
NEWTON RAPHSON METHOD
#include <iostream>
#include <cmath>
using namespace std;
double f(double x) {
return x*x*x - x - 2;
double df(double x) {
return 3*x*x - 1;
double newtonRaphson(double x0, double tolerance) {
double x1;
while (true) {
double fx = f(x0);
double dfx = df(x0);
if (dfx == 0) {
cout << "Derivative is zero. Cannot continue." << endl;
return x0;
x1 = x0 - fx / dfx;
if (fabs(x1 - x0) < tolerance)
break;
x0 = x1;
}
return x1;
int main() {
double x0, tolerance;
cout << "Enter initial guess: ";
cin >> x0;
cout << "Enter tolerance (e.g. 0.0001): ";
cin >> tolerance;
double root = newtonRaphson(x0, tolerance);
cout << "Approximate root = " << root << endl;
cout << "f(root) = " << f(root) << endl;
return 0;
TRAPEZOIDAL METHOD
#include <iostream>
#include <cmath>
#include <iomanip>
using namespace std;
double f(double x) {
return x * x; // Example: f(x) = x^2
int main() {
double a, b;
int n;
cout << "Enter lower limit a: ";
cin >> a;
cout << "Enter upper limit b: ";
cin >> b;
cout << "Enter number of subintervals n: ";
cin >> n;
double h = (b - a) / n;
double sum = f(a) + f(b);
for (int i = 1; i < n; i++) {
double x = a + i * h;
sum += 2 * f(x);
double result = (h / 2) * sum;
cout << fixed << setprecision(6);
cout << "\nApproximate integral = " << result << endl;
return 0;
}
Simpsons rule [1/3]
#include <iostream>
#include <cmath>
#include <iomanip>
using namespace std;
double f(double x) {
return x * x; // Example: f(x) = x^2
int main() {
double a, b;
int n;
cout << "Enter lower limit a: ";
cin >> a;
cout << "Enter upper limit b: ";
cin >> b;
cout << "Enter number of subintervals n (even number): ";
cin >> n;
if (n % 2 != 0) {
cout << "Number of subintervals must be even for Simpson's 1/3 rule." << endl;
return 1;
}
double h = (b - a) / n;
double sum = f(a) + f(b);
for (int i = 1; i < n; i++) {
double x = a + i * h;
if (i % 2 == 0) {
sum += 2 * f(x);
} else {
sum += 4 * f(x);
double integral = (h / 3) * sum;
cout << fixed << setprecision(10);
cout << "\nApproximate integral = " << integral << endl;
return 0;
Gauss quadrature method
#include <iostream>
#include <cmath>
using namespace std;
double f(double x) {
return x * x;
}
int main() {
double a, b;
cout << "Enter lower limit a: ";
cin >> a;
cout << "Enter upper limit b: ";
cin >> b;
double w1 = 1.0, w2 = 1.0;
double x1 = -1.0 / sqrt(3);
double x2 = 1.0 / sqrt(3);
double xm = 0.5 * (b + a);
double xr = 0.5 * (b - a);
double integral = xr * (w1 * f(xm + xr * x1) + w2 * f(xm + xr * x2));
cout << "Approximate integral using 2-point Gauss Quadrature = " << integral
<< endl;
return 0;
}
Gauss seidel method
#include <iostream>
#include <cmath>
using namespace std;
int main() {
int n;
cout << "Enter number of unknowns: ";
cin >> n;
double A[10][10], b[10], x[10];
cout << "Enter the coe icients of the matrix A row-wise:" << endl;
for (int i = 0; i < n; i++)
for (int j = 0; j < n; j++)
cin >> A[i][j];
cout << "Enter the constants vector b:" << endl;
for (int i = 0; i < n; i++) {
cin >> b[i];
x[i] = 0;
double tol;
cout << "Enter tolerance (e.g. 0.0001): ";
cin >> tol;
int maxIter = 100;
for (int iter = 0; iter < maxIter; iter++) {
double maxDi = 0;
for (int i = 0; i < n; i++) {
double sum = 0;
for (int j = 0; j < n; j++)
if (j != i)
sum += A[i][j] * x[j];
double x_new = (b[i] - sum) / A[i][i];
maxDi = max(maxDi , fabs(x_new - x[i]));
x[i] = x_new;
if (maxDi < tol)
break;
cout << "Solution:" << endl;
for (int i = 0; i < n; i++)
cout << "x" << i + 1 << " = " << x[i] << endl;
return 0;
Guass elimination
#include <iostream>
#include <cmath>
using namespace std;
int main() {
int n;
cout << "Enter number of unknowns: ";
cin >> n;
double A[10][10], b[10], x[10];
cout << "Enter the coe icients of matrix A row-wise:" << endl;
for (int i = 0; i < n; i++)
for (int j = 0; j < n; j++)
cin >> A[i][j];
cout << "Enter constants vector b:" << endl;
for (int i = 0; i < n; i++)
cin >> b[i];
for (int k = 0; k < n - 1; k++) {
for (int i = k + 1; i < n; i++) {
double factor = A[i][k] / A[k][k];
for (int j = k; j < n; j++)
A[i][j] -= factor * A[k][j];
b[i] -= factor * b[k];
for (int i = n - 1; i >= 0; i--) {
x[i] = b[i];
for (int j = i + 1; j < n; j++)
x[i] -= A[i][j] * x[j];
x[i] /= A[i][i];
}
cout << "Solution:" << endl;
for (int i = 0; i < n; i++)
cout << "x" << i + 1 << " = " << x[i] << endl;
return 0;
Runga kutta method
#include <iostream>
#include <iomanip>
using namespace std;
double f(double x, double y) {
return x + y; // Example: dy/dx = x + y
int main() {
double x0, y0, xn, h;
cout << "Enter initial x (x0): ";
cin >> x0;
cout << "Enter initial y (y0): ";
cin >> y0;
cout << "Enter final x (xn): ";
cin >> xn;
cout << "Enter step size h: ";
cin >> h;
double n = (xn - x0) / h;
double x = x0, y = y0;
cout << fixed << setprecision(6);
cout << "x\t\ty" << endl;
cout << x << "\t" << y << endl;
for (int i = 0; i < n; i++) {
double k1 = h * f(x, y);
double k2 = h * f(x + h/2.0, y + k1/2.0);
double k3 = h * f(x + h/2.0, y + k2/2.0);
double k4 = h * f(x + h, y + k3);
y = y + (k1 + 2*k2 + 2*k3 + k4) / 6.0;
x = x + h;
cout << x << "\t" << y << endl;
return 0;
Tylour series method
#include <iostream>
#include <iomanip>
using namespace std;
double f(double x, double y) {
return x + y;
double f_prime(double x, double y) {
return 1 + f(x, y);
int main() {
double x0, y0, xn, h;
cout << "Enter initial x (x0): ";
cin >> x0;
cout << "Enter initial y (y0): ";
cin >> y0;
cout << "Enter final x (xn): ";
cin >> xn;
cout << "Enter step size h: ";
cin >> h;
double n = (xn - x0) / h;
double x = x0, y = y0;
cout << fixed << setprecision(6);
cout << "x\t\ty" << endl;
cout << x << "\t" << y << endl;
for (int i = 0; i < n; i++) {
double y_new = y + h * f(x, y) + (h*h/2.0) * f_prime(x, y);
x = x + h;
y = y_new;
cout << x << "\t" << y << endl;
return 0;
Matrix multiplication
#include <iostream>
using namespace std;
int main() {
int r1, c1, r2, c2;
cout << "Enter rows and columns of first matrix: ";
cin >> r1 >> c1;
cout << "Enter rows and columns of second matrix: ";
cin >> r2 >> c2;
if (c1 != r2) {
cout << "Matrix multiplication not possible (columns of first must equal rows of
second)." << endl;
return 0;
int A[10][10], B[10][10], C[10][10] = {0};
cout << "Enter elements of first matrix:" << endl;
for (int i = 0; i < r1; i++)
for (int j = 0; j < c1; j++)
cin >> A[i][j];
cout << "Enter elements of second matrix:" << endl;
for (int i = 0; i < r2; i++)
for (int j = 0; j < c2; j++)
cin >> B[i][j];
// Matrix multiplication
for (int i = 0; i < r1; i++) {
for (int j = 0; j < c2; j++) {
C[i][j] = 0;
for (int k = 0; k < c1; k++) {
C[i][j] += A[i][k] * B[k][j];
}
cout << "Resultant matrix:" << endl;
for (int i = 0; i < r1; i++) {
for (int j = 0; j < c2; j++)
cout << C[i][j] << " ";
cout << endl;
return 0;
Matrix addition
#include <iostream>
using namespace std;
int main() {
int r, c;
cout << "Enter number of rows and columns: ";
cin >> r >> c;
int A[10][10], B[10][10], C[10][10];
cout << "Enter elements of first matrix:" << endl;
for (int i = 0; i < r; i++)
for (int j = 0; j < c; j++)
cin >> A[i][j];
cout << "Enter elements of second matrix:" << endl;
for (int i = 0; i < r; i++)
for (int j = 0; j < c; j++)
cin >> B[i][j];
// Matrix addition
for (int i = 0; i < r; i++)
for (int j = 0; j < c; j++)
C[i][j] = A[i][j] + B[i][j];
cout << "Resultant matrix after addition:" << endl;
for (int i = 0; i < r; i++) {
for (int j = 0; j < c; j++)
cout << C[i][j] << " ";
cout << endl;
return 0;
Transpose
#include <iostream>
using namespace std;
int main() {
int r, c;
cout << "Enter number of rows and columns: ";
cin >> r >> c;
int A[10][10], T[10][10];
cout << "Enter elements of the matrix:" << endl;
for (int i = 0; i < r; i++)
for (int j = 0; j < c; j++)
cin >> A[i][j];
for (int i = 0; i < r; i++)
for (int j = 0; j < c; j++)
T[j][i] = A[i][j];
cout << "Transpose of the matrix:" << endl;
for (int i = 0; i < c; i++) {
for (int j = 0; j < r; j++)
cout << T[i][j] << " ";
cout << endl;
return 0;
Substraction
#include <iostream>
using namespace std;
int main() {
int r, c;
cout << "Enter number of rows and columns: ";
cin >> r >> c;
int A[10][10], B[10][10], C[10][10];
cout << "Enter elements of first matrix:" << endl;
for (int i = 0; i < r; i++)
for (int j = 0; j < c; j++)
cin >> A[i][j];
cout << "Enter elements of second matrix:" << endl;
for (int i = 0; i < r; i++)
for (int j = 0; j < c; j++)
cin >> B[i][j];
// Matrix subtraction
for (int i = 0; i < r; i++)
for (int j = 0; j < c; j++)
C[i][j] = A[i][j] - B[i][j];
cout << "Resultant matrix after subtraction:" << endl;
for (int i = 0; i < r; i++) {
for (int j = 0; j < c; j++)
cout << C[i][j] << " ";
cout << endl;
return 0;
Vector (addition, dot product cross product)
#include <iostream>
using namespace std;
int main() {
int n;
cout << "Enter dimension of vectors: ";
cin >> n;
double A[10], B[10], sum[10];
cout << "Enter elements of first vector:" << endl;
for (int i = 0; i < n; i++)
cin >> A[i];
cout << "Enter elements of second vector:" << endl;
for (int i = 0; i < n; i++)
cin >> B[i];
for (int i = 0; i < n; i++)
sum[i] = A[i] + B[i];
cout << "Vector addition result: ";
for (int i = 0; i < n; i++)
cout << sum[i] << " ";
cout << endl;
double dot = 0;
for (int i = 0; i < n; i++)
dot += A[i] * B[i];
cout << "Dot product: " << dot << endl;
if (n == 3) {
double cross[3];
cross[0] = A[1]*B[2] - A[2]*B[1];
cross[1] = A[2]*B[0] - A[0]*B[2];
cross[2] = A[0]*B[1] - A[1]*B[0];
cout << "Cross product: ";
for (int i = 0; i < 3; i++)
cout << cross[i] << " ";
cout << endl;
} else {
cout << "Cross product is only defined for 3D vectors." << endl;
return 0;