Question about Integrating a 2D matrix that depends on 3 variables

I have a question about numeric integration specfically trying to compute the integral form below. I have a function in a general form of

where is a first order Bessel function of the 0th order and i is the complex number. If I understand this correctly - to integrate this function numerically over two variables k and z I put this function, f(r,L) into a for loop for every value of r and compute the double integral. The issue is I don't think I did this correctly code wise since I am trying to avoid any use of symbolic equations which most other posts use.

for NN = 1:length(r)
    Val = r(NN)
    for MM = 1:length(z)
      Temp(MM) = (integral(@(K) K.*MAS(K).*f(K,z(MM),Val),0,1));
    end
    Intval(NN) = sum(Temp);
end

My idea is to loop through all values of r, in this case a 1D vector, and then do the same for each value of z and then sum up all values for each value of z at a particular r which is the second integral. Is this the correct way to go about a numerical double integral?

If anyone has insight or a way to do this using the function integral2 it would be greatly appriciated since I am really having a problem solving this integral correctly.

Answer 1

[darova]

I used numerical calculations

kk = linspace(0,1,20);
zz = linspace(0,10,20);
rr = linspace(0,15,20);
dk = kk(2)-kk(1);
dz = zz(2)-zz(1);
[K,Z,R] = meshgrid(kk,zz,rr);
b = 0.2;
J0 = 1.5;
C = 1.2;
A = -2;
G = 1.1;
f1 = 1 - exp(-b.*(K.*Z).^2).*J0.*([1-(1-b).*Z-1i*C*Z].*K.*R);
F = A*exp(-b*K.*G.*f1);

RES = sum(sum(F,1),2);
RES = RES(:)*dk*dz;
plot(rr,RES)

This how course integration works

x = linspace(0,10,20);
y = sin(x);
dx = x(2)-x(1);
disp('course integration')
sum(y)*dx
disp('trapz integration')
trapz(x,y)