Я пытаюсь реализовать 2D-образы Лапласа с помощью рядов Фурье в R. Мой код работал в случае единственной инверсии, но не для двойной инверсии. Любая помощь будет очень высоко ценится. Спасибо!
Ниже приведен код в R для простой функции f(x,t) = exp(-3(x+t))
с двойным преобразованием Лапласа F(p,s) = 1/((p+3)*(s+3))
.
# General implementation of 2D Laplace Inversion
# How to use this?
# 1. Change the Laplace Transform function (F)
# 2. Change the Original function(f_1)
# Double Laplace Transform function, parameters p and s
F <- function(p,s) {1/((p+3)*(s+3))}
# Adjust parameters here
a <- 0.1
n <- 1000
# Note that the approximation works accurately up to x=T/2 only
T <- 10
# First inversion - p
F_1 <- function(x,s,k) {Re(F(complex(real=a, imaginary=k*pi/T),s))*cos(k*pi*x/T)}
# Summing over k from 1 to n
F_2 <- function(x,s) {sapply(1:n, F_1, x=x, s=s)}
# Approximate Function using Fourier Series
G <- function(x,s) {2*exp(a*x)/T*(0.5*Re(F(a,s))+sum(F_2(x,s)))}
# Second inversion - s
G_1 <- function(x,t,k) {Re(G(x,complex(real=a, imaginary=k*pi/T)))*cos(k*pi*t/T)}
# Summing over l from 1 to n
G_2 <- function(x,t) {sapply(1:n, G_1, x=x, t=t)}
# Approximate Function using Fourier Series
f2 <- function(x,t) {2*exp(a*t)/T*(0.5*Re(G(x,a))+sum(G_2(x,t)))}
# The other way round.
# First inversion - s
F_1 <- function(t,p,k) {Re(F(p, complex(real=a, imaginary=k*pi/T)))*cos(k*pi*t/T)}
# Summing over k from 1 to n
F_2 <- function(t,p) {sapply(1:n, F_1, t=t, p=p)}
# Approximate Function using Fourier Series
G <- function(t,p) {2*exp(a*t)/T*(0.5*Re(F(p,a))+sum(F_2(t,p)))}
# Second inversion - p
G_1 <- function(x,t,k) {Re(G(t,complex(real=a, imaginary=k*pi/T)))*cos(k*pi*x/T)}
# Summing over l from 1 to n
G_2 <- function(x,t) {sapply(1:n, G_1, x=x, t=t)}
# Approximate Function using Fourier Series
f2 <- function(x,t) {2*exp(a*x)/T*(0.5*Re(G(a,t))+sum(G_2(x,t)))}
# Original Function
f1 <- function(x,t) {exp(-3*(x+t))}