convergence_analysis.m File Reference

Functions
id	i: (id folders)
	fprintf ('Processing folder:%s\n', folder_name)
	fprintf (' Reading file %d:%s\n', j, file_name)
Extract mesh size from file	h_values (i, j)
	fprintf (' h=%.6e\n', h_values(i, j))
Read data matrix from	file (skip comment lines starting with '#') data
Compute energy error	errors (i, j)
	fprintf (' h=%.6e, error=%.6e\n', h_values(i, j), errors(i, j))
end	fprintf ('\n')
end Sort data by mesh	size (finest to coarsest or vice versa) for i
	errors (i,:)
	convergence_rates (i)
	fprintf ('Folder %s:Convergence rate=%.3f\n', folders{i}, convergence_rates(i))
end Plot convergence results	figure ('Position', [100, 100, 800, 600])
Position for	annotation (midpoint in log space) h_mid
Add text annotation	text (h_mid, error_mid, sprintf('%.2f', local_slope),... 'FontSize', 9, 'Color', colors(i,:),... 'HorizontalAlignment', 'center',... 'BackgroundColor', 'white', 'EdgeColor', colors(i,:),... 'Margin', 2)
end end	xlabel ('Mesh size h(log scale)', 'FontSize', 12)
	ylabel ('Energy Error(log scale)', 'FontSize', 12)
	title ('Convergence Analysis', 'FontSize', 14)
	legend ('Location', 'best', 'FontSize', 10)
	set (gca, 'FontSize', 11)
Print summary table	fprintf ('\n===Convergence Summary===\n')
	fprintf ('%-15s %-15s %-15s\n', 'Folder', 'Conv. Rate', 'Error Range')
	fprintf ('%-15s %-15s %-15s\n', '------', '----------', '-----------')
otherwise	error ('Unknown folder name:%s', folder_name)

Variables
Convergence Analysis for Triangle and Quadrilateral Meshes	clear
close	all
	clc
Setup	folders = {'ca_tri', 'ca_quad'}
	Nfolders = 2
	Nfiles = 5
Initialize storage arrays	h_values = zeros(Nfolders, Nfiles)
mesh sizes	errors = zeros(Nfolders, Nfiles)
for	j
end Compute convergence rates	convergence_rates = zeros(Nfolders, 1)
for	i
hold	on
Define	markers = {'^-', 's-'}
triangle and square markers	colors = lines(Nfolders)
	error_mid = exp((log(errors(i,j)) + log(errors(i,j+1))) / 2)
hold	off
	end
	dt = dt_values(file_index)
Determine mesh name based on folder switch folder_name case ca_tri For triangular	meshes
Determine mesh name based on folder switch folder_name case ca_tri For triangular	mesh1_2
Determine mesh name based on folder switch folder_name case ca_tri For triangular etc	mesh_name = sprintf('mesh1_%d', file_index)
case ca_quad For quadrilateral	quad10x10
case ca_quad For quadrilateral etc	mesh_size = 5 * 2^(file_index-1)
end Construct full filename	filename
end function	h
	h_cell = regexp(text, 'MeshSize:\s*([0-9.e+-]+)', 'tokens', 'once')
end function	error
	H1velocity = data(:,5).^2
	L2density = data(:,7).^2
	UPWdensity = data(:,9).^2
	nu = 1
	rhomin = 1
	err = max(L2density) + rhomin * max(L2velocity) + dt * sum(UPWdensity) + (nu * dt / 2) * sum(H1velocity)

Function Documentation

annotation()

Position for annotation

(

midpoint in log

space

)

convergence_rates()

convergence_rates

(

i

)

error()

otherwise error	(	'Unknown folder name:%s'	,
		folder_name
	)

errors() [1/2]

Compute energy error errors	(	i	,
		j
	)

errors() [2/2]

errors	(	i	,
		:
	)

figure()

end Plot convergence results figure

(

'Position'

)

file()

Read data matrix from file

(

skip comment lines starting with '#'

)

fprintf() [1/9]

fprintf	(	'	h = %.6e,
		error	= %.6e\n',
		h_values(i, j)	,
		errors(i, j)
	)

fprintf() [2/9]

fprintf	(	'	h = %.6e\n',
		h_values(i, j)
	)

fprintf() [3/9]

fprintf	(	' Reading file %d:%s\n'	,
		j	,
		file_name
	)

fprintf() [4/9]

fprintf	(	'%-15s %-15s %-15s\n'	,
		'------'	,
		'----------'	,
		'-----------'
	)

fprintf() [5/9]

fprintf	(	'%-15s %-15s %-15s\n'	,
		'Folder'	,
		'Conv. Rate'	,
		'Error Range'
	)

fprintf() [6/9]

end fprintf

(

'\n'

)

fprintf() [7/9]

Print summary table fprintf

(

'\

n = ==Convergence Summary===\n'

)

fprintf() [8/9]

fprintf	(	'Folder %s:Convergence	rate = %.3f\n',
		folders{i}	,
		convergence_rates(i)
	)

fprintf() [9/9]

fprintf	(	'Processing folder:%s\n'	,
		folder_name
	)

h_values()

Extract mesh size from file h_values	(	i	,
		j
	)

i:()

id i: ( id  folders )

virtual

legend()

legend	(	'Location'	,
		'best'	,
		'FontSize'	,
		10
	)

set()

set	(	gca	,
		'FontSize'	,
		11
	)

size()

end Sort data by mesh size

(

finest to coarsest or vice

versa

)

text()

Add text annotation text	(	h_mid	,
		error_mid	,
		sprintf('%.2f', local_slope)	,
		... 'FontSize'	,
		9	,
		'Color'	,
		colors(i,:)	,
		... 'HorizontalAlignment'	,
		'center'	,
		... 'BackgroundColor'	,
		'white'	,
		'EdgeColor'	,
		colors(i,:)	,
		... 'Margin'	,
		2
	)

title()

title	(	'Convergence Analysis'	,
		'FontSize'	,
		14
	)

xlabel()

end end xlabel	(	'Mesh size h(log scale)'	,
		'FontSize'	,
		12
	)

ylabel()

ylabel	(	'Energy Error(log scale)'	,
		'FontSize'	,
		12
	)

Variable Documentation

all

close all

clc

clc

clear

Convergence Analysis for Triangle and Quadrilateral Meshes clear

colors

triangle and square markers colors = lines(Nfolders)

convergence_rates

end Compute convergence rates convergence_rates = zeros(Nfolders, 1)

dt

dt = dt_values(file_index)

end

end

Initial value:

========================================
%% PLACEHOLDER FUNCTIONS (to be implemented)
%% ========================================
 
function filename = construct_filename(folder_name, file_index)
    % Construct filename based on folder and file index
    % Format: folder/timehistory_GuermondQuartapelle_k0_meshname_dt=value.txt
    
    % Define dt values for each file index
    dt_values = [0.001, 0.0005, 0.00025, 0.000125, 0.000063]

err

err = max(L2density) + rhomin * max(L2velocity) + dt * sum(UPWdensity) + (nu * dt / 2) * sum(H1velocity)

error

error

Initial value:

= compute_energy_error(data)
    % TODO: Compute energy error from data
    % Based on your example: data(end,7) is L2_rho
    % You might need to compute energy norm differently
 
    L2velocity = data(:,4).^2

error_mid

error_mid = exp((log(errors(i,j)) + log(errors(i,j+1))) / 2)

errors

mesh sizes errors = zeros(Nfolders, Nfiles)

filename

else filename

Initial value:

= sprintf('%s/GuermondQuartapelle_k0_%s_dt=%.6f.txt', ...
                       folder_name, mesh_name, dt)

folders

Setup folders = {'ca_tri', 'ca_quad'}

h

h

Initial value:

= extract_mesh_size(filename)
    % Extract mesh size from file header
    % Looks for pattern: MeshSize: <value>
    text = fileread(filename)

H1velocity

H1velocity = data(:,5).^2

h_cell

h_cell = regexp(text, 'MeshSize:\s*([0-9.e+-]+)', 'tokens', 'once')

h_values

Initialize storage arrays h_values = zeros(Nfolders, Nfiles)

i

for i

Initial value:

= 1:Nfolders
    % Fit: log(error) = slope * log(h) + intercept
    % Using polyfit on log-log scale
    p = polyfit(log(h_values(i,:)), log(errors(i,:)), 1)

j

end for j

Initial value:

= 1:Nfiles
        % Construct file name (you'll need to adjust the pattern)
        file_name = construct_filename(folder_name, j)

L2density

L2density = data(:,7).^2

markers

Define square for quad markers = {'^-', 's-'}

mesh1_2

Determine mesh name based on folder switch folder_name case ca_tri For triangular mesh1_2

mesh_name

mesh_name = sprintf('mesh1_%d', file_index)

mesh_size

case ca_quad For quadrilateral etc mesh_size = 5 * 2^(file_index-1)

meshes

case ca_quad For quadrilateral meshes

Nfiles

Nfiles = 5

Nfolders

Nfolders = 2

nu

nu = 1

off

hold off

on

grid on

quad10x10

case ca_quad For quadrilateral quad10x10

rhomin

rhomin = 1

UPWdensity

UPWdensity = data(:,9).^2