WIP: Sc/flux ruedaramirez etal normal direction

[SimonCan]

Added fluxes using normal_direction instead of orientation to IdealGlmMhdMultiIonEquations2D.

[github-actions]

Review checklist

This checklist is meant to assist creators of PRs (to let them know what reviewers will typically look for) and reviewers (to guide them in a structured review process). Items do not need to be checked explicitly for a PR to be eligible for merging.

Purpose and scope

• The PR has a single goal that is clear from the PR title and/or description.
• All code changes represent a single set of modifications that logically belong together.
• No more than 500 lines of code are changed or there is no obvious way to split the PR into multiple PRs.

Code quality

• The code can be understood easily.
• Newly introduced names for variables etc. are self-descriptive and consistent with existing naming conventions.
• There are no redundancies that can be removed by simple modularization/refactoring.
• There are no leftover debug statements or commented code sections.
• The code adheres to our conventions and style guide, and to the Julia guidelines.

Documentation

• New functions and types are documented with a docstring or top-level comment.
• Relevant publications are referenced in docstrings (see example for formatting).
• Inline comments are used to document longer or unusual code sections.
• Comments describe intent ("why?") and not just functionality ("what?").
• If the PR introduces a significant change or new feature, it is documented in NEWS.md with its PR number.

Testing

• The PR passes all tests.
• New or modified lines of code are covered by tests.
• New or modified tests run in less then 10 seconds.

Performance

• There are no type instabilities or memory allocations in performance-critical parts.
• If the PR intent is to improve performance, before/after time measurements are posted in the PR.

Verification

• The correctness of the code was verified using appropriate tests.
• If new equations/methods are added, a convergence test has been run and the results
  are posted in the PR.

Created with ❤️ by the Trixi.jl community.

[github-actions]

Remaining comments which cannot be posted as a review comment to avoid GitHub Rate Limit

JuliaFormatter

---

[JuliaFormatter] _{reported by reviewdog 🐶}

Trixi.jl/src/equations/ideal_glm_mhd_multiion_2d.jl

Line 1176 in a29b2e2

(rho_e_rr - 0.5f0 * rho_rr * vel_norm_rr - 0.5f0 * mag_norm_rr -

---

[JuliaFormatter] _{reported by reviewdog 🐶}

Trixi.jl/src/equations/ideal_glm_mhd_multiion_2d.jl

Line 1182 in a29b2e2

vk3_plus_ll[k] * B3_ll

---

[JuliaFormatter] _{reported by reviewdog 🐶}

Trixi.jl/src/equations/ideal_glm_mhd_multiion_2d.jl

Line 1184 in a29b2e2

vk3_plus_rr[k] * B3_rr

---

[JuliaFormatter] _{reported by reviewdog 🐶}

Trixi.jl/src/equations/ideal_glm_mhd_multiion_2d.jl

Lines 1213 to 1214 in a29b2e2

	f2 = (f1 * v1_avg + p_mean) * normal_direction[1] + f1 * v1_avg * normal_direction[2]
	f3 = f1 * v2_avg * normal_direction[1] + (f1 * v2_avg + p_mean) * normal_direction[2]

---

[JuliaFormatter] _{reported by reviewdog 🐶}

Trixi.jl/src/equations/ideal_glm_mhd_multiion_2d.jl

Line 1219 in a29b2e2

vk1_plus_rr[k] * mag_norm_rr)

---

[JuliaFormatter] _{reported by reviewdog 🐶}

Trixi.jl/src/equations/ideal_glm_mhd_multiion_2d.jl

Line 1221 in a29b2e2

vk2_plus_rr[k] * mag_norm_rr)

---

[JuliaFormatter] _{reported by reviewdog 🐶}

Trixi.jl/src/equations/ideal_glm_mhd_multiion_2d.jl

Line 1224 in a29b2e2

f2 * v1_avg + f3 * v2_avg + f4 * v3_avg

---

[JuliaFormatter] _{reported by reviewdog 🐶}

Trixi.jl/src/equations/ideal_glm_mhd_multiion_2d.jl

Lines 1227 to 1239 in a29b2e2

	0.5f0 * v1_plus_mag_avg * normal_direction[1] - 0.5f0 * v2_plus_mag_avg + normal_direction[2] +
	(B1_avg * normal_direction[1] + B2_avg * normal_direction[2]) * vel_dot_mag_avg # Same terms as in Derigs (but with v_plus)
	+ f9 * psi_avg - equations.c_h * (psi_B1_avg * normal_direction[1] + psi_B1_avg * normal_direction[2]) # GLM term
	+
	0.5f0 * (vk1_plus_avg * normal_direction[1] + vk2_plus_avg * normal_direction[2]) * mag_norm_avg -
	vk1_plus_avg * (B1_avg * normal_direction[1] + B2_avg * normal_direction[2]) * B1_avg -
	vk2_plus_avg * (B1_avg * normal_direction[1] + B2_avg * normal_direction[2]) * B2_avg -
	vk3_plus_avg * (B1_avg * normal_direction[1] + B2_avg * normal_direction[2]) * B3_avg # Additional terms related to the Lorentz non-conservative term (momentum eqs)
	-
	B2_avg * (vk1_minus_avg * B2_avg - vk2_minus_avg * B1_avg) * normal_direction[1] -
	B3_avg * (vk1_minus_avg * B3_avg - vk3_minus_avg * B1_avg) * normal_direction[1] -
	B1_avg * (vk2_minus_avg * B1_avg - vk1_minus_avg * B2_avg) * normal_direction[2] -
	B3_avg * (vk2_minus_avg * B3_avg - vk3_minus_avg * B2_avg) * normal_direction[2]) # Terms related to the multi-ion non-conservative term (induction equation!)

---

[JuliaFormatter] _{reported by reviewdog 🐶}

Trixi.jl/src/equations/ideal_glm_mhd_multiion_2d.jl

Line 1293 in a29b2e2

v_ll = max(v_ll, abs((rho_v1 * normal_direction[1] + rho_v2 * normal_direction[2]) / rho))

---

[JuliaFormatter] _{reported by reviewdog 🐶}

Trixi.jl/src/equations/ideal_glm_mhd_multiion_2d.jl

Line 1295 in a29b2e2

v_rr = max(v_rr, abs((rho_v1 * normal_direction[1] + rho_v2 * normal_direction[2]) / rho))

---

[JuliaFormatter] _{reported by reviewdog 🐶}

Trixi.jl/src/equations/ideal_glm_mhd_multiion_2d.jl

Lines 1418 to 1420 in a29b2e2

	sqrt(0.5f0 * (a_square + b_square) +
	0.5f0 *
	sqrt((a_square + b_square)^2 - 4 * a_square * (b1 *normal_direction[1] + b2 * normal_direction[2])^2)))

[codecov]

Codecov Report

Attention: Patch coverage is 0% with 248 lines in your changes missing coverage. Please review.

Project coverage is 96.26%. Comparing base (e33e645) to head (16015ab).

Files with missing lines	Patch %	Lines
src/equations/ideal_glm_mhd_multiion_2d.jl	0.00%	248 Missing ⚠️

Additional details and impacted files

@@            Coverage Diff             @@
##             main    #2446      +/-   ##
==========================================
- Coverage   96.83%   96.26%   -0.57%     
==========================================
  Files         504      504              
  Lines       41813    42061     +248     
==========================================
  Hits        40487    40487              
- Misses       1326     1574     +248     

Flag	Coverage Δ
unittests	96.26% <0.00%> (-0.57%)	⬇️

Flags with carried forward coverage won't be shown. Click here to find out more.

☔ View full report in Codecov by Sentry.
📢 Have feedback on the report? Share it here.

🚀 New features to boost your workflow:

• ❄️ Test Analytics: Detect flaky tests, report on failures, and find test suite problems.

[github-actions]

Remaining comments which cannot be posted as a review comment to avoid GitHub Rate Limit

JuliaFormatter

---

[JuliaFormatter] _{reported by reviewdog 🐶}

Trixi.jl/src/equations/ideal_glm_mhd_multiion_2d.jl

Line 1220 in 31cc07d

vk1_plus_rr[k] * mag_norm_rr)

---

[JuliaFormatter] _{reported by reviewdog 🐶}

Trixi.jl/src/equations/ideal_glm_mhd_multiion_2d.jl

Line 1222 in 31cc07d

vk2_plus_rr[k] * mag_norm_rr)

---

[JuliaFormatter] _{reported by reviewdog 🐶}

Trixi.jl/src/equations/ideal_glm_mhd_multiion_2d.jl

Line 1225 in 31cc07d

f2 * v1_avg + f3 * v2_avg + f4 * v3_avg

---

[JuliaFormatter] _{reported by reviewdog 🐶}

Trixi.jl/src/equations/ideal_glm_mhd_multiion_2d.jl

Lines 1228 to 1240 in 31cc07d

	0.5f0 * v1_plus_mag_avg * normal_direction[1] - 0.5f0 * v2_plus_mag_avg + normal_direction[2] +
	(B1_avg * normal_direction[1] + B2_avg * normal_direction[2]) * vel_dot_mag_avg # Same terms as in Derigs (but with v_plus)
	+ f9 * psi_avg - equations.c_h * (psi_B1_avg * normal_direction[1] + psi_B1_avg * normal_direction[2]) # GLM term
	+
	0.5f0 * (vk1_plus_avg * normal_direction[1] + vk2_plus_avg * normal_direction[2]) * mag_norm_avg -
	vk1_plus_avg * (B1_avg * normal_direction[1] + B2_avg * normal_direction[2]) * B1_avg -
	vk2_plus_avg * (B1_avg * normal_direction[1] + B2_avg * normal_direction[2]) * B2_avg -
	vk3_plus_avg * (B1_avg * normal_direction[1] + B2_avg * normal_direction[2]) * B3_avg # Additional terms related to the Lorentz non-conservative term (momentum eqs)
	-
	B2_avg * (vk1_minus_avg * B2_avg - vk2_minus_avg * B1_avg) * normal_direction[1] -
	B3_avg * (vk1_minus_avg * B3_avg - vk3_minus_avg * B1_avg) * normal_direction[1] -
	B1_avg * (vk2_minus_avg * B1_avg - vk1_minus_avg * B2_avg) * normal_direction[2] -
	B3_avg * (vk2_minus_avg * B3_avg - vk3_minus_avg * B2_avg) * normal_direction[2]) # Terms related to the multi-ion non-conservative term (induction equation!)

---

[JuliaFormatter] _{reported by reviewdog 🐶}

Trixi.jl/src/equations/ideal_glm_mhd_multiion_2d.jl

Line 1294 in 31cc07d

v_ll = max(v_ll, abs((rho_v1 * normal_direction[1] + rho_v2 * normal_direction[2]) / rho))

---

[JuliaFormatter] _{reported by reviewdog 🐶}

Trixi.jl/src/equations/ideal_glm_mhd_multiion_2d.jl

Line 1296 in 31cc07d

v_rr = max(v_rr, abs((rho_v1 * normal_direction[1] + rho_v2 * normal_direction[2]) / rho))

---

[JuliaFormatter] _{reported by reviewdog 🐶}

Trixi.jl/src/equations/ideal_glm_mhd_multiion_2d.jl

Lines 1419 to 1421 in 31cc07d

	sqrt(0.5f0 * (a_square + b_square) +
	0.5f0 *
	sqrt((a_square + b_square)^2 - 4 * a_square * (b1 *normal_direction[1] + b2 * normal_direction[2])^2)))

[amrueda]

Thanks for adding this @SimonCan!
Besides the formatting, I think that you can improve the implementation by precomputing the inner product of vector quantities with the normal direction. Here some examples:

[amrueda]

This looks correct (up to formatting), but you can get a cleaner and better performing implementation by precomputing the inner product of the vector quantities with the normal direction. This is what is done in other equations.

I already implemented this flux using normal_direction instead of orientation for IdealGlmMhdMultiIonEquations3D. It's already in main, so you can have a look there to see how it's implemented. For instance, here is this part:

Trixi.jl/src/equations/ideal_glm_mhd_multiion_3d.jl

Lines 281 to 298 in fa9d6b1

	v_normal = v1 * normal_direction[1] + v2 * normal_direction[2] +
	v3 * normal_direction[3]
	rho_v_normal = rho * v_normal
	vk_plus_normal = vk1_plus[k] * normal_direction[1] +
	vk2_plus[k] * normal_direction[2] +
	vk3_plus[k] * normal_direction[3]
	gamma = equations.gammas[k]
	p = (gamma - 1) * (rho_e - kin_en - mag_en - div_clean_energy)
	f1 = rho_v_normal
	f2 = rho_v_normal * v1 + p * normal_direction[1]
	f3 = rho_v_normal * v2 + p * normal_direction[2]
	f4 = rho_v_normal * v3 + p * normal_direction[3]
	f5 = (kin_en + gamma * p / (gamma - 1)) * v_normal +
	2 * mag_en * vk_plus_normal -
	B_normal * (vk1_plus[k] * B1 + vk2_plus[k] * B2 + vk3_plus[k] * B3) +
	equations.c_h * psi * B_normal

It should look pretty much identical for IdealGlmMhdMultiIonEquations2D, but removing the terms that depend on normal_direction[3]

[amrueda]

Same as above. Here the 3D implementation:

Trixi.jl/src/equations/ideal_glm_mhd_multiion_3d.jl

Lines 576 to 587 in fa9d6b1

	f2 = charge_ratio_ll[k] *
	((0.5f0 * mag_norm_avg + pe_mean) * normal_direction[1] -
	B_dot_n_avg * B1_avg)
	f3 = charge_ratio_ll[k] *
	((0.5f0 * mag_norm_avg + pe_mean) * normal_direction[2] -
	B_dot_n_avg * B2_avg)
	f4 = charge_ratio_ll[k] *
	((0.5f0 * mag_norm_avg + pe_mean) * normal_direction[3] -
	B_dot_n_avg * B3_avg)
	f5 = (vk1_plus_ll[k] * normal_direction[1] +
	vk2_plus_ll[k] * normal_direction[2] +
	vk3_plus_ll[k] * normal_direction[3]) * pe_mean

[SimonCan]

Done!

[github-actions]

Remaining comments which cannot be posted as a review comment to avoid GitHub Rate Limit

JuliaFormatter

---

[JuliaFormatter] _{reported by reviewdog 🐶}

Trixi.jl/src/equations/ideal_glm_mhd_multiion_2d.jl

Line 1302 in c461f93

v_rr = max(v_rr, abs((rho_v1 * normal_direction[1] + rho_v2 * normal_direction[2]) / rho))

---

[JuliaFormatter] _{reported by reviewdog 🐶}

Trixi.jl/src/equations/ideal_glm_mhd_multiion_2d.jl

Lines 1425 to 1427 in c461f93

	sqrt(0.5f0 * (a_square + b_square) +
	0.5f0 *
	sqrt((a_square + b_square)^2 - 4 * a_square * (b1 *normal_direction[1] + b2 * normal_direction[2])^2)))

[SimonCan]

Currently the test simulation in src/callbacks_step/analysis_dg2d.jl runs into instabilities. This examples should be very similar to the tree mesh example with a comparable resolution.

[amrueda]

This would be more efficient if you computed B_dot_n_avg at the beginning of the function, as in the 3D version. This quantity is the same for all species.

See #2446 (comment).

[amrueda]

Here, again, B_dot_n_avg is computed three more times for each species.

[amrueda]

I had already suggested changing this in #2446 (comment).

This would be more efficient if you computed B_dot_n_avg at the beginning of the function, as in the 3D version. This quantity is the same for all species.

[github-actions]

[JuliaFormatter] _{reported by reviewdog 🐶}

Suggested change

	# f2 = charge_ratio_ll[k] *
	# ((0.5f0 * mag_norm_avg - B1_avg * B1_avg + pe_mean) * normal_direction[1] +
	# (-B2_avg * B1_avg) * normal_direction[2])
	# f3 = charge_ratio_ll[k] * ((-B1_avg * B2_avg) * normal_direction[1] +
	# (-B2_avg * B2_avg + 0.5f0 * mag_norm_avg + pe_mean) * normal_direction[2])
	# f4 = charge_ratio_ll[k] * (-B1_avg * B3_avg * normal_direction[1] -
	# B2_avg * B3_avg * normal_direction[2])
	# f2 = charge_ratio_ll[k] *
	# ((0.5f0 * mag_norm_avg - B1_avg * B1_avg + pe_mean) * normal_direction[1] +
	# (-B2_avg * B1_avg) * normal_direction[2])
	# f3 = charge_ratio_ll[k] * ((-B1_avg * B2_avg) * normal_direction[1] +
	# (-B2_avg * B2_avg + 0.5f0 * mag_norm_avg + pe_mean) * normal_direction[2])
	# f4 = charge_ratio_ll[k] * (-B1_avg * B3_avg * normal_direction[1] -
	# B2_avg * B3_avg * normal_direction[2])

[github-actions]

[JuliaFormatter] _{reported by reviewdog 🐶}

Suggested change

	# f5 += dot(B_ll, cross(vk_minus_avg, B_avg)) ⋅ normal_direction
	f5 += dot(B_ll, cross(vk_minus_avg, B_avg)) * dot(normal_direction, normal_direction)
	# f5 += (B2_ll * (vk1_minus_avg * B2_avg - vk2_minus_avg * B1_avg) +
	# B3_ll * (vk1_minus_avg * B3_avg - vk3_minus_avg * B1_avg)) *
	# normal_direction[1]
	# f5 += (B1_ll * (vk2_minus_avg * B1_avg - vk1_minus_avg * B2_avg) +
	# B3_ll * (vk2_minus_avg * B3_avg - vk3_minus_avg * B2_avg)) *
	# normal_direction[2]
	# f5 += dot(B_ll, cross(vk_minus_avg, B_avg)) ⋅ normal_direction
	f5 += dot(B_ll, cross(vk_minus_avg, B_avg)) *
	dot(normal_direction, normal_direction)
	# f5 += (B2_ll * (vk1_minus_avg * B2_avg - vk2_minus_avg * B1_avg) +
	# B3_ll * (vk1_minus_avg * B3_avg - vk3_minus_avg * B1_avg)) *
	# normal_direction[1]
	# f5 += (B1_ll * (vk2_minus_avg * B1_avg - vk1_minus_avg * B2_avg) +
	# B3_ll * (vk2_minus_avg * B3_avg - vk3_minus_avg * B2_avg)) *
	# normal_direction[2]

[github-actions]

[JuliaFormatter] _{reported by reviewdog 🐶}

Suggested change

	GP_term = charge_ratio_ll[k] * (B_ll * dot(B_avg, SVector((normal_direction..., 0))))
	GP_term = charge_ratio_ll[k] *
	(B_ll * dot(B_avg, SVector((normal_direction..., 0))))

[github-actions]

[JuliaFormatter] _{reported by reviewdog 🐶}

Suggested change

	# f2 += charge_ratio_ll[k] * B1_ll *
	# (B1_avg * normal_direction[1] + B2_avg * normal_direction[2])
	# f3 += charge_ratio_ll[k] * B2_ll *
	# (B1_avg * normal_direction[1] + B2_avg * normal_direction[2])
	# f4 += charge_ratio_ll[k] * B3_ll *
	# (B1_avg * normal_direction[1] + B2_avg * normal_direction[2])
	# f2 += charge_ratio_ll[k] * B1_ll *
	# (B1_avg * normal_direction[1] + B2_avg * normal_direction[2])
	# f3 += charge_ratio_ll[k] * B2_ll *
	# (B1_avg * normal_direction[1] + B2_avg * normal_direction[2])
	# f4 += charge_ratio_ll[k] * B3_ll *
	# (B1_avg * normal_direction[1] + B2_avg * normal_direction[2])

[github-actions]

[JuliaFormatter] _{reported by reviewdog 🐶}

Suggested change

	vk1_plus_rr[k] * mag_norm_rr)
	vk1_plus_rr[k] * mag_norm_rr)