// zeros.cpp : Defines the entry point for the console application.
//

#include "stdafx.h"
#include "math.h"
#include "conio.h"

// definição de um ponteiro de função real
typedef float (*real_function)(float x);
float zero=0.0;

// nossa função
float f(float x)
{
	return pow(sin(x),2)-x*sqrt(2)/2-1;
}

// derivade aproximativa de uma funçao
#define H 0.0078125
float derive(real_function f, float a)
{
	return (f(a+H)-f(a-H))/(2*H);
}

// determinação de uma raiz da função f no intervalo [a,b] com a tolerância tol
float raiz_bissecao(real_function f, float a, float b, float tol)
{
	// f não troca de sinal entre a e b, a bisseção não tem sentido, resultado NaN
	if(f(a)*f(b)>0) return zero/zero;
	// f é já egal a 0 em a (resp. b), resultado a (resp. b)
	if(f(a)==0) return a;
	if(f(b)==0) return b;

	// enquanto o tamanho do intervalo [a,b] é maior que a tolerância, continuamos
	while(fabs(a-b)>tol) {
		// m, media de a e b
		float m=(a+b)/2;

		cprintf("Intervalo [%f, %f], ",a,b);
		// se f(m) é egal a zero, m é raiz, resultado m
		if(f(m)==0) return m;
		// se a media de a e b é egal a a ou b, significa que o intervalo é pequena de mais 
		// para ser representado com o computador
		if(m==a || m==b) {cprintf("tolerancia pequena de mais\r\n");return m;}
		// se f(a) e f(b) têm mesma valor, a tolerancia é pequena de mais para calcular f
		if(f(a)==f(b)) {cprintf("precisão pequena de mais para calcular f\r\n");return m;}
		// se f troca de sinal entre a e m, o novo intervalo é [a,m]
		if(f(a)*f(m)<0) b=m;
		// senão, o novo intervalo é [m,b]
		else a=m;

		cprintf("media: %f\r\n",m);
	}
	cprintf("conseguimos com a tolerencia pedida\r\n");
	cprintf("Intervalo [%f, %f]\r\n",a,b);
	return (a+b)/2;
}

float raiz_secante(real_function f, float a, float b, int iter)
{
	int i;
	for(i=0;i<iter;i++) {
		a=(f(a)*b-f(b)*a)/(f(a)-f(b));
		cprintf("passo:%d: %f\r\n",i+1,a);
	}
	return a;
}

float raiz_falsaposicao(real_function f, float a, float b, int iter)
{
	if(f(a)*f(b)>0) return zero/zero;
	else {
		int i;
		float x=zero/zero;
	
		for(i=0;i<iter;i++) {
			x=(f(a)*b-f(b)*a)/(f(a)-f(b));
			if(f(a)*f(x)<0) b=x;
			else a=x;
			cprintf("passo:%d: %f\r\n",i+1,x);
		}
		return x;
	}
}

float raiz_newtonraphson(real_function f, float a, int iter)
{
	int i;
	for(i=0;i<iter;i++) {
		a=a-f(a)/derive(f,a);
		cprintf("passo:%d: %f\r\n",i+1,a);
	}
	return a;
}

int main()
{
	float x;
	x=raiz_bissecao(f,-1,0,H);
	cprintf("result bissecao:%f, valor em f:%f\r\n",x,f(x));
	x=raiz_secante(f,-1,0,8);
	cprintf("result secante:%f, valor em f:%f\r\n",x,f(x));
	x=raiz_falsaposicao(f,-1,0,8);
	cprintf("result falsa posicao:%f, valor em f:%f\r\n",x,f(x));
	x=raiz_newtonraphson(f,-0,8);
	cprintf("result Newton Raphson:%f, valor em f:%f\r\n",x,f(x));
	while(!getch());
}