What would be right approach to compute this 1D convolution ?

I am trying to generate some data which is a result of 1D convolution. These equations are related to acoustics of a moving medium.

This is the Green's function:

This is the Green's function differentiated w.r.t to time "t":

I am interested in calculating the pressure which is given by this convolution product:

I am interested in plotting the results on x-axis bounded as . I have included the script which has given invalid values. I am wondering what could be going wrong. Here is the plot which I want to generate.

Can someone please provide some tips to troubleshoot this?

This my script:

close all
clear 
clc

% frequency 
freq = 1/30;

% angular frequency
omega = 2*pi*(1/30);
 
% time
t=270;

% speed of sound
c = 290;

% constants 
alpha = log(2)/2;
epsilon = 0.5;

% Mach number
M = 0.5;

% wavenumber
k = omega/c;

% discretize the x axis bounded from -100 to 100
x = -100:0.5:100;

% hankel function
hankel_function = besselh(0, k*x/(1-M^2));

% exponential function
exp_function = -(1j*M*k*x/(1-M^2)) - (1j*omega*t);

% Green's function - complex quantity
G = (1j/4) * (1/( c^2*sqrt(1-M^2) )) * hankel_function .* exp_function;

% Green's function differentiated w.r.t time "t" - complex quantity
diffG = -1j*omega * G;

% spatial portion of the source term = f(x)
f = epsilon * exp(-alpha*x);

% convolution of "f" and real part of "diffG" should give pressure p(x,t)
p = conv( f, real(diffG),'same')*0.5;       

[Abhaya]

".. the script which has given invalid values."

The hankel_function has value NaN at x = 0, as besselh(0,0) returns NaN.

[Walter Roberson]

The convolution then spreads the NaN out to the entire result.

[ishan]

Yes that's right. I am wondering then how should this convolution be approached because these expressions are from a well-known paper in the field of acoustics.

[Walter Roberson]

There is not much you can do. limit(besselh(0,x)) is -1-1i*inf from the left and +1-1i*inf from the right. They have the same complex component, but different real components.

It is true that you are multiplying by exp_function and exp_function involves x so there is the potential for the values canceling out. However, exp_function has an additional constant term, and that constant term times +/-1 - 1i*inf is still going to be disjoint.

[Torsten]

Maybe a link to the "well-known paper in the field of acoustics" can help.