Taesan Kim 2024. 9. 4. 23:17

Practice[Integral]

 

Intro

Today we will calculate integral by c programming.

 

Header File

/*
================================================
Handong Global University
------------------------------------------------
Name:		Taesan Kim
ID:			22300203
Create:		2024.07.19
Modifier:	2024.07.19
------------------------------------------------
전역변수를 이용하여 적분을 지원한다.
================================================
*/

#ifndef MYINTEGRAL_H
#define MYINTEGRAL_H
#include <stdio.h>
#include <math.h>

/**
* @breif
* 
* @f: 적분할 함수
* @a: 적분 시작값
* @b: 적분 끝값
* @n: 적분구간 분할개수
* 
* @return: 적분값 
*/

double rect(double(*f)(double), double a, double b, int n);

double trpz(double (*f)(double), double a, double b, int n);

double smps(double (*f) (double), double a, double b, int n);

double adapt(double(*f)(double), double a, double b, double tolerance);

#endif //MYINTEGRAL_H

 

Source File

#include "include/myIntegral.h"

double rect(double(*f)(double), double a, double b, int n)
{
	int i;

	double dt = (b - a) / n, sum = 0;
	for (i = 0; i < n; i++) sum += f(a + dt / 2 + i * dt) * dt;
	return sum;
}


double trpz(double (*f)(double), double a, double b, int n)
{
	int i;

	double dt = (b - a) / n, sum = 0;
	for (i = 0; i < n; i++) sum += ((f(a + i * dt) + f(a + (i + 1) * dt)) / 2) * dt;

	return sum;
}


double smps(double (*f) (double), double a, double b, int n)
{
	int i;
	double dt = (b - a) / n, sum = 0;

	for (i = 0; i < n; i++)sum += (f(a+i*dt) + 4 * f(a + dt / 2 + i*dt) + f(i*dt + a + dt)) * dt / 6;

	return sum;
}


double adapt(double(*f)(double), double a, double b, double tolerance)
{
	double x = smps(f, a, b, 4);
	double x2 = smps(f, a, b, 2);

	if (fabs(x - x2) > tolerance)
	{
		double tegral_1 = adapt(f, a, (a + b) / 2, tolerance);
		double tegral_2 = adapt(f, (a + b) / 2, b, tolerance);
		return tegral_1 + tegral_2;
	}
	else if (fabs(x - x2) <= tolerance)
	{
		return x;
	}
}

 

Main

#include "../include/myIntegral.h"

#define f(x) (4*x*x + 6*x + 9)
//매크로는 단순한 대체이므로 함수 포인터로 사용할 수 없다.

double f_function(double x)
{
	return powf(x, 5.0) - powf(x, 4.0);
}

int main() {

	printf("rectangular Rule : %f\n\n", rect(f_function, 0.0, 10.0, 10));
	printf("trapezoid Rule : %f\n\n", trpz(f_function, 0.0, 10.0, 10));
	printf("smpson Rule : %f\n\n", smps(f_function, 0.0, 10.0, 10));
	printf("adapted Rule : %f\n\n", adapt(f_function, 0.0, 10.0, 0.01));
	

	double input[11] = { 0.0, 1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0, 9.0, 10.0 };
	double dt = 0.1;
	double n;

	double integral;
	double pre_val = 0.0;

	for (int i = 0; i < 10; i++)
	{
		n = (input[i+1] - input[i]) / dt;
		integral = smps(f_function, input[i], input[i + 1], n);
		pre_val += integral;
		printf("integral[%d] : %f\n", i, pre_val);
	}

	return 0;
}