Новые математические возможности Java: Часть 1. Вещественные числа

В этой серии из двух статей Эллиотта Расти Гарольда рассматриваются "старые-новые" возможности хорошо известного класса java.lang.Math. Первая часть посвящена исключительно математическим функциям, а во второй речь пойдет о функциях, созданных для работы с числами с плавающей точкой.

Элиот Гарольд, адъюнкт-профессор, компания

Элиот Расти Гарольд (Elliotte Rusty Harold) родом из Нового Орлеана, куда он периодически возвращается в поисках достойной тарелки супа из стручков бамии. Однако он постоянно проживает в Проспект Хайт, пригороде Бруклина, со своей женой Бет и кошками Шарм (названной в честь кварка) и Марджори (названной в честь тещи). Он занимает должность адъюнкт-профессора теории вычислительных машин и систем в Политехническом университете, где он преподает Java и объектно-ориентированное программирование. Его сайт Cafe au Lait стал одним из самых популярных независимых ресурсов Интернет, посвященных Java, а сопутствующий сайт Cafe con Leche стал одним из самых посещаемых сайтов об XML. Его самая свежая книга – Java I/O, 2-ое издание. В последнее время он работает над XOM, прикладным программным интерфейсом для обработки XML, а также виртуальной XPath машиной Jaxen и Jester инструментальным средством для тестирования эффективности текста. Вы можете связаться с Элиотом по адресу elharo@metalab.unc.edu.



07.07.2010

Иногда вам кажется, что вы настолько хорошо знакомы с некоторым классом, что перестаете обращать на него внимание. Например, станете ли вы читать Javadoc по классу, который знаете настолько, что сами готовы написать по нему документацию, а также учитывая возможности автодополнения в Eclipse? Я думал точно так же в отношении класса java.lang.Math, который, как мне казалось, я знал досконально. Представьте себе мое недавнее удивление, когда, обратившись к его Javadoc впервые за последние пять лет, я обнаружил, что он увеличился как минимум вдвое и содержит 20 новых методов, о которых я никогда не слышал. Очевидно, что пришло время освежить свое представление об этом классе.

В спецификации Java™ 5 в класс java.lang.Math (а также в его порочный брат-близнец java.lang.StrictMath) было добавлено 10 новых методов, а в Java 6 - еще 10. Эта статья будет посвящена чистым математическим функциям, таким как log10 и cosh, а в следующей статье будут рассматриваться функции, предназначенные для работы с типами данных с плавающей точкой, а не просто с абстрактными числами.

Существует важное различие между абстрактными вещественными числами, такими как π или 0.2, и типом данных double в Java. Во-первых, платонически-идеальное представление вещественных чисел является бесконечным, в то время как представление в Java ограничено числом бит. Это очень важно при работе с очень большими или очень малыми числами. Например, число 2,000,000,001 (два миллиарда один) помещается в формат int, но не float. Ближайшее к нему число типа float – это 2.0E09, т.е. два миллиарда. Тип double имеет преимущество, так как его представление длиннее (именно поэтому double практически всегда используется вместо float), однако все равно, как правило, существуют ограничения на точность представления.

Второй недостаток машинной арифметики в Java и других языках заключается в использовании двоичного, а не десятичного представления чисел. Такие дроби как, например, 1/5 и 7/50, которые точно представимы в десятичной нотации (0.2 и 0.14 соответственно) оказываются периодическими в двоичном коде. Эта проблема аналогична трудностям в представлении периодических дробей, например 1/3 в десятичном коде (0.3333333...). Любая дробь, в которой единственными простыми множителями делителя являются 2 и 5, всегда точно представима в системе счисления с основанием 10. В системе с основанием 2 таковыми дробями являются только те, у которых делитель представляет собой степень двойки, например 1/2, 1/4, 1/8, 1/16 и т.д.

Подобные неточности в представлении чисел являются одной из главных причин необходимости специального класса для математических операций. Разумеется, вы можете реализовать тригонометрические и другие функции при помощи разложения в ряд Тейлора, используя лишь простой цикл, а также операции сложения и умножения, как показано в листинге 1.

Листинг 1. Вычисление синуса путем разложения в ряд Тейлора
public class SineTaylor {

    public static void main(String[] args) {
        for (double angle = 0; angle <= 4*Math.PI; angle += Math.PI/8) {
            System.out.println(degrees(angle) + "\t" + taylorSeriesSine(angle) 
                    + "\t" + Math.sin(angle));
        }
    }
    
    public static double degrees(double radians) {
        return 180 * radians/ Math.PI;
    }
    
    public static double taylorSeriesSine(double radians) {       
        double sine = 0;
        int sign = 1;
        for (int i = 1; i < 40; i+=2) {
            sine += Math.pow(radians, i) * sign / factorial(i);
            sign *= -1;
        }
        return sine;
    }

    private static double factorial(int i) {
        double result = 1;
        for (int j = 2; j <= i; j++)  {
            result *= j;
        }
        return result;
    }
}

Для небольших углов этот подход работает очень неплохо. Расхождение между вычисленным и точным значением в худшем случае начинается лишь в последней десятичной цифре.

0.0       0.0       0.0
22.5       0.3826834323650897       0.3826834323650898
45.0       0.7071067811865475       0.7071067811865475
67.5       0.923879532511287        0.9238795325112867
90.0       1.0000000000000002       1.0

Однако с увеличением угла ошибки начинают накапливаться, поэтому метод начинает показывать худшие результаты.

630.0000000000003       -1.0000001371557132       -1.0
652.5000000000005       -0.9238801080153761       -0.9238795325112841
675.0000000000005       -0.7071090807463408       -0.7071067811865422
697.5000000000006       -0.3826922100671368       -0.3826834323650824

На самом деле результаты вычислений путем разложения в ряд Тейлора оказались даже несколько лучше, чем я ожидал. Тем не менее при увеличении угла до 360, 720 градусов (2 и 4 π радиан) и более, для обеспечения точности необходимо считать все больше членов в ряду. Класс java.lang.Math использует более сложные алгоритмы, позволяющие этого избежать.

Кроме того, разложение в ряд Тейлора очень неэффективно по сравнению с функцией sin для вычисления синуса, встроенной в современные процессоры. Вычисление синуса и других функций достаточно быстро и с необходимой точностью требует сложных алгоритмов, специально спроектированных таким образом, чтобы не допускать накопления погрешностей. Зачастую эти алгоритмы реализуются непосредственно в электронных схемах, что еще повышает производительность. В частности, практически каждый чип x86, произведенный за последние 10 лет, включал в себя реализацию функций синуса и косинуса, которые могут непосредственно вызываться виртуальной машиной x86 вместо программного вычисления на основе примитивных арифметических операций. Благодаря использованию этих команд в HotSpot удалось существенно ускорить вычисление тригонометрических функций.

Прямоугольные треугольники и Евклидовы нормы

Любой учащийся старших классов школы знает теорему Пифагора: квадрат гипотенузы прямоугольного треугольника равен сумме квадратов катетов (c2 = a2 + b2).

Те из вас, кто изучал физику и высшую математику в ВУЗах, знают, что это соотношение встречается далеко не только в прямоугольных треугольниках. Например, это также квадрат Евклидовой нормы R2, длина двухмерного вектора, часть неравенства треугольника и т.д. На самом деле все это представляет собой примеры различных взглядов на одно и то же соотношение, но смысл в том, что теорема Пифагора значительно важнее, чем может показаться на первый взгляд.

В Java 5 был добавлен метод Math.hypot, выполняющий именно эту операцию (это хороший пример того, как может быть полезна данная библиотека). Наивный вариант вычисления суммы квадратов мог бы выглядеть примерно так:

public static double hypot(double x, double y){
  return Math.sqrt (x*x + y*y);
}

В реальности же код несколько более сложен (листинг 2). Первое, что бросается в глаза – это то, что метод реализован на языке С для максимальной производительности. Кроме того, видно, что были предприняты все меры для минимизации погрешности вычисления. Более того, в зависимости от относительного размера переменных x и y выбираются разные алгоритмы.

Листинг 2. Использующаяся на практике реализация метода Math.hypot
/*
 * ====================================================
 * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
 *
 * Developed at SunSoft, a Sun Microsystems, Inc. business.
 * Permission to use, copy, modify, and distribute this
 * software is freely granted, provided that this notice 
 * is preserved.
 * ====================================================
 */

#include "fdlibm.h"

#ifdef __STDC__
       double __ieee754_hypot(double x, double y)
#else
       double __ieee754_hypot(x,y)
       double x, y;
#endif
{
       double a=x,b=y,t1,t2,y1,y2,w;
       int j,k,ha,hb;

       ha = __HI(x)&0x7fffffff;       /* high word of  x */
       hb = __HI(y)&0x7fffffff;       /* high word of  y */
       if(hb > ha) {a=y;b=x;j=ha; ha=hb;hb=j;} else {a=x;b=y;}
       __HI(a) = ha;       /* a <- |a| */
       __HI(b) = hb;       /* b <- |b| */
       if((ha-hb)>0x3c00000) {return a+b;} /* x/y > 2**60 */
       k=0;
       if(ha > 0x5f300000) {       /* a>2**500 */
          if(ha >= 0x7ff00000) {       /* Inf or NaN */
              w = a+b;                     /* for sNaN */
              if(((ha&0xfffff)|__LO(a))==0) w = a;
              if(((hb^0x7ff00000)|__LO(b))==0) w = b;
              return w;
          }
          /* scale a and b by 2**-600 */
          ha -= 0x25800000; hb -= 0x25800000;       k += 600;
          __HI(a) = ha;
          __HI(b) = hb;
       }
       if(hb < 0x20b00000) {       /* b < 2**-500 */
           if(hb <= 0x000fffff) {       /* subnormal b or 0 */       
              if((hb|(__LO(b)))==0) return a;
              t1=0;
              __HI(t1) = 0x7fd00000;       /* t1=2^1022 */
              b *= t1;
              a *= t1;
              k -= 1022;
           } else {              /* scale a and b by 2^600 */
               ha += 0x25800000;        /* a *= 2^600 */
              hb += 0x25800000;       /* b *= 2^600 */
              k -= 600;
                 __HI(a) = ha;
                 __HI(b) = hb;
           }
       }
    /* если a и b среднего размера */
       w = a-b;
       if (w>b) {
           t1 = 0;
           __HI(t1) = ha;
           t2 = a-t1;
           w  = sqrt(t1*t1-(b*(-b)-t2*(a+t1)));
       } else {
           a  = a+a;
           y1 = 0;
           __HI(y1) = hb;
           y2 = b - y1;
           t1 = 0;
           __HI(t1) = ha+0x00100000;
           t2 = a - t1;
           w  = sqrt(t1*y1-(w*(-w)-(t1*y2+t2*b)));
       }
       if(k!=0) {
           t1 = 1.0;
           __HI(t1) += (k<<20);
           return t1*w;
       } else return w;
}

На самом деле, будете вы использовать именно эту реализацию или одну из нескольких похожих, зависит от используемой JVM для вашей платформы. Однако, скорее всего, именно эта функция будет использоваться при работе со стандартной JDK от Sun (другие поставщики JDK могут попытаться разработать более быстродействующую реализацию этой функции).

Этот код (а также большая часть нативного кода реализации математических функций из библиотеки Java Development Library от Sun) был позаимствована из открытой библиотеки fdlibm, созданной в Sun примерно 15 лет назад. Эта библиотека была специально спроектирована для реализации стандарта IEE754, а также обеспечения точности вычислений, возможно даже ценой некоторой потери производительности.


Логарифмы по основанию 10

Значением логарифма является число, в степень которого должно быть возведено основание, чтобы получить заданный результат. Другими словами, функция логарифма является обратной по отношению к Math.pow(). В инженерных расчетах чаще встречаются логарифмы по основанию 10, а логарифмы по основанию e (натуральные логарифмы) – при вычислении сложных процентов, а также в большом числе научных и математических приложений. Логарифмы по основанию 2, как правило, приходится вычислять при анализе алгоритмов.

Метод для вычисления натурального логарифма присутствовал в классе Math начиная с Java 1.0. Получая на вход параметр х, данный метод возвращает степень, в которую должно быть возведено число e, чтобы получить значение х. К сожалению, в Java аналогично C, Fortran и Basic, функция натурального алгоритма ошибочно названа log(), в то время как во всех книгах по математике, которые мне встречались, так обозначается логарифм по основанию 10, тогда как стандартными обозначениями для натурального логарифма и логарифма по основанию 2 являются ln и lg соответственно. Теперь эту ошибку исправлять уже поздно, но в Java 5 появился новый метод log10(), вычисляющий значение логарифма по основанию 10.

В листинге 3 показана простая программа, печатающая значения логарифмов по основанию 2, 10 и e для всех целых чисел от 1 до 100.

Листинг 3. Вычисление логарифмов по различным основаниям для целых чисел от 1 до 100
public class Logarithms {
    public static void main(String[] args) {
        for (int i = 1; i <= 100; i++) {
            System.out.println(i + "\t" +
                              Math.log10(i) + "\t" +
                              Math.log(i) + "\t"  +
                              lg(i));
        }
    }
    
    public static double lg(double x) {
        return Math.log(x)/Math.log(2.0);
    }
}

Ниже приведены первые 10 результатов.

1    0.0                                0.0                                 0.0
2    0.3010299956639812  0.6931471805599453   1.0
3    0.47712125471966244 1.0986122886681096  1.584962500721156
4    0.6020599913279624  1.3862943611198906   2.0
5    0.6989700043360189  1.6094379124341003   2.321928094887362
6    0.7781512503836436  1.791759469228055     2.584962500721156
7    0.8450980400142568  1.9459101490553132   2.807354922057604
8    0.9030899869919435  2.0794415416798357   3.0
9    0.9542425094393249  2.1972245773362196   3.1699250014423126
10  1.0                                2.302585092994046     3.3219280948873626

Метод Math.log10() обладает обычными свойствами функции логарифмирования. В частности, при попытке взять логарифм нуля или отрицательного числа будет возвращено значение NaN.


Кубические корни

Я не могу припомнить ни одного случая в своей жизни, когда мне пришлось бы брать кубический корень, хотя я и отношусь к тем немногим людям, которые регулярно сталкиваются с алгеброй и геометрией в повседневной жизни, не говоря уже об отдельных случаях использования математического анализа, дифференциальных уравнений и даже абстрактной алгебры. Таким образом, полезность этого метода от меня ускользает. Тем не менее, если вам все же потребовалось взять кубический корень, то в Java 5 есть метод Math.cbrt(). Пример его использования для чисел в диапазоне от -5 до 5 приведен в листинге 4.

Листинг 4. Вычисление кубических корней для чисел от -5 до 5
public class CubeRoots {
    public static void main(String[] args) {
        for (int i = -5; i <= 5; i++) {
            System.out.println(Math.cbrt(i));
        }
    }
}

Вот как выглядит результат:

-1.709975946676697
-1.5874010519681996
-1.4422495703074083
-1.2599210498948732
-1.0
0.0
1.0
1.2599210498948732
1.4422495703074083
1.5874010519681996
1.709975946676697

Этот пример демонстрирует одну отличительную особенность кубических корней по сравнению с квадратными: каждое вещественное число имеет ровно один кубический корень. Данный метод возвращает NaN только в том случае, когда NaN был передан в качестве параметра.


Гиперболические тригонометрические функции

Гиперболические тригонометрические функции имеют такое же отношение к гиперболам, как обычные тригонометрические функции к окружностям. Представьте себе, что значения этих функций нанесены на Декартову плоскость для всех возможных значений аргумента t.

x = r cos(t)
y = r sin(t)

В результате у вас получится окружность радиуса r. Теперь представьте, что вместо обычных синуса и косинуса используются функции sinh и cosh.

x = r cosh(t)
y = r sinh(t)

В этом случае вы нарисуете равнобочную гиперболу, в которой точкой, наиболее близкой к началу координат, является r.

О гиперболических функциях можно рассуждать следующим образом: известно, что sin(x) равен (eix - e-ix, а cos(x) = (eix + e-ix)/2. Функции sinh и cosh можно получить, убрав мнимую единицу из этих формул, т.е. sinh(x) = (ex - e-x)/2, а cosh(x) = (ex + e-x)/2.

В Java 5 были добавлены все три метода: Math.cosh(), Math.sinh() и Math.tanh(). При этом в состав класса пока не были включены обратные гиперболические функции: acosh, asinh и atanh.

График функции cosh(z) имеет вид свисающего каната, закрепленного с обоих концов. Он известен под именем цепная линия (catenary). В листинге 5 показан текст простой программы, рисующей цепную линию при помощи функции Math.cosh.

Листинг 5. Рисование цепной линии при помощи функции Math.cosh()
import java.awt.*;

public class Catenary extends Frame {

    private static final int WIDTH = 200;
    private static final int HEIGHT = 200;
    private static final double MIN_X = -3.0;
    private static final double MAX_X = 3.0;
    private static final double MAX_Y = 8.0;

    private Polygon catenary = new Polygon();

    public Catenary(String title) {
        super(title);
        setSize(WIDTH, HEIGHT);
        for (double x = MIN_X; x <= MAX_X; x += 0.1) {
            double y = Math.cosh(x);
            int scaledX = (int) (x * WIDTH/(MAX_X - MIN_X) + WIDTH/2.0);
            int scaledY = (int) (y * HEIGHT/MAX_Y);
            // В компьютерной графике ось y направлена вниз, а не вверх
            // (как в Декартовых координатах), поэтому необходимо изменить знак
            scaledY = HEIGHT - scaledY;
            catenary.addPoint(scaledX, scaledY);
        }
    }

    public static void main(String[] args) {
        Frame f = new Catenary("Catenary");
        f.setVisible(true);
    }

    public void paint(Graphics g) {
        g.drawPolygon(catenary);
    }

}

График показан на рисунке 1.

Рисунок 1. График цепной линии на Декартовой плоскости
График цепной линии на Декартовой плоскости

Функции sinh, cosh и tanh фигурируют в различных вычислениях в общей и специальной теориях относительности.


Знак

Метод Math.signum преобразует положительные числа в 1.0, отрицательные – в -1.0, а ноль оставляет без изменений. По сути, она просто возвращает знак числа. Это может быть полезно при реализации интерфейса Comparable.

Существуют две версии этого метода: для типов float и double. Метод может показаться тривиальным, однако необходим для корректной обработки специальных случаев, встречающихся при вычислениях с плавающей точкой, например NaN, а также положительные и отрицательные нули. NaN считается равным нулю, а знаком положительного и отрицательного нуля считается положительный и отрицательный ноль соответственно. Например, допустим, что вам потребовалось реализовать данный метод, и вы выбрали наивное решение, подобное показанному в листинге 6.

Листинг 6. Ошибочная реализация метода Math.signum
public static double signum(double x) {
  if (x == 0.0) return 0;
  else if (x < 0.0) return -1.0;
  else return 1.0;
}

Во-первых, этот метод преобразует положительные нули в отрицательные (отрицательные нули могут выглядеть несколько странно, однако они являются неотъемлемой частью спецификации IEEE 754). Во-вторых, в соответствии с этим методом значение NaN является положительным. Поэтому на практике используется реализация, показанная в листинге 7, которая несколько сложнее, но зато корректно обрабатывает специальные случаи.

Листинг 7. Корректная реализация метода Math.signum, использующаяся на практике
public static double signum(double d) {
    return (d == 0.0 || isNaN(d))?d:copySign(1.0, d);
}

public static double copySign(double magnitude, double sign) {
    return rawCopySign(magnitude, (isNaN(sign)?1.0d:sign));
}

public static double rawCopySign(double magnitude, double sign) {
    return Double.longBitsToDouble((Double.doubleToRawLongBits(sign) &
                                   (DoubleConsts.SIGN_BIT_MASK)) |
                                   (Double.doubleToRawLongBits(magnitude) &
                                   (DoubleConsts.EXP_BIT_MASK |
                                    DoubleConsts.SIGNIF_BIT_MASK)));
}

Добивайтесь большего меньшими усилиями

Как говорится, наиболее эффективным является код, который попросту не пришлось писать. Поэтому не стоит делать то, что уже сделано экспертами до вас. Код, использующий как новые, так и старые методы класса java.lang.Math, будет быстрее, эффективнее и точнее, чем тот, что вы напишете сами. Просто используйте его.

Ресурсы

Научиться

Получить продукты и технологии

  • fdlibm: данную библиотеку, написанную на языке С и реализующую арифметику с плавающей точкой в соответствии со стандартом IEEE 745, можно загрузить из репозитория математического программного обеспечения Netlib. (EN)
  • OpenJDK: взгляните на исходный код математических классов в этой открытой реализации Java SE. (EN)

Обсудить

Комментарии

developerWorks: Войти

Обязательные поля отмечены звездочкой (*).


Нужен IBM ID?
Забыли Ваш IBM ID?


Забыли Ваш пароль?
Изменить пароль

Нажимая Отправить, Вы принимаете Условия использования developerWorks.

 


Профиль создается, когда вы первый раз заходите в developerWorks. Информация в вашем профиле (имя, страна / регион, название компании) отображается для всех пользователей и будет сопровождать любой опубликованный вами контент пока вы специально не укажите скрыть название вашей компании. Вы можете обновить ваш IBM аккаунт в любое время.

Вся введенная информация защищена.

Выберите имя, которое будет отображаться на экране



При первом входе в developerWorks для Вас будет создан профиль и Вам нужно будет выбрать Отображаемое имя. Оно будет выводиться рядом с контентом, опубликованным Вами в developerWorks.

Отображаемое имя должно иметь длину от 3 символов до 31 символа. Ваше Имя в системе должно быть уникальным. В качестве имени по соображениям приватности нельзя использовать контактный e-mail.

Обязательные поля отмечены звездочкой (*).

(Отображаемое имя должно иметь длину от 3 символов до 31 символа.)

Нажимая Отправить, Вы принимаете Условия использования developerWorks.

 


Вся введенная информация защищена.


static.content.url=http://www.ibm.com/developerworks/js/artrating/
SITE_ID=40
Zone=Технология Java
ArticleID=499436
ArticleTitle=Новые математические возможности Java: Часть 1. Вещественные числа
publish-date=07072010