diff options
Diffstat (limited to 'ML_BSSN')
-rw-r--r-- | ML_BSSN/param.ccl | 6 | ||||
-rw-r--r-- | ML_BSSN/src/ML_BSSN_RHS1.c | 9 | ||||
-rw-r--r-- | ML_BSSN/src/ML_BSSN_convertFromADMBaseGamma.c | 65 | ||||
-rw-r--r-- | ML_BSSN/src/ML_BSSN_convertToADMBaseDtLapseShift.c | 10 | ||||
-rw-r--r-- | ML_BSSN/src/ML_BSSN_convertToADMBaseDtLapseShiftBoundary.c | 10 | ||||
-rw-r--r-- | ML_BSSN/src/ML_BSSN_convertToADMBaseFakeDtLapseShift.c | 10 |
6 files changed, 72 insertions, 38 deletions
diff --git a/ML_BSSN/param.ccl b/ML_BSSN/param.ccl index 8fb0583..188e0f1 100644 --- a/ML_BSSN/param.ccl +++ b/ML_BSSN/param.ccl @@ -103,6 +103,12 @@ CCTK_REAL SpatialBetaDriverRadius "Radius at which the BetaDriver starts to be r } 1000000000000 restricted: +CCTK_REAL SpatialShiftGammaCoeffRadius "Radius at which the ShiftGammaCoefficient starts to be reduced" +{ + "*:*" :: "" +} 1000000000000 + +restricted: CCTK_INT harmonicN "d/dt alpha = - f alpha^n K (harmonic=2, 1+log=1)" { *:* :: "" diff --git a/ML_BSSN/src/ML_BSSN_RHS1.c b/ML_BSSN/src/ML_BSSN_RHS1.c index 98f904a..932cdc1 100644 --- a/ML_BSSN/src/ML_BSSN_RHS1.c +++ b/ML_BSSN/src/ML_BSSN_RHS1.c @@ -486,17 +486,20 @@ void ML_BSSN_RHS1_Body(cGH const * restrict const cctkGH, int const dir, int con CCTK_REAL eta = BetaDriver*IfThen(rL > SpatialBetaDriverRadius,SpatialBetaDriverRadius*INV(rL),1); + CCTK_REAL theta = IfThen(rL > SpatialShiftGammaCoeffRadius,exp(1 - + rL*INV(SpatialShiftGammaCoeffRadius)),1); + CCTK_REAL beta1rhsL = (PDupwindNth1(beta1, i, j, k)*beta1L + PDupwindNth2(beta1, i, j, k)*beta2L + PDupwindNth3(beta1, i, j, - k)*beta3L)*ShiftAdvectionCoeff + B1L*ShiftGammaCoeff; + k)*beta3L)*ShiftAdvectionCoeff + B1L*ShiftGammaCoeff*theta; CCTK_REAL beta2rhsL = (PDupwindNth1(beta2, i, j, k)*beta1L + PDupwindNth2(beta2, i, j, k)*beta2L + PDupwindNth3(beta2, i, j, - k)*beta3L)*ShiftAdvectionCoeff + B2L*ShiftGammaCoeff; + k)*beta3L)*ShiftAdvectionCoeff + B2L*ShiftGammaCoeff*theta; CCTK_REAL beta3rhsL = (PDupwindNth1(beta3, i, j, k)*beta1L + PDupwindNth2(beta3, i, j, k)*beta2L + PDupwindNth3(beta3, i, j, - k)*beta3L)*ShiftAdvectionCoeff + B3L*ShiftGammaCoeff; + k)*beta3L)*ShiftAdvectionCoeff + B3L*ShiftGammaCoeff*theta; CCTK_REAL B1rhsL = -(B1L*eta) + ((PDupwindNth1(B1, i, j, k) - PDupwindNth1(Xt1, i, j, k))*beta1L + (PDupwindNth2(B1, i, j, k) - diff --git a/ML_BSSN/src/ML_BSSN_convertFromADMBaseGamma.c b/ML_BSSN/src/ML_BSSN_convertFromADMBaseGamma.c index 7d9e3aa..1a6159a 100644 --- a/ML_BSSN/src/ML_BSSN_convertFromADMBaseGamma.c +++ b/ML_BSSN/src/ML_BSSN_convertFromADMBaseGamma.c @@ -119,6 +119,7 @@ void ML_BSSN_convertFromADMBaseGamma_Body(cGH const * restrict const cctkGH, int CCTK_REAL gt22L = gt22[index]; CCTK_REAL gt23L = gt23[index]; CCTK_REAL gt33L = gt33[index]; + CCTK_REAL rL = r[index]; /* Include user supplied include files */ @@ -247,32 +248,44 @@ void ML_BSSN_convertFromADMBaseGamma_Body(cGH const * restrict const cctkGH, int CCTK_REAL AL = -(dtalpL*(-1 + LapseAdvectionCoeff)*INV(harmonicF)*pow(alphaL,-harmonicN)); - CCTK_REAL B1L = 6*IfThen(ShiftGammaCoeff != - 0,dtbetaxL*INV(ShiftGammaCoeff),0) + IfThen(ShiftGammaCoeff != - 0,(dtbetaxL - PDupwindNth1(betax, i, j, - k)*beta1L*ShiftAdvectionCoeff)*INV(ShiftGammaCoeff),0) + - IfThen(ShiftGammaCoeff != 0,(dtbetaxL - PDupwindNth2(betax, i, j, - k)*beta2L*ShiftAdvectionCoeff)*INV(ShiftGammaCoeff),0) + - IfThen(ShiftGammaCoeff != 0,(dtbetaxL - PDupwindNth3(betax, i, j, - k)*beta3L*ShiftAdvectionCoeff)*INV(ShiftGammaCoeff),0); - - CCTK_REAL B2L = 6*IfThen(ShiftGammaCoeff != - 0,dtbetayL*INV(ShiftGammaCoeff),0) + IfThen(ShiftGammaCoeff != - 0,(dtbetayL - PDupwindNth1(betay, i, j, - k)*beta1L*ShiftAdvectionCoeff)*INV(ShiftGammaCoeff),0) + - IfThen(ShiftGammaCoeff != 0,(dtbetayL - PDupwindNth2(betay, i, j, - k)*beta2L*ShiftAdvectionCoeff)*INV(ShiftGammaCoeff),0) + - IfThen(ShiftGammaCoeff != 0,(dtbetayL - PDupwindNth3(betay, i, j, - k)*beta3L*ShiftAdvectionCoeff)*INV(ShiftGammaCoeff),0); - - CCTK_REAL B3L = 6*IfThen(ShiftGammaCoeff != - 0,dtbetazL*INV(ShiftGammaCoeff),0) + IfThen(ShiftGammaCoeff != - 0,(dtbetazL - PDupwindNth1(betaz, i, j, - k)*beta1L*ShiftAdvectionCoeff)*INV(ShiftGammaCoeff),0) + - IfThen(ShiftGammaCoeff != 0,(dtbetazL - PDupwindNth2(betaz, i, j, - k)*beta2L*ShiftAdvectionCoeff)*INV(ShiftGammaCoeff),0) + - IfThen(ShiftGammaCoeff != 0,(dtbetazL - PDupwindNth3(betaz, i, j, - k)*beta3L*ShiftAdvectionCoeff)*INV(ShiftGammaCoeff),0); + CCTK_REAL theta = IfThen(rL > SpatialShiftGammaCoeffRadius,exp(1 - + rL*INV(SpatialShiftGammaCoeffRadius)),1); + + CCTK_REAL B1L = 6*IfThen(ShiftGammaCoeff*theta != + 0,dtbetaxL*INV(ShiftGammaCoeff)*INV(theta),0) + + IfThen(ShiftGammaCoeff*theta != 0,(dtbetaxL - PDupwindNth1(betax, i, + j, + k)*beta1L*ShiftAdvectionCoeff*theta)*INV(ShiftGammaCoeff)*INV(theta),0) + + IfThen(ShiftGammaCoeff*theta != 0,(dtbetaxL - PDupwindNth2(betax, i, + j, + k)*beta2L*ShiftAdvectionCoeff*theta)*INV(ShiftGammaCoeff)*INV(theta),0) + + IfThen(ShiftGammaCoeff*theta != 0,(dtbetaxL - PDupwindNth3(betax, i, + j, + k)*beta3L*ShiftAdvectionCoeff*theta)*INV(ShiftGammaCoeff)*INV(theta),0); + + CCTK_REAL B2L = 6*IfThen(ShiftGammaCoeff*theta != + 0,dtbetayL*INV(ShiftGammaCoeff)*INV(theta),0) + + IfThen(ShiftGammaCoeff*theta != 0,(dtbetayL - PDupwindNth1(betay, i, + j, + k)*beta1L*ShiftAdvectionCoeff*theta)*INV(ShiftGammaCoeff)*INV(theta),0) + + IfThen(ShiftGammaCoeff*theta != 0,(dtbetayL - PDupwindNth2(betay, i, + j, + k)*beta2L*ShiftAdvectionCoeff*theta)*INV(ShiftGammaCoeff)*INV(theta),0) + + IfThen(ShiftGammaCoeff*theta != 0,(dtbetayL - PDupwindNth3(betay, i, + j, + k)*beta3L*ShiftAdvectionCoeff*theta)*INV(ShiftGammaCoeff)*INV(theta),0); + + CCTK_REAL B3L = 6*IfThen(ShiftGammaCoeff*theta != + 0,dtbetazL*INV(ShiftGammaCoeff)*INV(theta),0) + + IfThen(ShiftGammaCoeff*theta != 0,(dtbetazL - PDupwindNth1(betaz, i, + j, + k)*beta1L*ShiftAdvectionCoeff*theta)*INV(ShiftGammaCoeff)*INV(theta),0) + + IfThen(ShiftGammaCoeff*theta != 0,(dtbetazL - PDupwindNth2(betaz, i, + j, + k)*beta2L*ShiftAdvectionCoeff*theta)*INV(ShiftGammaCoeff)*INV(theta),0) + + IfThen(ShiftGammaCoeff*theta != 0,(dtbetazL - PDupwindNth3(betaz, i, + j, + k)*beta3L*ShiftAdvectionCoeff*theta)*INV(ShiftGammaCoeff)*INV(theta),0); /* Copy local copies back to grid functions */ diff --git a/ML_BSSN/src/ML_BSSN_convertToADMBaseDtLapseShift.c b/ML_BSSN/src/ML_BSSN_convertToADMBaseDtLapseShift.c index 15dbd8a..48f6554 100644 --- a/ML_BSSN/src/ML_BSSN_convertToADMBaseDtLapseShift.c +++ b/ML_BSSN/src/ML_BSSN_convertToADMBaseDtLapseShift.c @@ -92,6 +92,7 @@ void ML_BSSN_convertToADMBaseDtLapseShift_Body(cGH const * restrict const cctkGH CCTK_REAL beta1L = beta1[index]; CCTK_REAL beta2L = beta2[index]; CCTK_REAL beta3L = beta3[index]; + CCTK_REAL rL = r[index]; CCTK_REAL trKL = trK[index]; /* Include user supplied include files */ @@ -105,6 +106,9 @@ void ML_BSSN_convertToADMBaseDtLapseShift_Body(cGH const * restrict const cctkGH int dir3 = Sign(beta3L); + CCTK_REAL theta = IfThen(rL > SpatialShiftGammaCoeffRadius,exp(1 - + rL*INV(SpatialShiftGammaCoeffRadius)),1); + CCTK_REAL dtalpL = (PDupwindNth1(alpha, i, j, k)*beta1L + PDupwindNth2(alpha, i, j, k)*beta2L + PDupwindNth3(alpha, i, j, k)*beta3L)*LapseAdvectionCoeff + harmonicF*(AL*(-1 + @@ -112,15 +116,15 @@ void ML_BSSN_convertToADMBaseDtLapseShift_Body(cGH const * restrict const cctkGH CCTK_REAL dtbetaxL = (PDupwindNth1(beta1, i, j, k)*beta1L + PDupwindNth2(beta1, i, j, k)*beta2L + PDupwindNth3(beta1, i, j, - k)*beta3L)*ShiftAdvectionCoeff + B1L*ShiftGammaCoeff; + k)*beta3L)*ShiftAdvectionCoeff + B1L*ShiftGammaCoeff*theta; CCTK_REAL dtbetayL = (PDupwindNth1(beta2, i, j, k)*beta1L + PDupwindNth2(beta2, i, j, k)*beta2L + PDupwindNth3(beta2, i, j, - k)*beta3L)*ShiftAdvectionCoeff + B2L*ShiftGammaCoeff; + k)*beta3L)*ShiftAdvectionCoeff + B2L*ShiftGammaCoeff*theta; CCTK_REAL dtbetazL = (PDupwindNth1(beta3, i, j, k)*beta1L + PDupwindNth2(beta3, i, j, k)*beta2L + PDupwindNth3(beta3, i, j, - k)*beta3L)*ShiftAdvectionCoeff + B3L*ShiftGammaCoeff; + k)*beta3L)*ShiftAdvectionCoeff + B3L*ShiftGammaCoeff*theta; /* Copy local copies back to grid functions */ diff --git a/ML_BSSN/src/ML_BSSN_convertToADMBaseDtLapseShiftBoundary.c b/ML_BSSN/src/ML_BSSN_convertToADMBaseDtLapseShiftBoundary.c index 44874d5..40804f6 100644 --- a/ML_BSSN/src/ML_BSSN_convertToADMBaseDtLapseShiftBoundary.c +++ b/ML_BSSN/src/ML_BSSN_convertToADMBaseDtLapseShiftBoundary.c @@ -89,6 +89,7 @@ void ML_BSSN_convertToADMBaseDtLapseShiftBoundary_Body(cGH const * restrict cons CCTK_REAL B1L = B1[index]; CCTK_REAL B2L = B2[index]; CCTK_REAL B3L = B3[index]; + CCTK_REAL rL = r[index]; CCTK_REAL trKL = trK[index]; /* Include user supplied include files */ @@ -96,14 +97,17 @@ void ML_BSSN_convertToADMBaseDtLapseShiftBoundary_Body(cGH const * restrict cons /* Precompute derivatives */ /* Calculate temporaries and grid functions */ + CCTK_REAL theta = IfThen(rL > SpatialShiftGammaCoeffRadius,exp(1 - + rL*INV(SpatialShiftGammaCoeffRadius)),1); + CCTK_REAL dtalpL = harmonicF*(AL*(-1 + LapseAdvectionCoeff) - LapseAdvectionCoeff*trKL)*pow(alphaL,harmonicN); - CCTK_REAL dtbetaxL = B1L*ShiftGammaCoeff; + CCTK_REAL dtbetaxL = B1L*ShiftGammaCoeff*theta; - CCTK_REAL dtbetayL = B2L*ShiftGammaCoeff; + CCTK_REAL dtbetayL = B2L*ShiftGammaCoeff*theta; - CCTK_REAL dtbetazL = B3L*ShiftGammaCoeff; + CCTK_REAL dtbetazL = B3L*ShiftGammaCoeff*theta; /* Copy local copies back to grid functions */ diff --git a/ML_BSSN/src/ML_BSSN_convertToADMBaseFakeDtLapseShift.c b/ML_BSSN/src/ML_BSSN_convertToADMBaseFakeDtLapseShift.c index b0f579b..828f254 100644 --- a/ML_BSSN/src/ML_BSSN_convertToADMBaseFakeDtLapseShift.c +++ b/ML_BSSN/src/ML_BSSN_convertToADMBaseFakeDtLapseShift.c @@ -89,6 +89,7 @@ void ML_BSSN_convertToADMBaseFakeDtLapseShift_Body(cGH const * restrict const cc CCTK_REAL B1L = B1[index]; CCTK_REAL B2L = B2[index]; CCTK_REAL B3L = B3[index]; + CCTK_REAL rL = r[index]; CCTK_REAL trKL = trK[index]; /* Include user supplied include files */ @@ -96,14 +97,17 @@ void ML_BSSN_convertToADMBaseFakeDtLapseShift_Body(cGH const * restrict const cc /* Precompute derivatives */ /* Calculate temporaries and grid functions */ + CCTK_REAL theta = IfThen(rL > SpatialShiftGammaCoeffRadius,exp(1 - + rL*INV(SpatialShiftGammaCoeffRadius)),1); + CCTK_REAL dtalpL = harmonicF*(AL*(-1 + LapseAdvectionCoeff) - LapseAdvectionCoeff*trKL)*pow(alphaL,harmonicN); - CCTK_REAL dtbetaxL = B1L*ShiftGammaCoeff; + CCTK_REAL dtbetaxL = B1L*ShiftGammaCoeff*theta; - CCTK_REAL dtbetayL = B2L*ShiftGammaCoeff; + CCTK_REAL dtbetayL = B2L*ShiftGammaCoeff*theta; - CCTK_REAL dtbetazL = B3L*ShiftGammaCoeff; + CCTK_REAL dtbetazL = B3L*ShiftGammaCoeff*theta; /* Copy local copies back to grid functions */ |