Используя три самостоятельно реализованные функции iPow(x, n), Ln(x) и Exp(x), я могу вычислить fPow(x, a), x и существо , удваивающее . Ни одна из функций ниже не использует библиотечные функции, а только итерацию.
Некоторые пояснения о реализованных функциях:
(1) iPow(x, n): x равен double, n равен int. Это простая итерация, так как n - целое число.
(2) Ln(x): эта функция использует итерацию ряда Тейлора. Ряд, используемый в итерации, равен Σ (from int i = 0 to n) {(1 / (2 * i + 1)) * ((x - 1) / (x + 1)) ^ (2 * n + 1)}. Символ ^ обозначает степенную функцию Pow(x, n), реализованную в 1-й функции, которая использует простую итерацию.
(3) Exp(x): Эта функция снова использует итерацию серии Тейлора. Ряд, используемый в итерации, равен Σ (from int i = 0 to n) {x^i / i!}. Здесь ^ обозначает степенную функцию, но она не вычисляется путем вызова 1-й Pow(x, n) функции; вместо этого он реализован в 3-й функции одновременно с факториалом, используя d *= x / i. Я чувствовал , что мне пришлось использовать этот трюк , потому что в этой функции итерация делает еще несколько шагов относительно других функций, а факториал (i!) переполняется большую часть времени. Чтобы убедиться, что итерация не переполняется, степенная функция в этой части повторяется одновременно с факториалом. Таким образом, я преодолел переполнение.
(4) fPow(x, a): x и a оба являются двойными . Эта функция не делает ничего, кроме как просто вызывает другие три функции, реализованные выше. Основная идея этой функции зависит от некоторого исчисления: fPow(x, a) = Exp(a * Ln(x)). И теперь у меня есть все функции iPow, Ln и Exp с итерацией.
n.b. Я использовал constant MAX_DELTA_DOUBLE, чтобы решить, на каком шаге остановить итерацию. Я установил его на 1.0E-15, что кажется разумным для двойников. Таким образом, итерация останавливается, если (delta < MAX_DELTA_DOUBLE) Если вам нужна дополнительная точность, вы можете использовать long double и уменьшить постоянное значение для MAX_DELTA_DOUBLE, например, до 1.0E-18 (1,0E-18 будет минимальным).
Вот код, который работает для меня.
#define MAX_DELTA_DOUBLE 1.0E-15
#define EULERS_NUMBER 2.718281828459045
double MathAbs_Double (double x) {
return ((x >= 0) ? x : -x);
}
int MathAbs_Int (int x) {
return ((x >= 0) ? x : -x);
}
double MathPow_Double_Int(double x, int n) {
double ret;
if ((x == 1.0) || (n == 1)) {
ret = x;
} else if (n < 0) {
ret = 1.0 / MathPow_Double_Int(x, -n);
} else {
ret = 1.0;
while (n--) {
ret *= x;
}
}
return (ret);
}
double MathLn_Double(double x) {
double ret = 0.0, d;
if (x > 0) {
int n = 0;
do {
int a = 2 * n + 1;
d = (1.0 / a) * MathPow_Double_Int((x - 1) / (x + 1), a);
ret += d;
n++;
} while (MathAbs_Double(d) > MAX_DELTA_DOUBLE);
} else {
printf("\nerror: x < 0 in ln(x)\n");
exit(-1);
}
return (ret * 2);
}
double MathExp_Double(double x) {
double ret;
if (x == 1.0) {
ret = EULERS_NUMBER;
} else if (x < 0) {
ret = 1.0 / MathExp_Double(-x);
} else {
int n = 2;
double d;
ret = 1.0 + x;
do {
d = x;
for (int i = 2; i <= n; i++) {
d *= x / i;
}
ret += d;
n++;
} while (d > MAX_DELTA_DOUBLE);
}
return (ret);
}
double MathPow_Double_Double(double x, double a) {
double ret;
if ((x == 1.0) || (a == 1.0)) {
ret = x;
} else if (a < 0) {
ret = 1.0 / MathPow_Double_Double(x, -a);
} else {
ret = MathExp_Double(a * MathLn_Double(x));
}
return (ret);
}