그냥 하는 노트와 메모장

이분법(Bisection method) 본문

Algorithms

이분법(Bisection method)

coloredrabbit 2018. 6. 6. 16:06

* 이분법(Bisection method)
  1차원 데이터에 대해 [low, high] 구간을 지정하고, 원하는 데이터를 이분적으로 찾아가는 방식을 말한다.
여기서 low와 high는 사용자 정의 f()함수에 대해 f(low) > f(high)인 경우 low와 high를 swap 해주는 것이 더욱 깔끔하게 소스를 짤 수 있다.
  즉, low = 10, high = 15, f(low) = 2, f(high) = 1 인 경우, 결과적으로 f(low)가 f(high)보다 크기 때문에 구간 [low, high]를 [15, 10]로 잡는 것이다.

  방식은 iteration(순차 접근, for문)을 이용하는 방법과 while로 근사치를 추정해나가는 방식이 있다. 두 방식 중 안전한 방법은 iteration한 방법이다. 이유는 while로 근사치를 추정하는 방식은 프로그래머에게 직접 데이터를 설정을 강요하기 때문이며, 이 근사치가 실제 해와 다른 근사치로 설정될 수 있기 때문이다.

  iteration을 이용하는 경우엔 약 100번 정도 돌면 2^100 만큼 데이터를 나누고 찾아가게 되는데, 이정도면 소수값 답을 매우 근사하게 다가갈 수 있기 때문이다.

(※ 이 게시글에선 두 방식을 섞여 있다.)


* 문제

1. 2017 아주대학교 프로그래밍 경시대회 (APC) - 구분구적법 (Large, Small)

   (https://www.acmicpc.net/problem/14609, https://www.acmicpc.net/problem/14608)


  구간 low, high를 해당 x축의 가로 크기로 두고 이분법을 사용하여 적당한(오차의 범위에 벗어나지 않는) 값을 찾아간다.

문제에서 주어진 오차값을 1e-4를 이용하여 답을 도출해낼 수 있다.


#include <cstdio>
#include <cmath>
double diff = 0.0001, m = 1e-15;
double abs_v(double a) { return a < 0 ? -a : a; }
double ing(int c, int n, int a) {
	return c*(pow(a, n + 1.) / (double)(n + 1));
}
double gv(int c, int n, double a) {
	return c*pow(a, n);
}
int main() {
	int K, c[10], a, b,N,i,j;
	double integral=0, ans=-1,l,r,mid,fx,seg,off;
	scanf("%d", &K);
	for (i = 0; i <= K; i++)
		scanf("%d", &c[i]);
	scanf("%d%d%d", &a, &b, &N);
	for (i = 0; i <= K; i++)
		integral += ing(c[i], K - i, b) - ing(c[i], K - i, a);
	l = 0, seg = r = (double)(b - a)/N;
	while (l <= r) {
		mid = (l + r) / 2;
		fx = 0;
		for (i = 0; i < N; i++)
			for (j = 0; j <= K; j++)
				fx += gv(c[j], K - j, a + i*seg + mid)*seg;

		if (abs_v(integral - fx) <= diff) {
			ans = mid;
			break;
		}

		if (integral > fx) l = mid + m;
		else r = mid - m;
	}
	printf("%lf", ans);
	return 0;
}



2. ROOTS

   (https://algospot.com/judge/problem/read/ROOTS)


  재밌는 문제 !

  해를 구하기 위해 극점 추출을 한다. 극점은 f(x)에 대해 미분을 한 g(x)의 해가 곧 f(x)의 극점 위치가 된다. 이 점을 활용해야 한당

  2차 이하 방정식에 대해서 해는 근의 공식으로 구할 수 있고, 그 이상에 대해서는 미분을 하고, 극점을 추출한 다음, 극점 간의 데이터를 확인하여 해를 구하는 방식을 취한다.(재귀)

  아래 그림에서 e1,e2,e3가 극점이고 x0,x1,x2,x3가 해일 때, 양 끝과 극점 사이에 답이 있음을 알 수 있다. 단 f(e_i)*f(e_i+1) <= 0 이어야만 해가 존재한다. 극점 사이에 해가 존재한다면, 양 끝 극점의 f(x)값은 서로 부호가 달라야하기 때문이다.




  극점말고 양끝 구간에 대해서는 문제에 보이듯 절대값이 10 이하이므로 그 구간을 -11 ~ 11으로 대충 잡아주자.


  이제 이분법으로 답을 찾아가야 한다.

  이분법은 x1, x2에 대해 가운데 값을 구하고(mid = (x1+x2)/2) 그 값의 f(mid)를 구한다. 그 값이 0에 근접하다면 이를 답으로 내놓는다. 근접의 기준은 문제에서 주어진 값을 이용하며 최대한 근사한 치로 구하도록 한다.


- 코드

#include <iostream>
#include <vector>
#include <cmath>
#include <algorithm>
using namespace std;
typedef double DT;
const int L = 10;
const DT diff = 1e-9;
DT abs_v(DT a) { return a < 0 ? -a : a; }

vector<DT> getF0(vector<DT> coefs){
	vector<DT> ret;
	DT a=coefs[0],b=coefs[1], ans;
	if(coefs.size() == 3) {
		DT c=coefs[2], t;
		ans = -b/(2*a);
		if(b*b-4*a*c < 0) return ret; // 2-degree equation doesn't have any ans.
		t = sqrt(b*b-4*a*c)/(2*a);
		if(t == 0)
			ret.push_back(ans);
		else{
			ret.push_back(ans + t);
			ret.push_back(ans - t);
		}
	} else
		ret.push_back(-b/a);
	return ret;
}

vector<DT> getDifferential(vector<DT> coefs){
	vector<DT> ret;
	int n = coefs.size(),i;
	for(i=0;i<n-1;i++)
		ret.push_back((n-1-i)*coefs[i]);
	return ret;
}

DT getFx(vector<DT> coefs, DT x){
	DT ret = 0;
	for(int i=0;i<coefs.size();i++){
		ret *= x;
		ret += coefs[i];
	}
	return ret;
}

vector<DT> getAnswer(vector<DT> coefs){
	int n = (int)coefs.size()-1;
	if(n <= 2)
		return getF0(coefs);
	vector<DT> extremumP = getAnswer(getDifferential(coefs));
	extremumP.insert(extremumP.begin(), -L-1);
	extremumP.insert(extremumP.end(), L+1);
	sort(extremumP.begin(), extremumP.end());
	DT x1, x2, y1, y2, low, high,mid, y;
	vector<DT> ret;
	for(int i=0;i<extremumP.size()-1;i++){
		x1 = extremumP[i], x2 = extremumP[i+1];
		y1 = getFx(coefs, x1), y2 = getFx(coefs, x2);
		if(y1*y2 > 0) continue;
		low = y1 >= 0 ? x2 : x1;
		high = y1 >= 0 ? x1 : x2;
		do {
			mid = (low + high) / 2;
			y = getFx(coefs, mid);
			if(y < 0) low = mid;
			else high = mid;
		} while(abs_v(y) > diff);
		ret.push_back(mid);
	}
	return ret;
}

int main() {
	int C,n,i;
	DT d;
	vector<DT> coefs, ans;
	scanf("%d",&C);
	while(C--){
		scanf("%d\n",&n);
		for(i=0;i<=n;i++){
			scanf("%lf",&d);
			coefs.push_back(d);
		}
		ans = getAnswer(coefs);
		sort(ans.begin(),ans.end());
		for(i=0;i<ans.size();i++)
			printf("%.12lf ", ans[i]);
		puts("");
		coefs.clear();
	}
	return 0;
}



3. LOAN

   (https://algospot.com/judge/problem/read/LOAN)


  low와 high를 제대로 설정만 한다면 충분히 쉽게 접근할 수 있는 문제.

  high의 경우 C로 처리해야하는 최대값을 설정한다.


#include <cstdio>
int main() {
	int T, M;
	double N, P, lo,hi,mid,i,TN;
	scanf("%d", &T);
	while (T--) {
		scanf("%lf %d %lf", &N, &M, &P);
		lo = 0;
		hi = N*(1.0 + P/1200.);

		for(int it=0;it<100;it++){
			TN = N;
			mid = (lo + hi) / 2.;
			for (i = 0; i < M; i++) {
				TN *= (1.0 + P / 1200.);
				TN -= mid;
			}
			if (TN <= 0) hi = mid;
			else lo = mid;
		}
		printf("%.12lf\n", hi);
	}
	return 0;
}


Comments