Блог

Главная - Вопросы по программированию - Решение с двойным интегралом с использованием scipy.integrate.nquad не соответствует integrate.dblquad

Решение с двойным интегралом с использованием scipy.integrate.nquad не соответствует integrate.dblquad

#python #math #scipy #integral #mismatch

#python #математика #scipy #интеграл #несоответствие

Вопрос:

Первая функция в приведенном ниже коде использует двойное интегрирование с scipy.integrate.dblquad для вычисления дифференциальной энтропии, c*np.log(c) , функции плотности связки c , которая имеет один параметр зависимости, theta , обычно положительный.

введите описание изображения здесь

Вторая функция в приведенном ниже коде пытается решить ту же проблему, что и выше, но с использованием решателя с несколькими интегралами scipy.integrate.nquad 

 from scipy import integrate
import numpy as np

def dblquad_(theta):
    "Double integration"
    c = lambda v, u: ((1 theta)*(u*v)**(-1-theta)) * (u**(-theta) v**(-theta)-1)**(-1/theta-2)
    return -integrate.dblquad(
        lambda u,v: c(v,u)*np.log(c(v,u)), 
        0, 1, lambda u: 0, lambda u: 1
        )[0]

def nquad_(n,theta):
    "Multiple integration"
    c = lambda *us: ((1 theta)*np.prod(us)**(-1-theta)) * (np.sum(np.power(us,-theta))-1)**(-1/theta-2)
    return -integrate.nquad(
        func   = lambda *us : c(*us)*np.log(c(*us)), 
        ranges = [(0,1) for i in range(n)],
        args   = (theta,) 
        )[0] 

n=2
theta = 1
print(dblquad_(theta))
print(nquad_(n,theta))

Функция на dblquad основе выдает ответ -0.7127 , в то время nquad как выдает -0.5823 и заметно занимает больше времени. Почему решения разные, хотя я установил оба для решения n=2 многомерной задачи?

Ответ №1:

С n theta учетом предоставленных вами значений и выходных данных вашего кода:

 -0.1931471805597395
0.17055845832017144,

не -0.7127 » и -0.5823 » .

Первое значение ( -0.1931471805597395 ) является правильным (вы можете проверить это здесь сами).

Проблема nquad_ заключается в обработке theta параметра. Спасибо @mikuszefski за объяснение; Я копирую его здесь для ясности:

nquad переходит lambda к функции по мере необходимости. Лямбда-выражение запрограммировано таким образом, что оно принимает произвольное количество аргументов, поэтому оно с радостью принимает его и помещает в список степеней и сумм. Следовательно, вы не получаете, например 1/u**t 1/v**t -1 , но 1/u**t 1/v**t 1/t**t -1 . Вызов функции просто не соответствует предполагаемому использованию функции. Если вы вместо этого напишете us[0]**() us[1]**() - 1 , это сработает.

Вот измененный код:

 from scipy import integrate
import numpy as np

def dblquad_(theta):
    "Double integration"
    c = lambda v, u: ((1 theta)*(u*v)**(-1-theta)) * (u**(-theta) v**(-theta)-1)**(-1/theta-2)
    return -integrate.dblquad(
        lambda u,v: c(v,u)*np.log(c(v,u)),
        0, 1, lambda u: 0, lambda u: 1
        )[0]

def nquad_(n,theta):
    "Multiple integration"
    c = lambda *us: ((1 theta)*np.prod((us[0], us[1]))**(-1-theta)) * (np.sum(np.power((us[0], us[1]),-theta))-1)**(-1/theta-2)
    return -integrate.nquad(
        func   = lambda *us : c(*us)*np.log(c(*us)),
        ranges = [(0,1) for i in range(n)],
        args=(theta,)
        )[0]

n=2
theta = 1
print(dblquad_(theta))
print(nquad_(n,theta))

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

 -0.1931471805597395
-0.1931471805597395

Комментарии:

1. Учитывая, что вы обнаружили, что существует определенное решение для двойного интеграла с использованием параметрического ввода c(theta) , не могли бы вы сказать, что это указывает на то, что решение в закрытой форме может быть получено для показанной математической формулы?

2. Кстати, параметр theta не используется неправильно. nquad передает его в функцию по мере необходимости. Он lambda запрограммирован таким образом, что принимает произвольное количество аргументов, поэтому с радостью принимает его и помещает в список степеней и сумм. Следовательно, вы не получаете, например 1/u**t 1/v**t -1 , но 1/u**t 1/v**t 1/t**t -1 . Вызов функции просто не соответствует предполагаемому использованию функции. Если вы вместо этого напишете us[0]**() us[1]**() - 1 , это сработает.

3. @develarist интеграл от c(u,v) на самом деле имеет замкнутую форму. c * log c скорее всего, нет.

4. в какой статье вы видели c решение в закрытой форме для его интеграла?

5. математически, почему nquad_(n=3,theta) дает точно такое же решение, что и двумерная связка nquad_(n=2,theta) , при этом вычисление занимает в два раза больше времени? Разве на многомерные связки не должен по-разному влиять один и тот же уровень theta ? Если нет, означает ли это, что решение для двумерных связок уже даст мне результаты связки более высокого размера?