From b3f6c47796f3b947e2cc6b2ecb37ba4ebf926aa6 Mon Sep 17 00:00:00 2001 From: Anton Khirnov Date: Sun, 8 Apr 2018 18:02:00 +0200 Subject: ML_BSSN: implement GH lapse --- ML_BSSN/param.ccl | 21 +++++++++++++++++++++ ML_BSSN/src/ML_BSSN_lapse_evol.cc | 16 +++++++++++++--- 2 files changed, 34 insertions(+), 3 deletions(-) diff --git a/ML_BSSN/param.ccl b/ML_BSSN/param.ccl index db572db..a0fcb37 100644 --- a/ML_BSSN/param.ccl +++ b/ML_BSSN/param.ccl @@ -159,6 +159,27 @@ CCTK_INT harmonicN "d/dt alpha = - f alpha^n K (harmonic=2, 1+log=1)" *:* :: "" } 2 +restricted: +CCTK_INT lapse_gh "" +{ + *:* :: "" +} 0 + +CCTK_REAL gh_p "" +{ + *:* :: "" +} 1.0 + +CCTK_REAL gh_q "" +{ + *:* :: "" +} 0.0 + +CCTK_REAL gh_eta_l_bar "" +{ + *:* :: "" +} 0.0 + restricted: CCTK_INT ShiftAlphaPower "ShiftAlphaPower" { diff --git a/ML_BSSN/src/ML_BSSN_lapse_evol.cc b/ML_BSSN/src/ML_BSSN_lapse_evol.cc index 040cbb4..988b586 100644 --- a/ML_BSSN/src/ML_BSSN_lapse_evol.cc +++ b/ML_BSSN/src/ML_BSSN_lapse_evol.cc @@ -159,6 +159,7 @@ static void ML_BSSN_lapse_evol_Body(const cGH* restrict const cctkGH, const int CCTK_REAL_VEC beta1L CCTK_ATTRIBUTE_UNUSED = vec_load(beta1[index]); CCTK_REAL_VEC beta2L CCTK_ATTRIBUTE_UNUSED = vec_load(beta2[index]); CCTK_REAL_VEC beta3L CCTK_ATTRIBUTE_UNUSED = vec_load(beta3[index]); + CCTK_REAL_VEC phiL CCTK_ATTRIBUTE_UNUSED = vec_load(phi[index]); CCTK_REAL_VEC PDupwindNthAnti1alpha CCTK_ATTRIBUTE_UNUSED; CCTK_REAL_VEC PDupwindNthSymm1alpha CCTK_ATTRIBUTE_UNUSED; @@ -213,9 +214,18 @@ static void ML_BSSN_lapse_evol_Body(const cGH* restrict const cctkGH, const int PDupwindNthAnti2alpha = 0.0; PDupwindNthSymm2alpha = 0.0; - alpharhs[index] = -harmonicF * pow(alphaL, harmonicN) * trKL + WFactor * WL - + LapseAdvectionCoeff * (beta1L * PDupwindNthAnti1alpha + beta3L * PDupwindNthAnti3alpha + - fabs(beta1L) * PDupwindNthSymm1alpha + fabs(beta3L) * PDupwindNthSymm3alpha); + if (lapse_gh) { + const double e4phi = conformalMethod == 1 ? 1.0 / (phiL * phiL) : exp(4.0 * phiL); + const double det_gamma = pow(e4phi, 3.0); + + alpharhs[index] = -(alphaL * alphaL) * trKL + gh_eta_l_bar * pow(alphaL, gh_q) * (alphaL * alphaL) * log(pow(det_gamma, 0.5 * gh_p) / alphaL) + + LapseAdvectionCoeff * (beta1L * PDupwindNthAnti1alpha + beta3L * PDupwindNthAnti3alpha + + fabs(beta1L) * PDupwindNthSymm1alpha + fabs(beta3L) * PDupwindNthSymm3alpha); + } else { + alpharhs[index] = -harmonicF * pow(alphaL, harmonicN) * trKL + WFactor * WL + + LapseAdvectionCoeff * (beta1L * PDupwindNthAnti1alpha + beta3L * PDupwindNthAnti3alpha + + fabs(beta1L) * PDupwindNthSymm1alpha + fabs(beta3L) * PDupwindNthSymm3alpha); + } } CCTK_ENDLOOP3STR(ML_BSSN_lapse_evol); -- cgit v1.2.3