Подтвердить что ты не робот

Производные в C/С++?

У меня есть некоторые выражения, такие как x^2+y^2, которые я бы хотел использовать для математических вычислений. Одна из вещей, которые я хотел бы сделать, - принять частные производные от выражений.

Итак, если f(x,y) = x^2 + y^2, то частичная f по отношению к x будет 2x, частичная относительно y будет 2y.

Я написал функцию dinky с использованием метода конечных разностей, но у меня много проблем с точностью с плавающей запятой. Например, я заканчиваю 1.99234 вместо 2. Существуют ли библиотеки, которые поддерживают символическую дифференциацию? Любые другие предложения?

4b9b3361

Ответ 1

Я реализовал такие библиотеки на нескольких разных языках, но, к сожалению, не C. Если вы имеете дело только с полиномами (суммами, продуктами, переменными, константами и степенями), это довольно легко. Не слишком плохие функции Trig. Что-нибудь более сложное, и вам, вероятно, будет лучше не торопиться, чтобы овладеть чужой библиотекой.

Если вы решили сворачивать свои собственные, у меня есть несколько советов, которые упростят вашу жизнь:

  • Используйте неизменяемые структуры данных (чисто функциональные структуры данных) для представления выражений.

  • Используйте сборщик мусора Hans Boehm, чтобы управлять памятью для вас.

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

Если вы захотите вставить Lua в свой код C и выполнять свои вычисления там, я поместил код Lua в http://www.cs.tufts.edu/~nr/drop/lua. Одна из лучших особенностей заключается в том, что он может принимать символическое выражение, дифференцировать его и компилировать результаты в Lua. Вы, конечно, не найдете никакой документации: - (

Ответ 2

Получение числового дифференцирования "вправо" (в смысле минимизации ошибок) может быть довольно сложным. Чтобы начать работу, вы можете взглянуть на раздел Numerical Recipies на числовые производные.

Для бесплатных символических математических пакетов вы должны посмотреть GiNaC. Вы также можете проверить SymPy, автономный, символический математический пакет с чистым питоном. Вы обнаружите, что SymPy намного легче изучить, поскольку вы можете использовать его в интерактивном режиме из командной строки Python.

В коммерческом конце, как Mathematica, так и Maple имеют C API. Вам нужна установленная/лицензионная версия программы для использования библиотек, но оба могут легко выполнять тип символического дифференциации, который вы после.

Ответ 3

Если вы выполняете численное дифференцирование ( "оцените производную от f (x) при x = x0" ), и вы знаете, что вы заранее задаете уравнения (т.е. не пользовательский ввод), тогда я бы рекомендовал FADBAD ++. Это библиотека шаблонов С++ для решения числовых производных, используя Автоматическое дифференцирование. Это очень быстро и точно.

Ответ 4

Вы можете улучшить точность своего численного дифференцирования двумя способами:

  • Используйте меньшую дельту. Кажется, вы использовали значение около 1e-2. Начните с 1e-8 и испытайте, если вы получаете меньше болей или помогает. Очевидно, вы не можете приблизиться к точности машины - около 1e-16 для двойного.

  • Используйте центральные различия, а не прямые (или наоборот) различия. т.е. df_dx =(f(x+delta) - f(x-delta)) / (2.0*delta) По причинам, связанным с отменой более высоких членов усечения, ошибка в оценке центральных различий имеет порядок delta^2, а не дельта переходов. См. http://en.wikipedia.org/wiki/Finite_difference

Ответ 5

Если это действительно такая функция, которую вы хотите использовать, было бы достаточно просто написать библиотеку классов. Начните с одного термина, с коэффициентом и показателем. Имейте Полиномиальный, который будет состоять из Собрания Условий.

Если вы определяете интерфейс для интересующих математических методов (например, add/sub/mul/div/diffiate/integrate), вы просматриваете шаблон GoF Composite. И Term, и Polynomial будут реализовывать этот интерфейс. Полином просто перебирал каждый термин в своей коллекции.

Ответ 6

Конечно, было бы легче использовать существующий пакет, чем писать самостоятельно, но если вы решите написать свой собственный, и вы готовы потратить некоторое время на изучение некоторых темных уголков С++, вы можете использовать Boost.Proto от Boost, чтобы спроектировать ваш собственной библиотеки.

В принципе, Boost.Proto позволяет вам преобразовать любое допустимое выражение С++, например x * x + y * y в шаблон выражения - в основном представление дерева разбора этого выражения с использованием вложенных struct - и затем выполнить любое произвольное вычисление над этим деревом разбора в более позднее время, вызвав на нем proto::eval(). По умолчанию proto::eval() используется для оценки дерева, как если бы он выполнялся напрямую, хотя нет причин, по которым вы не могли бы изменить поведение каждой функции или оператора вместо символической производной.

Хотя это было бы чрезвычайно сложным решением вашей проблемы, тем не менее было бы намного проще, чем пытаться перевернуть ваши собственные шаблоны выражений, используя методы метапрограммирования шаблонов С++.

Ответ 7

Я создал такую ​​библиотеку на С++, здесь вы можете увидеть образец: https://github.com/MartinKosicky/auto_diff/blob/master/sample/main.cpp

Использование довольно простое, оно также поддерживает матрицы. Я использовал его для создания повторяющейся нейронной сети... Мне нужно было бы добавить умножение матрицы GPU, хотя

Ответ 8

Это отчасти отложено, поскольку оно относится к Lisp, а не к C/С++, но может помочь другим, которые ищут похожие задачи, или вы можете получить некоторые идеи о реализации чего-то подобного в C/С++ самостоятельно. В SICP есть несколько лекций по теме для lisp:

  • производные правила 3b
  • алгебраические правила 4a

В Lisp он довольно прямолинейный (и на других функциональных языках с мощным сопоставлением шаблонов и полиморфными типами). В C вам, вероятно, придется сильно использовать перечисления и структуры для получения той же мощности (не говоря уже о распределении/освобождении). Можно было однозначно закодировать то, что вам нужно в ocaml в течение часа - я бы сказал, что скорость печати является ограничивающим фактором. Если вам нужно C, вы можете называть ocaml из C (и наоборот).

Ответ 9

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

  • значение
  • массив производных значений по сравнению с независимыми переменными

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

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

усеченный ряд Тейлора - для этого доступны две библиотеки:

http://code.google.com/p/libtaylor/

http://www.foelsche.com/ctaylor

Ответ 10

Посмотрите на Theano, он поддерживает символическое дифференцирование (в контексте нейронных сетей). Проект является открытым исходным кодом, поэтому вы должны уметь видеть, как они это делают.

Ответ 11

Извините за то, что принесли это через 6 лет. Тем не менее, я искал такую ​​библиотеку для своего проекта, и я видел, что @eduffy предлагает FADBAD ++. Я прочитал документацию и вернулся к вашему вопросу. Я чувствую, что мой ответ будет полезен, поэтому следующий код для вашего дела.

#include <iostream>
#include "fadiff.h"

using namespace fadbad;

F<double> func(const F<double>& x, const F<double>& y)
{
    return x*x + y*y;
}

int main()
{
    F<double> x,y,f;     // Declare variables x,y,f
    x=1;                 // Initialize variable x
    x.diff(0,2);         // Differentiate with respect to x (index 0 of 2)
    y=1;                 // Initialize variable y
    y.diff(1,2);         // Differentiate with respect to y (index 1 of 2)
    f=func(x,y);         // Evaluate function and derivatives

    double fval=f.x();   // Value of function
    double dfdx=f.d(0);  // Value of df/dx (index 0 of 2)
    double dfdy=f.d(1);  // Value of df/dy (index 1 of 2)

    std::cout << "    f(x,y) = " << fval << std::endl;
    std::cout << "df/dx(x,y) = " << dfdx << std::endl;
    std::cout << "df/dy(x,y) = " << dfdy << std::endl;

    return 0;
}

Выходной сигнал

    f(x,y) = 2
df/dx(x,y) = 2
df/dy(x,y) = 2

Другой пример, скажем, что нас интересует первая производная от sin(). Аналитически это cos. Это здорово, потому что нам нужно сравнить истинную производную от данной функции и ее числовую копию, чтобы вычислить истинную ошибку.

#include <iostream>
#include "fadiff.h"

using namespace fadbad;

F<double> func(const F<double>& x)
{
    return sin(x);
}



int main()
{
    F<double> f,x;
    double dfdx;
    x = 0.0;
    x.diff(0,1);
    f = func(x);
    dfdx=f.d(0);


    for (int i(0); i < 8; ++i ){
        std::cout << "       x: " << x.val()        << "\n"
                  << "    f(x): " << f.x()          << "\n" 
                  << " fadDfdx: " << dfdx           << "\n"
                  << "trueDfdx: " << cos(x.val())   << std::endl;
        std::cout << "=========================="   << std::endl;

        x += 0.1;
        f = func(x);
        dfdx=f.d(0);
    }


    return 0;
}

Результат

       x: 0
    f(x): 0
 fadDfdx: 1
trueDfdx: 1
==========================
       x: 0.1
    f(x): 0.0998334
 fadDfdx: 0.995004
trueDfdx: 0.995004
==========================
       x: 0.2
    f(x): 0.198669
 fadDfdx: 0.980067
trueDfdx: 0.980067
==========================
       x: 0.3
    f(x): 0.29552
 fadDfdx: 0.955336
trueDfdx: 0.955336
==========================
       x: 0.4
    f(x): 0.389418
 fadDfdx: 0.921061
trueDfdx: 0.921061
==========================
       x: 0.5
    f(x): 0.479426
 fadDfdx: 0.877583
trueDfdx: 0.877583
==========================
       x: 0.6
    f(x): 0.564642
 fadDfdx: 0.825336
trueDfdx: 0.825336
==========================
       x: 0.7
    f(x): 0.644218
 fadDfdx: 0.764842
trueDfdx: 0.764842
==========================