Untitled

bzeros = zeros(num_roots, num_functions); % define the matrix bzeros, size(bzeros) = [5, 6] 
% We just need to finding the root -> J_n(x) don't need to be divided by pi -> omit pi
integrand = @(theta, n, x) cos(x.*sin(theta)-n*theta); % define the function in integral
for n = 0:num_functions-1
    J_n = @(x) integral(@(theta) integrand(theta, n, x), 0, pi); % define J_n(x) without pi
    for i = 1:num_roots
        bzeros(i, n+1) = fzero(J_n, zeros_guess(i, n+1)); % note the index in matlab
    end
end