-
Notifications
You must be signed in to change notification settings - Fork 55
add tdc uq into HLLC #883
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
base: main
Are you sure you want to change the base?
add tdc uq into HLLC #883
Conversation
Co-authored-by: Warrick Ball <[email protected]>
star/private/hydro_energy.f90
Outdated
| if (k < s% nz) then | ||
| TDC_eturb_cell_start = 0.5d0*(pow2(s% mlt_vc_old(k)/sqrt_2_div_3) + & | ||
| pow2(s% mlt_vc_old(k+1)/sqrt_2_div_3)) | ||
| TDC_eturb_cell = 0.5d0*(pow2(s% mlt_vc_ad(k)/sqrt_2_div_3) + & | ||
| pow2(shift_p1(s% mlt_vc_ad(k+1))/sqrt_2_div_3)) | ||
| else ! center cell averaged with 0 for inner face | ||
| TDC_eturb_cell_start = 0.5d0*pow2(s% mlt_vc_old(k)/sqrt_2_div_3) | ||
| TDC_eturb_cell = 0.5d0*pow2(s% mlt_vc(k)/sqrt_2_div_3) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I feel like there's some simple factorisation of the constants available here. We've got (√(3/2))² from the pow2 expressions and ½ out front, so is this not the same as
| if (k < s% nz) then | |
| TDC_eturb_cell_start = 0.5d0*(pow2(s% mlt_vc_old(k)/sqrt_2_div_3) + & | |
| pow2(s% mlt_vc_old(k+1)/sqrt_2_div_3)) | |
| TDC_eturb_cell = 0.5d0*(pow2(s% mlt_vc_ad(k)/sqrt_2_div_3) + & | |
| pow2(shift_p1(s% mlt_vc_ad(k+1))/sqrt_2_div_3)) | |
| else ! center cell averaged with 0 for inner face | |
| TDC_eturb_cell_start = 0.5d0*pow2(s% mlt_vc_old(k)/sqrt_2_div_3) | |
| TDC_eturb_cell = 0.5d0*pow2(s% mlt_vc(k)/sqrt_2_div_3) | |
| if (k < s% nz) then | |
| TDC_eturb_cell_start = 0.75d0*(pow2(s% mlt_vc_old(k)) + & | |
| pow2(s% mlt_vc_old(k+1))) | |
| TDC_eturb_cell = 0.75d0*(pow2(s% mlt_vc_ad(k)) + & | |
| pow2(shift_p1(s% mlt_vc_ad(k+1)))) | |
| else ! center cell averaged with 0 for inner face | |
| TDC_eturb_cell_start = 0.75d0*pow2(s% mlt_vc_old(k)) | |
| TDC_eturb_cell = 0.75d0*pow2(s% mlt_vc(k)) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Ahh yes, I didn't even consider reducing this. Good catch. (3/2) * (1/2) -> 3/4
| integer :: k, op_err | ||
| include 'formats' | ||
| ierr = 0 | ||
| op_err = 0 |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Do we actually use op_err? It looks like both it and ierr and initialised as zero, then ierr is changed to match op_err whenever the latter is non-zero. But that seems the same as just using ierr instead of op_err.
Correctly includes tdc Uq into HLLC as a non-HSE source in the du/dt eqn, and fixes other bugs. Needs further polish before merge.