Peter.alexander.au

#include <iostream>
#include <vector>
#include <math.h>
#include <complex>
#include <algorithm>
#include <iterator>

using namespace std;

typedef complex<double> C;

C fast_exp(C x, int y)
{
	if (y == 0) return 1.0;
	if (y == 1) return x;
	int h = y / 2;
	return fast_exp(x, h) * fast_exp(x, y-h);
}

typedef vector<double> Poly;

C eval_poly(const Poly& p, C x)
{
	C y = 0.0;
	for (int i = 0; i < p.size(); ++i)
		y += p[i] * fast_exp(x, i);
	return y;
}

int main()
{
	int n;
	cin >> n;
	Poly p(n, 0.0);
	for (int i = 0; i < n; ++i)
		cin >> p[i];
	
	vector<C> q(n-1, C(0.4, 0.9));
	for (int i = 0; i < n-1; ++i)
		q[i] = fast_exp(C(0.4, 0.9), i);

	for (int k = 0; k < 10; ++k)
	{
		for (int i = 0; i < n-1; ++i)
		{
			C denom(1.0, 0.0);
			for (int j = 0; j < n-1; ++j)
				if (i != j)
					denom *= q[i] - q[j];
				
			q[i] -= eval_poly(p, q[i]) / denom;
		}
	}
	
	copy(q.begin(), q.end(), ostream_iterator<C>(cout));
	cout << endl;
}

Start a discussion with Peter.alexander.au

Talk pages are where people discuss how to make content on Wikipedia the best that it can be. Start a new discussion to connect and collaborate with Peter.alexander.au. What you say here will be public for others to see.

Start a discussion