Ошибка загрузки при попытке передать сложную функцию в правило Симпсона - PullRequest
Новая
       31

Ошибка загрузки при попытке передать сложную функцию в правило Симпсона

0 голосов
/ 17 апреля 2019

Я написал метод, который аппроксимирует определенный интеграл по правилу составного Симпсона. 

#=
f  integrand
a  lower integration bound
b  upper integration bound
n  number of iterations or panels

h     step size
=#
function simpson(f::Function, a::Number, b::Number, n::Number)
   n % 2 == 0 || error("`n` must be even")
   h = (b - a) / n
   s = f(a) + f(b)
   s += 4*sum(f(a .+ collect(1:2:n) .* h))
   s += 2*sum(f(a .+ collect(2:2:n-1) .* h))
   return h/3 * s
end

Для «простых» функций, таких как e^(-x^2), работает функция simpson. 

Input: simpson(x -> simpson(x -> exp.(-x.^2), 0, 5, 100)
Output: 0.8862269254513949

Однако для более сложной функции f(x)

gArgs(x) = (30 .+ x, 0)
f(x) = exp.(-x.^2) .* maximum(generator.(gArgs.(x)...)[1])

, где generator(θ, plotsol) - это функция, которая принимает дефект θ в процентах и ​​логическое значение plotsol (либо 0, либо 1), которое определяет, следует ли строить график генератора, и возвращает вектор с намагниченностью в определенных точках. в генераторе.

Когда я пытаюсь вычислить интеграл, запустив следующий код 

gArgs(x) = (30 .+ x, 0)
f(x) = exp.(-x.^2) .* maximum(generator.(gArgs.(x)...)[1])
println(simpson(x -> f(x), 0, 5, 10))

Я сталкиваюсь с ошибкой MethodError: no method matching generator(::Float64). С небольшими вариантами выражения для f(x) я сталкиваюсь с различными ошибками, такими как DimensionMismatch("array could not be broadcast to match destination") и InexactError: Bool(33.75). В конце концов, я думаю, что причина ошибки сводится к тому, что я не могу понять, как правильно ввести выражение для подынтегрального выражения f(x). Может ли кто-нибудь помочь мне разобраться, как правильно ввести f(x)? Дайте мне знать, если что-то неясно в моем вопросе.

1 Ответ

1 голос
/ 17 апреля 2019

Для массива x, gArgs.(x) возвращает массив Tuple s, и вы пытаетесь транслировать через массив кортежей.Но поведение вещания с кортежами немного отличается.Кортежи не рассматриваются как один элемент, и они сами транслируются. 

julia> println.(gArgs.([0.5, 1.5, 2.5, 3.5, 4.5])...)
30.531.532.533.534.5
00000

Это не то, что вы ожидали, не так ли?

Вы также можете увидеть проблему в следующем примере; 

julia> (2, 5) .!= [(2, 5)]
2-element BitArray{1}:
 true
 true

Я считаю, f - это функция, которая на самом деле берет скаляр и возвращает скаляр.Вместо того, чтобы f работать с массивами, вы должны оставить вещание вызывающей стороне.Скорее всего, вам будет лучше реализовать f поэлементно.Это более удобный для Джулии способ работы, который значительно облегчит вашу работу.

Тем не менее, я считаю, что ваша реализация должна работать со следующими модификациями , если у вас нет ошибки в generator. 

function simpson(f::Function, a::Number, b::Number, n::Number)
   n % 2 == 0 || error("`n` must be even")
   h = (b - a) / n
   s = f(a) + f(b)
   s += 4*sum(f.(a .+ collect(1:2:n) .* h)) # broadcast `f`
   s += 2*sum(f.(a .+ collect(2:2:n-1) .* h)) # broadcast `f`
   return h/3 * s
end

# define `gArg` and `f` element-wise and `generator`, too.
gArgs(x) = (30 + x, 0) # get rid of broadcasting dot. Shouldn't `0` be `false`?
f(x) = exp(-x^2) * maximum(generator(gArgs(x)...)[1]) # get rid of broadcasting dots
println(simpson(f, 0, 5, 10)) # you can just write `f`

Вы также должны определить функцию generator поэлементно.

...