2014 dxdy logo

Научный форум dxdy

Математика, Физика, Computer Science, Machine Learning, LaTeX, Механика и Техника, Химия,
Биология и Медицина, Экономика и Финансовая Математика, Гуманитарные науки




Начать новую тему Ответить на тему
 
 Определитель
Сообщение23.02.2017, 15:07 


11/11/12
172
Здравствуйте! Написал программу, приводящую матрицу к диагональной элементарными преобразованиями. На некоторых матрицах она работает корректно, а на других -- нет. Например, если ввести такую:
Код:
30 20 15 12
20 15 12 15
15 12 15 20
12 15 20 30

программа выдаёт следующее:
Код:
30.000000  0.000000  0.000000  0.000000
0.000000  0.000000  0.000000  0.000000
0.000000  0.000000  -38280596832649256.000000  0.000000
0.000000  0.000000  0.000000  -10.349020

Что интересно, определитель она считает верно, хотя видно в данном случае, что на диагонали есть нуль (а детерминант равен 2639). Помогите, пожалуйста, найти ошибку.
Код:
   
#include<stdio.h>
#include<stdlib.h>
#include<math.h>

void det(double **m, int *k, int n, double *D);

void det(double **m, int *k, int n, double *D)
{
   
    int i, j, p, u, s;
    double x, y;
   
    (*D)=1;
   
    for(j=0; j<n; j++)
    {
        for(i=0; i<n; i++)
        { 
          if((m[i][j]!=0)&&(j<=i))
          {
              for(p=0; p<n; p++)
              {
                  if(p!=i)
                  {
                    y=m[p][j]/m[i][j];
                    for(u=0; u<n; u++)
                    {
                        m[p][u]-=y*m[i][u];
                    }
                  }
                 
              }
             
             for(s=0; s<n; s++)
             {
               if(j!=i)
               {
                 x=m[j][s];
                 m[j][s]=m[i][s];
                 m[i][s]=x;
                 (*D)=(*D)*(-1);
               }
             }
         
         
          }
         
         
         
        }
    }
   
 


 Профиль  
                  
 
 Re: Определитель
Сообщение23.02.2017, 15:21 
Заслуженный участник
Аватара пользователя


06/10/08
6422
А как Вы делаете вывод матрицы? Может, это не нуль на диагонали, а маленькое число?

 Профиль  
                  
 
 Re: Определитель
Сообщение23.02.2017, 15:28 


11/11/12
172
Xaositect в сообщении #1194794 писал(а):
А как Вы делаете вывод матрицы? Может, это не нуль на диагонали, а маленькое число?



Код:
for(i=0; i<n; i++)
  {
    for(j=0; j<k[i]; j++)
    {
      fprintf(f," %lf ",m[i][j]);
      if(j==k[i]-1) fprintf(f,"\n");
    }
  }


-- 23.02.2017, 15:29 --

Откуда может взяться это маленькое число?

 Профиль  
                  
 
 Re: Определитель
Сообщение23.02.2017, 16:13 
Аватара пользователя


07/01/15
1145
function в сообщении #1194797 писал(а):
Откуда может взяться это маленькое число?

Цифр после запятой может не хватить $-$ даблы не мелочатся.

 Профиль  
                  
 
 Re: Определитель
Сообщение23.02.2017, 17:41 


11/11/12
172
Как это можно исправить?

 Профиль  
                  
 
 Re: Определитель
Сообщение23.02.2017, 19:23 
Заслуженный участник


27/04/09
28128
Заменить lf на другой спецификатор формата — например, e (или E, отличается регистром).

Справка по printf и спецификаторам формата.

 Профиль  
                  
 
 Re: Определитель
Сообщение23.02.2017, 19:46 


11/11/12
172
arseniiv в сообщении #1194828 писал(а):
Заменить lf на другой спецификатор формата — например, e (или E, отличается регистром).


Теперь, конечно, они отображаются, но проблема не ушла: почему эти числа возникают, непонятно. За массив я не вылезаю, вроде, да и деления на ноль нет... :?

 Профиль  
                  
 
 Re: Определитель
Сообщение23.02.2017, 20:02 
Заслуженный участник


20/08/14
11065
Россия, Москва
function
Мне думается дело в конечной точности деления. Т.е. код
Код:
y=m[p][j]/m[i][j]; y*m[i][u];
не обязательно (не для всех возможных аргументов) восстановит точное исходное значение m[p][j] при u=j, будет какая-то погрешность. И эта погрешность вылезет на следующей итерации внешних циклов и вместо точного равенства будет лишь приближённое. А потом будет деление на эту маленькую разницу - и получите большое число, которого на самом деле быть не должно.

 Профиль  
                  
 
 Re: Определитель
Сообщение23.02.2017, 20:24 


10/04/12
704
Напрашивается такой вариант: (1) выполнить все действия на калькуляторе (2) вывести определители после каждой итерации (3) сравнить результаты.

 Профиль  
                  
 
 Re: Определитель
Сообщение23.02.2017, 20:35 
Заслуженный участник


20/08/14
11065
Россия, Москва
Калькуляторы тоже не всегда правильно округляют, так что это не панацея.
Если входные данные в виде целых или рациональных чисел, то выходом было бы проводить все операции тоже в рациональных числах (с числителем и знаменателем дроби), при этом результат будет точным. Для действительных чисел решение в общем виде без переделки алгоритма вряд ли существует (впрочем могу и не знать).

 Профиль  
                  
 
 Re: Определитель
Сообщение23.02.2017, 20:45 


11/11/12
172
Dmitriy40 в сообщении #1194837 писал(а):
function
Мне думается дело в конечной точности деления. Т.е. код
Код:
y=m[p][j]/m[i][j]; y*m[i][u];
не обязательно (не для всех возможных аргументов) восстановит точное исходное значение m[p][j] при u=j, будет какая-то погрешность. И эта погрешность вылезет на следующей итерации внешних циклов и вместо точного равенства будет лишь приближённое. А потом будет деление на эту маленькую разницу - и получите большое число, которого на самом деле быть не должно.

Спасибо Вам большое! Я просто обнулил соответствующие позиции ($u=j$), и всё получилось!

 Профиль  
                  
 
 Re: Определитель
Сообщение23.02.2017, 20:51 
Заслуженный участник


20/08/14
11065
Россия, Москва
function
Я хотел сразу предложить этот метод, но потом подумал что погрешность при вычислении остальных элементов (при $u\ne j$) так не уберётся. И при других исходных данных может снова вылезти. Так что не считайте это окончательным решением.

 Профиль  
                  
 
 Re: Определитель
Сообщение23.02.2017, 20:56 
Заслуженный участник


27/04/09
28128
Известна ещё модификация этого метода с перестановкой строк так, чтобы делилось всегда на как можно большее по модулю число (кажется), и в результате выдаются кроме матрицы (если что-то нужно считать дальше) ещё и новые индексы строк.

 Профиль  
                  
 
 Re: Определитель
Сообщение23.02.2017, 21:16 
Заслуженный участник


20/08/14
11065
Россия, Москва
Да, известна, но такая модификация лишь облегчает проблему потери точности, не исключая её полностью.

 Профиль  
                  
 
 Re: Определитель
Сообщение23.02.2017, 21:20 
Заслуженный участник


27/04/09
28128
Ну, это понятно. :-) Я думаю, ваш совет насчёт рациональных чисел, если в все элементы матрицы предполагаются небольшими целыми, самый лучший. Правда, наверняка они здесь небольшие и целые просто для простоты (хотя ещё лучше не угадывать за ТС).

 Профиль  
                  
Показать сообщения за:  Поле сортировки  
Начать новую тему Ответить на тему  [ Сообщений: 15 ] 

Модераторы: Karan, Toucan, PAV, maxal, Супермодераторы



Кто сейчас на конференции

Сейчас этот форум просматривают: нет зарегистрированных пользователей


Вы не можете начинать темы
Вы не можете отвечать на сообщения
Вы не можете редактировать свои сообщения
Вы не можете удалять свои сообщения
Вы не можете добавлять вложения

Найти:
Powered by phpBB © 2000, 2002, 2005, 2007 phpBB Group