From 8fe9d89923061e59e7141808461f01aaecf7069a Mon Sep 17 00:00:00 2001 From: Zhi Date: Mon, 23 Sep 2024 18:06:17 -0400 Subject: [PATCH] update thermo source for spherical2d geometry, this adds -p divU source term to internal energy --- Source/sources/Castro_thermo.cpp | 14 ++++++++++++-- 1 file changed, 12 insertions(+), 2 deletions(-) diff --git a/Source/sources/Castro_thermo.cpp b/Source/sources/Castro_thermo.cpp index 21ef8b9034..0b1ee6c571 100644 --- a/Source/sources/Castro_thermo.cpp +++ b/Source/sources/Castro_thermo.cpp @@ -169,8 +169,18 @@ Castro::fill_thermo_source (MultiFab& state_in, MultiFab& thermo_src) } #if AMREX_SPACEDIM >= 2 - src(i,j,k,UEINT) += -0.5_rt*(U(i,j+1,k,UMY)/U(i,j+1,k,URHO) - - U(i,j-1,k,UMY)/U(i,j-1,k,URHO))/dx[1]; + if (coord == 2) { + Real thetap = prob_lo[1] + (static_cast(j) + 1.5_rt)*dx[1]; + Real thetam = prob_lo[1] + (static_cast(j) - 0.5_rt)*dx[1]; + Real theta = 0.5_rt*(thetam + thetap); + + src(i,j,k,UEINT) += -0.5_rt*(std::sin(thetap)*U(i,j+1,k,UMY)/U(i,j+1,k,URHO) - + std::sin(thetam)*U(i,j-1,k,UMY)/U(i,j-1,k,URHO)) / + (r*std::sin(theta)*dx[1]); + } else { + src(i,j,k,UEINT) += -0.5_rt*(U(i,j+1,k,UMY)/U(i,j+1,k,URHO) - + U(i,j-1,k,UMY)/U(i,j-1,k,URHO))/dx[1]; + } #endif #if AMREX_SPACEDIM == 3 src(i,j,k,UEINT) += -0.5_rt*(U(i,j,k+1,UMZ)/U(i,j,k+1,URHO) -