Is the computation of Beta correct ? Please could anyone correct the codes based on the equations ?

My code for Beta along with the values :

function beta_n_values = Beta(n, d)
    % Beta function to calculate Beta^n(k - d)
    % n: Degree of the polynomial
    % d: Interpolation point as a random variable with |d| < 0.5

    % Calculate the lower and upper bounds of k
    k0 = ceil(d - (n + 1) / 2);
    k1 = floor(d + (n + 1) / 2);
    
    % Define the range of k values
    k_values = k0:k1;

    % Pre-compute the coefficients a_p(k) for each k using the calculate_ap function
    a_p_values = calculate_ap(n, d); % Make sure the calculate_ap function is defined

    % Initialize an array to store the values of Beta^n(k - d)
    beta_n_values = zeros(size(k_values));

    % Loop over each k in k_values to compute Beta^n(k - d)
    for k_idx = 1:length(k_values)
        k = k_values(k_idx);
        
        % Initialize the sum for this k
        beta_n_k = 0;
        
        % Loop over p from 0 to n
        for p = 0:n
            % Get the coefficient a_p(k) from the precomputed matrix a_p_values
            a_p_k = a_p_values(p+1, k_idx); % MATLAB indexing starts from 1, so use p+1
            
            % Add the term a_p(k) * d^p to the sum
            beta_n_k = beta_n_k + a_p_k * d^p;
        end
        
        % Store the computed value of Beta^n(k - d)
        beta_n_values(k_idx) = beta_n_k;
    end

    % Display the computed Beta^n(k - d) values
    disp('Computed Beta^n(k - d) values:');
    disp(beta_n_values);
end

% Input the degree of the polynomial

n = input('Enter the degree (n): '); % Example: n = 3;

% Generate a random interpolation point d with |d| < 0.5

num_points = 1;

d = -0.5 + (0.5 - (-0.5)) * rand(num_points, 1);

a_p_values = calculate_ap(n, d);

disp('The computed a_p^k values are:');

disp(a_p_values);

beta_n_values = Beta(n, d);

disp('Computed Beta^n(k - d) values:');

disp(beta_n_values);

The solution for a_p and beta is :

>> SplineCalculation

Enter the degree (n): 3

Calculated a_p(k) values:

0.6667 0.6667 0.6667 0.6667

-1.0000 -1.0000 -1.0000 -1.0000

0.5000 0.5000 0.5000 0.5000

0 0 0 0

Computed Beta^n(k - d) values:

0.3432 0.3432 0.3432 0.3432

The code for function a_p:

function a_p_values = calculate_ap(n, d)
% n: Degree of the polynomial
% d: Interpolation point as a random variable with |d| < 0.5

% Define the range for p (0 to n)
p_values = 0:n;

% Calculate the lower and upper bounds of k
k0 = ceil(d - (n + 1) / 2);
k1 = floor(d + (n + 1) / 2);

% Define the range for k
k_values = k0:k1;

% Initialize a matrix to store a_p(k) values
a_p_values = zeros(length(p_values), length(k_values));

% Loop through each value of p to compute a_p for each k
for p_idx = 1:length(p_values)
    p = p_values(p_idx);
    
    % Calculate the binomial coefficient (n choose p)
    binomial_np = nchoosek(n, p);
    
    % Calculate the factorial term 1/n!
    n_factorial_term = 1 / factorial(n);
    
    % Loop through each value of k
    for k_idx = 1:length(k_values)
        k = k_values(k_idx);
        
        % Initialize the sum for this combination of p and k
        sum_result = 0;
        
        % Sum over m from 0 to n+1
        for m = 0:n+1
            % Calculate the binomial coefficient (n+1 choose m)
            binomial_n1m = nchoosek(n+1, m);
            
            % Calculate the term (-1)^m
            sign_term = (-1)^m;
            
            % Calculate the term (n+1)/2 - m
            difference_term = (n+1) / 2 - m;
            
            % Calculate the term (difference_term)^(n-p)
            power_term = difference_term^(n - p);
            
            % Calculate the Heaviside function mu(p - m + (n+1)/2)
            mu_term = heaviside(p - m + (n+1)/2);
            
            % Add the current term to the sum
            sum_result = sum_result + binomial_n1m * sign_term * power_term * mu_term;
        end
        
        % Calculate a_p(k) for the current p and k
        a_p_values(p_idx, k_idx) = n_factorial_term * binomial_np * sum_result;
    end
end

% Display the calculated a_p(k) values
disp('Calculated a_p(k) values:');
disp(a_p_values);
end

Is the code correct cause all are convolutions basically

[Torsten]

As already noted in a previous question, i don't see where a parameter k appears in your definition for a_p(k).

[Rohitashya]

Is it possible to make changes this is a definition from an existing paper .I'm trying to compute in code.Will yoube able to do it ?

[Torsten]

Will you be able to do it ?

How should I if you don't tell us what k and mu mean in your formula ?

[Rohitashya]

All the details are in the code

[Rohitashya]

You can modify the code it's given with relevant data