Burgers-FLRW Model with WENO and Lax-Friedrichs Flux

[Kadir Çar]

Hello everyone, as I have mentioned in the title section I want to implement FLRW equation to clawpack package and get the output from clawpack.  I already have the matlab code for this equation and it has f_plus and f_minus fluxes, for that example I will provide the matlab code and asking from you guys to direct me how can I write the riemann solver for this problem since it has two fluxes. I am a bit confused and if there is any example like this matlab code implemented on clawpack  I can check and try to implement it to my example but  I have checked and haven't seen an example like this.  My matlab code is provided below:

clear; clc;

%The best

%for k=0, k=1, k=-1

%FLRW model 16 april

% Parameters

J = 200; % Number of spatial points

N = 256; % Number of time steps

Rmin = 0; Rmax = 1; % Spatial domain

deltaR = (Rmax - Rmin) / (J - 1);

Tmax = 1;

deltaT = Tmax / N;

r = linspace(Rmin, Rmax, J); % Spatial grid

% Scale factor and curvature parameter

a = @(t) t^2; % Scale factor

k = -1; % Curvature parameter

% Initial condition (Riemann problem)

v = zeros(N+1, J);

threshold = 0.5; % Threshold for the Riemann problem

for i = 1:J

if r(i) < threshold

v(1, i) = 1; % Left state

else

v(1, i) = 0; % Right state

end

end

% Time evolution using a simplified WENO scheme and Lax-Friedrichs flux

for n = 1:N

vn = v(n, :);

vn_new = vn; % Initialize the new time step

for j = 3:J-2

at = a((n-1)*deltaT + 1);

dot_a = 2 * ((n-1)*deltaT + 1);

% Lax-Friedrichs Flux, simplified

f_minus = 0.5 * ((1 - k*r(j-1)^2)^0.5 / at) * (vn(j-1)^2);

f_plus = 0.5 * ((1 - k*r(j+1)^2)^0.5 / at) * (vn(j+1)^2);

F_j_minus_half = 0.5 * (f_plus + f_minus) - 0.5 * max(abs(vn(j+1)), abs(vn(j-1))) * (vn(j+1) - vn(j-1));

f_minus = 0.5 * ((1 - k*r(j)^2)^0.5 / at) * (vn(j)^2);

f_plus = 0.5 * ((1 - k*r(j+2)^2)^0.5 / at) * (vn(j+2)^2);

F_j_plus_half = 0.5 * (f_plus + f_minus) - 0.5 * max(abs(vn(j+2)), abs(vn(j))) * (vn(j+2) - vn(j));

% Source term

source_term = -dot_a / at * vn(j) * (1 - vn(j)^2) - k*r(j) / at * vn(j)^2 / 2 * (1 - k*r(j)^2)^(-0.5);

% Update rule

vn_new(j) = vn(j) - deltaT / deltaR * (F_j_plus_half - F_j_minus_half) + deltaT * source_term;

end

% Boundary conditions

vn_new(1) = vn_new(2);

vn_new(J) = vn_new(J-1);

v(n+1, :) = vn_new;

end

% Compute maximum velocity for CFL condition

max_velocity = max(max(abs(v)));

CFL_condition = max_velocity * deltaT / deltaR;

fprintf('CFL Condition = %.3f\n', CFL_condition);

% Visualization

figure;

hold on;

xlabel('r');

ylabel('v(r, t)');

title('Burgers-FLRW Model with WENO and Lax-Friedrichs Flux');

xlim([Rmin Rmax]);

colormap jet;

colors = jet(N+1);

intervals = round(linspace(1, N+1, min(10, N))); % Adjust the number of plots

for idx = intervals

plot(r, v(idx, :), 'Color', colors(idx, :), 'LineWidth', 2);

pause(0.3); % Optional: to see the plotting in progress

end

legend(arrayfun(@(x) sprintf('Time = %.2f', (x-1) * deltaT), intervals, 'UniformOutput', false), 'Location', 'best');

grid on;

ylim([-0.2, 1.2]);

hold off;

[David Ketcheson]

Hi,

This sounds like an interesting project!  Unfortunately, I think it's unlikely that anyone will do all the work for you in the way you have suggested.  It's also hard to answer your questions since we don't know the equations that you're solving.  I'd be happy to answer more specific questions, though, if you decide to do this yourself and can explain the equations you wish to solve.

David