Суммирование очень больших чисел без использования наборов инструментов

Я пытаюсь суммировать очень большие числа в MATLAB, такие как e^800 и e^1000, и получить ответ.

Я знаю, что в Double-Precision максимальное число, которое я могу представить, равно 1.8 * 10^308, в противном случае я получаю Inf, которое я получаю при попытке суммировать эти числа.

Мой вопрос в том, как мне оценить ответ для сумм очень, очень больших чисел, подобных этим, без использования vpa или какого-либо другого набора инструментов?

Должен ли я использовать строки? Можно ли это сделать с помощью логов? Могу ли я представить поплавки как m x 2^E, и если да, то как мне взять число, такое как e^700, и преобразовать его в это? Если число больше порогового значения для Inf, следует ли разделить его на два и сохранить в двух разных переменных?

Например, как я могу получить приблизительный ответ для:

e^700 + e^800 + e^900 + e^1000 ?


matlab approximation floating-point addition largenumber
person Sunden    schedule 09.09.2019    source источник
comment
Может быть, вы хотите использовать символьную математику? Невозможно легко закодировать такие большие числа, и еще сложнее закодировать сумму чисел с такой большой разницей в величине. Вы уверены, что вам действительно нужно это сделать?   -  person Cris Luengo    schedule 09.09.2019
comment
Почему вы не хотите использовать vpa ? Слишком медленно ? Это единственная доступная встроенная арифметическая машина переменной точности.   -  person obchardon    schedule 09.09.2019
comment
exp(1000) + exp(900) + exp(800) + exp(700) = exp(1000) * ( 1 + exp(-100) + exp(-200) + exp(-300)) поэтому с относительным ошибка порядка ехр(-100) сумма ехр(1000)   -  person dmuir    schedule 09.09.2019
comment
Это больше вопрос, смогу ли я сделать это без использования каких-либо дополнительных наборов инструментов. Я знаю, что могу использовать vpa или символический или эквивалентный.   -  person Sunden    schedule 11.09.2019

Ответы (2)


arrow_upward
2
arrow_downward

Возможным приближением является использование округленных значений этих чисел (лично я использовал Wolfram|Alpha), затем выполните «длинное сложение», как учат в начальной школе:

function sumStr = q57847408()
% store rounded values as string:
e700r = "10142320547350045094553295952312676152046795722430733487805362812493517025075236830454816031618297136953899163768858065865979600395888785678282243008887402599998988678389656623693619501668117889366505232839133350791146179734135738674857067797623379884901489612849999201100199130430066930357357609994944589";
e800r = "272637457211256656736477954636726975796659226578982795071066647118106329569950664167039352195586786006860427256761029240367497446044798868927677691427770056726553709171916768600252121000026950958713667265709829230666049302755903290190813628112360876270335261689183230096592218807453604259932239625718007773351636778976141601237086887204646030033802";
e900r = "7328814222307421705188664731793809962200803372470257400807463551580529988383143818044446210332341895120636693403927733397752413275206079839254190792861282973356634441244426690921723184222561912289431824879574706220963893719030715472100992004193705579194389741613195142957118770070062108395593116134031340597082860041712861324644992840377291211724061562384383156190256314590053986874606962229";
e1000r = "197007111401704699388887935224332312531693798532384578995280299138506385078244119347497807656302688993096381798752022693598298173054461289923262783660152825232320535169584566756192271567602788071422466826314006855168508653497941660316045367817938092905299728580132869945856470286534375900456564355589156220422320260518826112288638358372248724725214506150418881937494100871264232248436315760560377439930623959705844189509050047074217568";

% pad to the same length with zeros on the left:
padded = pad([e700r; e800r; e900r; e1000r], 'left', '0');

% convert the padded value to an array of digits:
dig = uint8(char(padded) - '0');

% some helpful computations for later:
colSum = [0 uint8(sum(dig, 1))]; % extra 0 is to prevent overflow of MSB
remainder = mod(colSum, 10);
carry = idivide(colSum, 10, 'floor');

while any(carry) % can also be a 'for' loop with nDigit iterations (at most)
  result = remainder + circshift(carry, -1);
  remainder = mod(result, 10);
  carry = idivide(result, 10, 'floor');
end

% remove leading zero (at most one):
if ~result(1)
  result = result(2:end);
end

% convert result back to string:
sumStr = string(char(result + '0'));

Это дает (округленный) результат:

197007111401704699388887935224332312531693805861198801302702004327171116872054081548301452764017301057216669857236647803717912876737392925607579016038517631441936559738211677036898431095605804172455718237264052427496060405708350697523284591075347592055157466708515626775854212347372496361426842057599220506613838622595904885345364347680768544809390466197511254544019946918140384750254735105245290662192955421993462796807599177706158188
person Dev-iL    schedule 12.09.2019

arrow_upward
0
arrow_downward

Опечатки исправлены ранее.

Десятичное приближение:

function [m, new_exponent] = base10_mantissa_exponent(base, exponent)
    exact_exp = exponent*log10(abs(base));
    new_exponent = floor(exact_exp);
    m = power(10, exact_exp - new_exponent);
end

Таким образом, значение e600 станет 3,7731 * 10260.

А значение 117150 станет 1,6899 * 10310.

Чтобы сложить эти два значения вместе, я взял разницу между двумя показателями и разделил на нее мантиссу меньшего члена. Тогда это так же просто, как сложить вместе мантиссы.

mantissaA = 3.7731;
exponentA = 260;
mantissaB = 1.6899;
exponentB = 310;

diff = abs(exponentA - exponentB);
if exponentA < exponentB
    mantissaA = mantissaA / (10^diff);
    finalExponent = exponentB;
elseif exponentB < exponentA
    mantissaB = mantissaB / (10^diff);
    finalExponent = exponentA;
end

finalMantissa = mantissaA + mantissaB;

Это было важно для меня, так как я выполнял такие суммы, как:

( ex) / ( xex)

От х=1 до х=1000.

person Sunden    schedule 13.09.2019