Работа алгоритмов на устройствах с разными процессорами

#java #android #алгоритм

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


Meizu M2 Note (MT6753 8 * 1.3Ghz Cortex A53, 64bit): 229 циклов;
ASUS Zenfon 5 (Intel Atom Z2560 2 * 1.6Ghz x86): 145 циклов;
Lenovo A850+ (MT6592): 149 циклов;
Samsung Galaxy A3 (Qualcomm Snapdragon 410 4*1.2Ghz Cortex A53 64bit): 1878 циклов.


Еще заметил, что количество циклов для всех устройств с qualcomm snapdragon одинаковое:
1878.
Код алгоритма:

ArrayList box = new ArrayList();
MyAdapter adapter;
ListView listView;
double precision = 1E-28;
double K1F1, K2F1, K3F1, K4F1, K1F2, K2F2, K3F2, K4F2, K1F3, K2F3, K3F3;
int flag=0;

 //инициализацию вытащил из обработчика кнопки, он сам не важен. А эти переменные
используются в countFx() и countWx()
 K1F1=16d;
           K2F1=16d;
           K3F1=1d;
           K4F1=16d;
           K1F2=1d;
           K2F2=1d;
           K3F2=1d;
           K4F2=3d;
           K1F3=1d;
           K2F3=1d;
           K3F3=0d;


 private double[] Countfx(double[] x){//вектор Fx
    double[] fx = new double[3];
    fx[0] = K1F1*Math.pow(x[0],4)+K2F1*Math.pow(x[1],4)+K3F1*Math.pow(x[2],4)-K4F1;
    fx[1] = K1F2*Math.pow(x[0],2)+K2F2*Math.pow(x[1],2)+K3F2*Math.pow(x[2],2)-K4F2;
    fx[2] = K1F3*Math.pow(x[0],3)-K2F3*x[1]-K3F3;

    return fx;
}
private double Countdet(double[][] ar){//определитель матрицы 3*3
    double detar;
    detar = ar[0][0]*ar[1][1]*ar[2][2] + ar[0][1]*ar[1][2]*ar[2][0] +
            ar[0][2]*ar[1][0]*ar[2][1] -
            ar[2][0]*ar[1][1]*ar[0][2] - ar[2][1]*ar[1][2]*ar[0][0] - ar[2][2]*ar[1][0]*ar[0][1];
    return detar;
}
private double[][] CountWx(double[] x){//матрица якоби
    double[][] Wx = new double[3][3];
    Wx[0][0] = K1F1*4d*Math.pow(x[0],3);
    Wx[1][0] = K2F1*4d*Math.pow(x[1],3);
    Wx[2][0] = K3F1*4d*Math.pow(x[2],3);
    Wx[0][1] = K1F2*2d*x[0];
    Wx[1][1] = K2F2*2d*x[1];
    Wx[2][1] = K3F2*2d*x[2];
    Wx[0][2] = K1F3*3d*Math.pow(x[0],2);
    Wx[1][2] = K2F3*x[1];
    Wx[2][2] = 0;

    return Wx;
}
private double[][] CountMinorMatrix(double[][] matrix){//матрица миноров (для нахождения
обратной матрицы)
    double[][] minor = new double[3][3];
    minor[0][0] = matrix[1][1]*matrix[2][2]-(matrix[1][2]*matrix[2][1]);
    minor[1][0] = matrix[0][1]*matrix[2][2]-(matrix[0][2]*matrix[2][1]);
    minor[2][0] = matrix[0][1]*matrix[1][2]-(matrix[1][1]*matrix[0][2]);
    minor[0][1] = matrix[1][0]*matrix[2][2]-(matrix[1][2]*matrix[2][0]);
    minor[1][1] = matrix[0][0]*matrix[2][2]-(matrix[0][2]*matrix[2][0]);
    minor[2][1] = matrix[0][0]*matrix[1][2]-(matrix[1][0]*matrix[0][2]);
    minor[0][2] = matrix[1][0]*matrix[2][1]-(matrix[1][1]*matrix[2][0]);
    minor[1][2] = matrix[0][0]*matrix[2][1]-(matrix[0][1]*matrix[2][0]);
    minor[2][2] = matrix[0][0]*matrix[1][1]-(matrix[0][1]*matrix[1][0]);

    return minor;
}
private double[][] getTransposedMatrix(double[][] matrix){//транспонирование матрицы
    double[][] tr = new double[3][3];
    for(int i=0; i<3; i++){
        for(int j=0; j<3; j++){
            tr[i][j]=matrix[j][i];
        }
    }
    return tr;
}
private double[][] MultiplyMatrixAtNumber(double[][] m, double num){//умножение матрицы
3*3 на число
    for(int i=0; i<3; i++){
        for(int j=0; j<3; j++){
            m[i][j]*=num;
        }
    }
    return m;
}
private double[][] getInverseMatrix(double[][] m){//нахождение обратной матрицы
    double det = Countdet(m);

    m=CountMinorMatrix(m);

    if(m[1][0]>0) m[1][0]=-Math.abs(m[1][0]);
    else m[1][0]=Math.abs(m[1][0]);
    if(m[0][1]>0) m[0][1]=-Math.abs(m[0][1]);
    else m[0][1]=Math.abs(m[0][1]);
    if(m[2][1]>0) m[2][1]=-Math.abs(m[2][1]);
    else m[2][1]=Math.abs(m[2][1]);
    if(m[1][2]>0) m[1][2]=-Math.abs(m[1][2]);
    else m[1][2]=Math.abs(m[1][2]);

    m=getTransposedMatrix(m);
    det = 1d/det;
    m=MultiplyMatrixAtNumber(m,det);

    return m;
}
private double[] MultiplyMatrix(double[][] m1, double[] m2){//умножение матрицы 3*3
и вектора 1*3
    double[] res = new double[3];
    for(int i=0; i<3; i++)
        res[i]=0;
    for(int j=0; j<3; j++){
        for(int i=0; i<3; i++){
            res[j]+=m1[i][j]*m2[i];
        }
    }
    return res;
}
private double[] CountDifMatrix(double[] m1, double[] m2){//разность матриц
    double[] sum = new double[3];
    for(int i=0; i<3; i++){
        sum[i] = m1[i]-m2[i];
    }
    return sum;
}
private boolean compare(double[] m1, double[] m2){//сравнение векторов
    int n=0;
    for(int i=0; i<3; i++){
        if(m1[i]==m2[i])
            n++;
    }
    if(n>=2)
        return true;
    else
        return false;
}
private void CountFunction(){
    double[] X0={1d,1d,1d};
    double[] Fx;
    double[][] Wx;
    double det;
    double[][] WxInverse;
    double[] multRes;
    double[] result;
    int w=0, it=0, col=0;
    while(w==0 && col<1){

        Fx=Countfx(X0);//нахождение вектора Fx
        Wx=CountWx(X0);//нахождение матрицы якоби
        det=Countdet(Wx);//нахождение определителя матрицы для проверки на 0
        if(det==0){
            break;
        }
        WxInverse=getInverseMatrix(Wx);//нахождение обратной матрицы
        multRes=MultiplyMatrix(WxInverse,Fx);
        result=CountDifMatrix(X0,multRes);
        Box b=new Box(result);//объект просто хранит массив (для ArrayList нужен
абстрактный тип)
        outputNumbersinLog(result);//лог
        box.add(b);// запись результата цикла в arraylist
        boolean compare=compare(X0,result);//сравнение результата текущего и предыдущего
циклов
        if(compare)
            col++;
        for(int i=0; i<3; i++){
            if(Math.abs(result[i])
Ответ 1

Для того, чтобы гарантировать правильность вычислений, необходимы две вещи.

Во-первых, необходимо отмаркировать классы, занимающиеся вычислениями, модификатором
strictfp. (Это отметил @zRrr в комментарии.) Без этого вычисления с плавающей точкой
выполняются оптимальным по скорости для текущего процессора способом, с возможным ущербом
в точности.

Во-вторых, необходимо вместо класса Math использовать StrictMath, который осуществляет
вычисления в режиме FP-strict.