Интеграл, непонятен один момент Опции

[Zelepuk]

Хочу найти функцию вычисления определённого интеграла для произвольного массива чисел.
В сети полно примеров взятия интеграла, например:

CODE

double integral(a, b, f, eps)
double a, b; /* концы отрезка */
double eps; /* требуемая точность */
double (*f)(); /* подынтегральная функция */
{
register long i;
double fab = (*f)(a) + (*f)(; /* сумма на краях */
double h, h2; /* шаг и удвоенный шаг */
long n, n2; /* число точек разбиения и оно же удвоенное */
double Sodd, Seven; /* сумма значений f в нечетных и в
четных точках */
double S, Sprev;/* значение интеграла на данной
и на предыдущей итерациях */
double x; /* текущая абсцисса */

niter = 0;
n = 10L; /* четное число */
n2 = n * 2;

h = fabs(b - a) / n2; h2 = h * 2.0;

/* Вычисляем первое приближение */
/* Сумма по нечетным точкам: */
for( Sodd = 0.0, x = a+h, i = 0;
i < n;
i++, x += h2 )
Sodd += (*f)(x);

/* Сумма по четным точкам: */
for( Seven = 0.0, x = a+h2, i = 0;
i < n-1;
i++, x += h2 )
Seven += f(x);

/* Предварительное значение интеграла: */
S = h / 3.0 * (fab + 4.0 * Sodd + 2.0 * Seven );
do{
niter++;
Sprev = S;

/* Вычисляем интеграл с половинным шагом */
h2 = h; h /= 2.0;
if( h == 0.0 ) break; /* потеря значимости */
n = n2; n2 *= 2;

Seven = Seven + Sodd;
/* Вычисляем сумму по новым точкам: */
for( Sodd = 0.0, x = a+h, i = 0;
i < n;
i++, x += h2 )
Sodd += (*f)(x);

/* Значение интеграла */
S = h / 3.0 * (fab + 4.0 * Sodd + 2.0 * Seven );

} while( niter < 31 && fabs(S - Sprev) / 15.0 >= eps );
/* Используем условие Рунге для окончания итераций */

return ( 16.0 * S - Sprev ) / 15.0 ;
/* Возвращаем уточненное по Ричардсону значение */
}

Непонятно следущее: f - указывает на подынтегральную функцию. Математически это понятно. Но в Си же это подпрограмма! Я наверное принципиально ен понимаю чего-то но как можно интеграл от подпрограммы взять? или программа должна возвращать указатель на массив?
я же надеялся найти реализацию где в качестве указателя будет указатель на массив, поскольку функция представлена в виде набора отсчётов (массива) от a до b.

Заранее спасибо за разъяснения.

Сообщение отредактировал Zelepuk - Dec 14 2012, 09:38

[MrYuran]

double f() - это функция, возвращающая значение в конкретной точке.
Все в соответствии с канонами.

[Zelepuk]

double f() - это функция, возвращающая значение в конкретной точке.
Все в соответствии с канонами.

непонятно как производится "перебор" этих точек. Нужно же посчитать для группы точек. Или в функции f() после каждого вызова должен быть инкремент указателя?
Думаю намного проще было бы сделать что-то вроде цикла где перебирались бы элементы массива.
Я полагаю конструкции вида (*f)(a) можно просто заменить на указатель на массив(с тем же индексом)?

Сообщение отредактировал Zelepuk - Dec 14 2012, 10:02

[xemul]

непонятно как производится "перебор" этих точек. Нужно же посчитать для группы точек. Или в функции f() после каждого вызова должен быть инкремент указателя?
Думаю намного проще было бы сделать что-то вроде цикла где перебирались бы элементы массива.

Вы в код, Вами же выложенный заглядывали? Там есть кучка обращений (*f)(x), и именно в циклах.

Я полагаю конструкции вида (*f)(a) можно просто заменить на указатель на массив(с тем же индексом)?

Ваша задача - нарисовать f(x), т.е. собственно подынтегральную функцию. Будет ли функция представлена аналитически или какой-либо интерполяцией - без разницы, но для любого x на [a, b] она должна вернуть осмысленное значение.

[DRUID3]

Я полагаю конструкции вида (*f)(a) можно просто заменить на указатель на массив...

Ну а если она не массивом задана?

(с тем же индексом)?

??

[Палыч]

непонятно как производится "перебор" этих точек. Нужно же посчитать для группы точек.

Для этого в программе присутствуют циклы. Аргумент вызываемой функции в каждом цикле увеличивается на величину шага

Код

x += h2

затем вычисляется значение функции в точке х и вычисляется сумма значений

Код

S... += (*f)(х)

[demiurg_spb]

Хочу найти функцию вычисления определённого интеграла для произвольного массива чисел.

Может я чего-то не понимаю?
В чём смысл (достоинство) приведенного вами кода для дискретной функции?
Для массива чисел очень просто найти площадь под графиком (интеграл) просто просуммировав площади между соседними точками (можно методом прямоугольников - очень быстро, а можно трапецией - чуть медленнее но чуть точнее).

[DRUID3]

Может я чего-то не понимаю?
В чём смысл (достоинство) приведенного вами кода для дискретной функции?

Точность... Ошибка вычисления... Функция на входе то ожидается аналитическая.

[demiurg_spb]

Точность... Ошибка вычисления... Функция на входе то ожидается аналитическая.

Это-то понятно, но у него не аналитическая функция, а массив. И мне думается что этот алгоритм вообще не применим для его случая.

[Xenia]

Для массива чисел очень просто найти площадь под графиком (интеграл) просто просуммировав площади между соседними точками (можно методом прямоугольников - очень быстро, а можно трапецией - чуть медленнее но чуть точнее).

Интересно, откуда всякий раз выползает предубеждение, что метод прямого суммирования дискрет менее точен, чем сумма трапеций, а последний более трудоёмок?

Ведь если точки идут по сетке через равные интервалы, то площадь отдельной трапеции равна полусумме соседних отсчетов (F_i-1 + F_i)/2 или (что тоже самое!) сумме половинок F_i-1/2 + F_i/2. Но когда начнём все эти трапеции друг с дружкой складывать, то соседние половинки объединятся в целое! В результате чего получится простая сумма дискрет.

Впрочем, небольшая разница тут всё же есть - на самых краях, там, где половинок не хватает. Поэтому метод трапеций отличается от метода прямого суммирования только величиной (F₀ + F_n)/2, где F₀ - начальная точка, а F_n - конечная. Разница в трудоёмкости, как видим, здесь совершенно ничтожна. А те, у кого она велика, занимаются избыточным суммированием и дележкой на 2 на каждом обороте цикла.

[demiurg_spb]

Так я прямо и написал об этом что чуть лучше и чуть трудозатратнее. Никакого заблуждения как и предубеждения нет:-)

Разница в трудоёмкости, как видим, здесь совершенно ничтожна.

Всё от задачи зависит. Вот если массив даблов да ещё и на 8-ми битнике.... я бы разницу точно ощутил

[_Ivana]

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

[Xenia]

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

Тогда demiurg_spb на своём 8-битнике просто зашьётся.

[_Ivana]

А кто его знает, с какой частотой у автора идут отсчеты его функции? Он об этом скромно умолчал. Может она близка к частоте Найквиста (при условии ограниченности спектра функции) и остается только использовать суровые методы аппроксимации между отсчетами Кстати, если "произвольные отсчеты" идут ещё реже, то задача становится ещё интереснее. Или, если у него отсчеты какой-нибудь функции Дирихле, к примеру

[Xenia]

А кстати, а чего это народ функцию от указателя так вычурно вычисляет?
S += (*f)(х);
На мой взгляд (а мои взгляды могут расходиться со взглядами компилятора ), такую функцию следует вызывать как обычную
S += f(х);
Тем более что прецедент уже есть - существует два разных способа объявления массива:
int array[100], *parray=array;
И тем не менее, обоими массивами пользуются одинаково, когда достают из них элементы:
S += array[x];
или
S += parray[x];
а вовсе не (*parray)[x];
И все только потому, что array - тоже указатель, только константный. А, значит, и оператор [] будет их разыменовывать правильно в обоих случаях.

По аналогии с этим, и функцию по указателю можно вызывать напрямую, как f(x), поскольку обычная функция (не указатель) тоже является указателем на функцию, только константным. Значит, и оператор () тоже должен работать в обоих случаях одинаково.
Или я не права?

[demiurg_spb]

Или я не права?

не совсем. попробуйте двумерный массив.
Ну или посмотрите что вам вернёт sizeof(array1) и что sizeof(array2) для таких вот случаев

Код

uint8_t array1[10];
uint8_t array2[2][5];

[Xenia]

не совсем. попробуйте двумерный массив.
Ну или посмотрите что вам вернёт sizeof(array1) и что sizeof(array2) для таких вот случаев

Код

uint8_t array1[10];
uint8_t array2[2][5];

Ну так я же не просила вас вызывать двухмерную функцию.
Речь шла о том, что оператор [] приложим к указателю, так же как и оператор (). Оба они требуют, чтобы слева от них стоял указатель, хоть бы и вычисляемый. Т.е. я не за идентичность указателей ратовала, а за единообразие выполнения операторов "скобочки". А именно, в выражении f() скобочки должны интерпретировать f, как указатель (в данном случае на функцию), а не как-то иначе. Зачем же тогда f предварительно разыновывать? Ведь если имя массива (без скобочек) это указатель, то и имя фукции без скобочек тоже должно быть указателем. А чем же ещё-то?

[ReAl]

А именно, в выражении f() скобочки должны интерпретировать f, как указатель (в данном случае на функцию), а не как-то иначе.

И это так и есть по стандарту.

6.3.2 Other operands
6.3.2.1 Lvalues, arrays, and function designators

...
3 Except when it is the operand of the sizeof operator or the unary & operator, or is a string literal used to initialize an array, an expression that has type ‘‘array of type’’ is converted to an expression with type ‘‘pointer to type’’ that points to the initial element of the array object and is not an lvalue.

4 A function designator is an expression that has function type. Except when it is the operand of the sizeof operator or the unary & operator, a function designator with type ‘‘function returning type’’ is converted to an expression that has type ‘‘pointer to function returning type’’.

Зачем же тогда f предварительно разыновывать?

Абсолютно незачем. Так как после этого разименованное ‘‘function returning type’’ моментально опять превращается в ‘‘pointer to function returning type’’ для того, чтобы скобочки () (которые в данном случае являются «function-call operator») могли этим воспользоваться.

6.5.2.2 Function calls
Constraints
1 The expression that denotes the called function shall have type pointer to function returning void or returning an object type other than an array type.

Как следствие, разыменовывать можно сколько угодно раз, это ничего не меняет.

Код

#include <stdio.h>

int (*pf)(const char *) = puts;

int main()
{
        pf("pf");
        // Разыменовываем, автоматически приводится к указателю, мы его опять разыменовываем и т.д.
        (****pf)("****pf");
        return 0;
}

С выборкой указателя на функцию из массива аналогично pfunc_array[index](arguments); и всё.
Но при большом желании можно и (*****pfunc_array[index])(arguments);. Эффект тот же.

И где-то на этом форуме об этом уже говорили. И вроде как не раз.

_____________________
Сегодня смотрел «День сурка». Опять…