МК61 - вычисление интеграла - вопрос/баг с точностью метода ТРАПЕЦИЙ

добрый день

расчет интеграла 2мя методами - "м.прямоугольников" и "м.трапеций"

оба метода/программы считают (вроде бы как :о), но "трапеции" - более сложный (для програмирования), и менее точных, хотя должно быть наоборот (если не изменяет память и учебники :о)

с математикой "прямоугольников", вроде как "понятно и так" (надеюсь), с "трапециями" чуть сложнее

// общая формула
In = Integral(sin(x)dx) dx={0...Pi}

// формула "прямоугольников
In = dx*[f(x[0]) + f(x[1]) + ... f(x[n-1])]

// формула трапеций (преобразовал)
I = dx/2 * [f(x[0])   + 2f(x[1] + 2f(x[2]) + ... + 2f(x[n-1]) + f(x[n])   ]
-->
I = dx   * [f(x[0])/2 +  f(x[1] +  f(x[2]) + ... +  f(x[n-1]) + f(x[n])/2 ]


таблица: время выполнения расчета
-------------------------------------------
dx   Integral  time/мин.сек
-------------------------------------------
0.1  1.9995479 02.59 - прямоугольники
0.1  1.9966292 05.09 - трапеции

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

дотошно проверил все, чЁ можно, (руки помыл, сделал зарядку, постоял на голове) - найти "жучка" не смог

далее более полный выкус из "дневников" и сканы "трапеций" из тетради, возможно поможет

спасибо

///////////////////////////////////////////////////////////////////////////////
### Рсчет интеграла синуса методом прямоугольников ###
///////////////////////////////////////////////////////////////////////////////

формулы брал тут - "астанин - прим ПМК инж. науч.djvu", стр. 72

In = Integral(sin(x)dx) dx={0...Pi}

In = dx*[f(x[0]) + f(x[1]) + ... f(x[n-1])]

///////////////////////////////////////////////////////////////////////////////
ver.1 (без подпрограммы)
///////////////////////////////////////////////////////////////////////////////

01 F sin - собственно вычисление y=f(x) или более точно в итерации: y[n]=f(x[n])
11 НОП   - удалил ненужный код (можно сдвинуть код на один шаг :о)

+--+------+--+------+--+------+--+
|  | 00   |xx| 10   |xx| 20   |xx|  
+--+------+--+------+--+------+--+
|00| ИП3   63  -     11  00   00
|01| F sin 1C  НОП   
|02| ИП4   64  F>=0  59
|03| +     10  00    00
|04| П4    44  ИП4   64
|05| ИП3   63  ИП2   62
|06| ИП2   62  *     12 
|07| +     10  П4    44
|08| П3    43  С/П   50
|09| ИП1   61  БП    51
+--+------+--+------+--+------+--+

// распределение регистров
a   -> R0 (мин.  граница - dx[min])
b   -> R1 (макс. граница - dx[max])
dx  -> R2 (шаг приращения интеграции/точность вычисления)
dxs -> R3 (изменяемы x[i] - в процессе итерации/"суммируемый шаг приращения")
In  -> R4 (интеграл)

// ввод исх. данных
0   П0  
FPi П1 // 3,1415...
0.1 П2 // шаг приращения dx (~ точность вычисления)
0   П3
0   П4

// запуск:
В/О С/П --> I(интеграл)

// повторный запуск (с другими параметрами):
изменение шага, В/О С/П 

таблица: время выполнения расчета, в зов-ти от заданного шага интеграции "dx"
---------------------------
dx   Integral  time/мин.сек
---------------------------
0.1  1.9995479 02.59
0.01 1.9999907 26.40 
---------------------------

///////////////////////////////////////////////////////////////////////////////
ver.2 (с исп. подпрограммы)
///////////////////////////////////////////////////////////////////////////////

A -
B - L
C - C
D - Г
E - E
F - ''

20 - адресс ППРГ

+--+------+--+------+--+------+--+
|  | 00   |xx| 10   |xx| 20   |xx|  
+--+------+--+------+--+------+--+
|00| ПП    53  -     11  ИП3   63
|01| 20    20  F>=0  59  F sin 1C 
|02| ИП4   64  00    00  В/О   52
|03| +     10  ИП4   64
|04| П4    44  ИП2   62
|05| ИП3   63  *     12
|06| ИП2   62  П4    44 
|07| +     10  С/П   50
|08| П3    43  БП    51
|09| ИП1   61  00    00
+--+------+--+------+--+------+--+

время выполнени то-же (собственно, там ни чего не изменилось, кроме вызова ППРГ)

///////////////////////////////////////////////////////////////////////////////
### Рсчет интеграла синуса методом трапеций ###
///////////////////////////////////////////////////////////////////////////////

формулы брал тут - "астанин - прим ПМК инж. науч.djvu", стр. 72

выкус:
при линейной интерполяции функ. f(x) получается т.н. формула трапеций , кот. 
апроксисимирует интеграл суммой площадей трапеций
- с основанием (x[i],x[i+1])
- и высотой, изменяющейся линейно от f(x[i]) -> f(x[i+1])

I = dx/2 * [f(x[0]) + 2f(x[1] + 2f(x[2]) + ... + 2f(x[n-1]) + f(x[n]) ]

точность интегрирования можно повысить, если заменить линейную интерполяцию квадратичной (формула симпсона)

I = dx/3 * [f(x[0]) + 4f(x[1]) + 2f(x[2]) + 4f(x[3]) + ... + 4f(x[n-1]) + f(x[n]) ]

в книге реализуется 33 и 19 шагов, правда, лажЁва - пограничные "нестандартные" значения функции вводятся мануально!!! 

//
// преобразовал немного метод трапеций
// 

I = dx/2 * [f(x[0])   + 2f(x[1] + 2f(x[2]) + ... + 2f(x[n-1]) + f(x[n])   ]
-->
I = dx   * [f(x[0])/2 +  f(x[1] +  f(x[2]) + ... +  f(x[n-1]) + f(x[n])/2 ]

A - '-'
B - 'L'
C - 'C'
D - 'Г'
E - 'E'
F - ' '

+--+------+--+------+--+------+--+------+--+------+--+------+--+------+--+
                                           |   ППРГ  | блок А  | блок В  |
+--+------+--+------+--+------+--+------+--+------+--+------+--+------+--+
|  | 00   |xx| 10   |xx| 20   |xx| 30   |xx| 40   |xx| 50   |xx| 60   |xx|
+--+------+--+------+--+------+--+------+--+------+--+------+--+------+--+
|00| Сх    0Г  60  бл.В  -     11  ИП2   62  ИП3   63  ИП4   64  ИП6   66 
|01| П3    43  ИП5   65  Fx>=0 59  *     12  F sin 1С  ИП7   67  ИП7   67
|02| П4    44  1     01  50  бл.А  П4    44  П7    47  +     10  +     10
|03| П5    45  +     10  ИП6   66  С/П   50  В/О   52  П4    44  П6    46
|04| П6    46  П5    45  ИП7   67  БП    51            БП    51  БП    51
|05| П7    47  ИП3   63  +     10  00    00            06   т.1  11   т.2
|06| ПП    53  ИП2   62  2     02  
|07| 40  ППРГ  +     10  /     13
|08| ИП5   65  П3    43  ИП4   64 
|09| Fx!=0 57  ИП1   61  +     10
+--+------+--+------+--+------+--+------+--+------+--+------+--+------+--+

// общ. комменты

00-05 - подг/очистака рабочих регистров
06    - т.1 начало основной программы
11    - т.2 
40    - ППРГ y[i]=f(x[i]) / подпрограмма вычисления значения функции в данной точке интегрирования
50    - блок А
60    - блок Б

// переходы

10 - блок.В
22 - блок.А
55 - т.1
65 - т.2

// исх.данные

a    - П0 // нач. значение х
b    - П1 // кон. значение х
dx   - П2 // приращение / шаг интегрирования
x[i] - П3 // тек. значение х
I    - П4 // интеграл
i    - П5 // счетчик итераций
TMP  - П6 // врем. переменная, для сохранения "нестандарт." членов полинома (крайние значения)
y[i] - П7 // тек. значение у

// подготовка данных

0    П0
F Pi П1
0.1  П2 изменение шага, С/П 


// запуск
В/О С/П --> интеграл

// повторный запуск (с другими параметрами):
изменение шага, С/П или В/О С/П 

таблица: время выполнения расчета, в зов-ти от заданного шага интеграции "dx"
---------------------------
dx   Integral  time/мин.сек
---------------------------
0.1  1.9966292 05.09
0.01 1.9999487 49.хх 
---------------------------
Прикрепленный файлРазмер
Plain text icon integral.txt7.47 KB
Image icon resized_scan-180112-0017.jpg194.83 KB
Image icon resized_scan-180112-0018.jpg199.16 KB

p.s.
уточняю:
- имеется расчет 2мя методами (1-м методом - две программы, 2-м методом - одна большая)
- второй метод должен быть более точный, на практике - все наоборот
- как итог вопрос, может кто прояснит?

если я ошибаюсь, то поправьте а не критикуйте :о)

Вычислять интегралы на МК-61 — дело хорошее. Каждый из участников сайта этим занимался, в своё время. Никаких особых аномалий на этом пути не обнаружено, математика работает как надо.

Конечно, здесь никто не будет специально искать ошибки в твоих программах. Это считается ответственностью программиста. Но если упёрся в тупик, можно дать пару советов.

1. Вычисление интегралов это типовая задача. Есть справочники Дьяконова и Цветкова, например, где такая работа уже произведена. Можно взять их проверенные программы, как эталон, и посмотреть, какая из твоих программ ошибается.

2. Попробовать интегрировать другую функцию. В математике не всегда бывает, что из двух приближений f1 и f2 к точному значению f ближе окажется «более точное». Главное, чтобы оба значения f1 и f2 были в пределах своей гарантированной погрешности от f. При этом любое из них может быть равно f или наоборот, максимально отличаться от f, в пределах своей погрешности.

### астанин
для начала я взял "астанин - прим ПМК инж. науч ..." стр. 72

там предлагается "пограничные" значение вводить вручную, предварительно посчитав на калькуляторе (хороша-же программа, которая без счет и палочек не работает) - закрыли, пошли смотреть другие справочники

(* пограничные) - те, которые не входят в стандартные условия суммирования/интегрирования

### цветков
есть след. методы - м.прямоугольников, симсона с заданной точностью, м.суммирования
м. трапеций нет, закрываем, идем дальше

### дьяконов 5.9 стр.243
- м.прямоугольников, трапеций, симпсона (парабол), симсона с задан. точностью, м.ньютона котеса, уэддля, чебышева, гауса, интегрирование табличных данных при постоянном шаге м. трапеций, ... просто кладезь материала

но! описание весьма и весьма мудрено

прг.5.63 / стр.244

дословно:
"Ввод: данные f(x), m, a и b в регистр X"

тут у меня голова накалилась! это чЁ, все данные в один регистр Х?! какой такой шашлык-машлык?!
ладно, может дальше будет уточнение, а далее и в самом деле упоминается про регистры:

дословно:
"текущее значение x->RA, регистры 0,А,В,С заняты"

фысе! ни какого упоминания про распределение (что, где, куда...)

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

поэтому я написал свою программу :о)
были подозрения на баг либо в блок-схеме, либо в реализации...

подумал, возможно, кто-ни будь заинтересуется, глянет... подскажет :о) да-да, нет-нет, без всяких притензий

спасибо

п.с.
к стати, у дьяконова - в таблице так же указано, что м.трапеций "проигрывает", хотя во многих справочниках говориться, что "трапеции" более точный метод...

про все остальные "вариации" и "биения" - в курсе, согласен, допускаю, политику партии понимаю и не осуждаю...

если я ошибаюсь, то поправьте а не критикуйте :о)

Научиться пользоваться справочником Дьяконова полезно. Несколько значений в RX означает, что между вводом этих значений нажимается С/П. Это должно быть ясно из текста программ. Также эти моменты объясняются в начале справочника.

Потратьте своё время, чтобы разобраться в программе профессора. Они проверены и составлены так хорошо, что идут даже на МК-161. Если вы считаете, что разбираться в чужих программах сложно – зачем просить других копаться в ваших?

P.S. Методы трапеций и прямоугольников одного порядка, точность примерно одинаковая. По моему опыту интегрирования на ПМК хорошую точность даёт метод Симпсона (парабол). Также я писал и выкладывал здесь программы для более точного метода численного интегрирования, но для МК-152/161.

Помню в молодости считал интеграл вероятностей в диапазоне от -10 до 0 методом трапеций и методом Симпсона. Симпсон оказался точнее при в 20 раз более грубом разбиении (единицу делил на 50 интервалов), чем трапеции (единицу делил на 1000 интервалов). С тех пор использую либо Симпсона, либо Гаусса, а прямоугольники и трапеции не использую.

Мои программируемые калькуляторы:
Б3-21, Б3-34, МК-61, МК-52, МК-85
CASIO: cfx-9850GB+, fx-9750G+, fx-9750GII, fx-9860G, Algebra fx-2.0, fx-5800P, fx-7400G+
HP: 50G, 48G, 35s
TI: Nspire-CAS, Voyage-200, 89Titanium
SHARP EL-9600G

спасибо

если я ошибаюсь, то поправьте а не критикуйте :о)