From 01c6d38d3b08d32c5fad17d987ae3575dc2e6212 Mon Sep 17 00:00:00 2001 From: knarf Date: Fri, 5 Mar 2010 21:01:28 +0000 Subject: crude proper distance steering for use with TwoPunctures, working but needs to be tidied up git-svn-id: http://svn.einsteintoolkit.org/cactus/EinsteinInitialData/TOVSolver/trunk@84 1bdb13ef-5d69-4035-bb54-08abeb3aa7f1 --- param.ccl | 4 ++++ schedule.ccl | 9 +++++++++ src/tov.c | 56 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++ src/utils.inc | 5 +++++ 4 files changed, 74 insertions(+) diff --git a/param.ccl b/param.ccl index 77dc77f..2f9ff3f 100644 --- a/param.ccl +++ b/param.ccl @@ -78,6 +78,10 @@ CCTK_REAL TOV_Velocity_z[10] "(fixed) Velocity of neutron star, z coordinate (ca *:* :: "real" } 0.0 +BOOLEAN TOV_ProperPosition "For use only with two NSs, atm only handles equal mass" +{ +} "no" + BOOLEAN TOV_Fast_Interpolation "Use faster interpolation algorithm? Default is yes." { } "yes" diff --git a/schedule.ccl b/schedule.ccl index 3de374f..34fc47f 100644 --- a/schedule.ccl +++ b/schedule.ccl @@ -51,6 +51,15 @@ if (!CCTK_Equals(TOV_save_to_datafile,"")) } "Save data to file and exit" } +if (TOV_ProperPosition) +{ + schedule TOV_Set_ProperPositions in TOV_Initial_Data after TOV_C_Integrate_RHS before TOV_C_Exact + { + LANG: C + Options: GLOBAL + } "Steer NS position parameters according to proper distance" +} + if (CCTK_Equals(initial_data, "tov") || CCTK_Equals(initial_hydro, "tov") || (TOV_Use_Old_Initial_Data > 0) || TOV_Enforce_Interpolation) { diff --git a/src/tov.c b/src/tov.c index d880558..8ca0883 100644 --- a/src/tov.c +++ b/src/tov.c @@ -41,6 +41,7 @@ CCTK_REAL * TOV_Surface=0; CCTK_REAL * TOV_R_Surface=0; +CCTK_REAL * TOV_RProp_Surface=0; CCTK_REAL * TOV_Atmosphere=0; CCTK_REAL * TOV_r_1d=0; @@ -75,6 +76,11 @@ void TOV_C_ParamCheck(CCTK_ARGUMENTS) CCTK_WARN(1, "TOV_Solve_for_TOVs is depreciated. " "Use TOV_Enforce_Interpolation instead.\n"); } + if (TOV_ProperPosition) + { + if (TOV_Num_TOVs != 2) + CCTK_WARN(0, "TOV_ProperPosition atm only works for TOV_Num_TOVs==2"); + } } /* centered differencing with one-sided differencing at the boundary */ @@ -195,6 +201,7 @@ void TOV_C_Integrate_RHS(CCTK_ARGUMENTS) assert(TOV_Surface!=0); assert(TOV_R_Surface!=0); + assert(TOV_RProp_Surface!=0); assert(TOV_Atmosphere!=0); assert(TOV_r_1d!=0); @@ -331,6 +338,7 @@ void TOV_C_Integrate_RHS(CCTK_ARGUMENTS) { TOV_Surface[star] = TOV_r_1d[i]; TOV_R_Surface[star] = TOV_rbar_1d[i]; + TOV_RProp_Surface[star] = TOV_rprop_1d[i]; TOV_Surface_Index = i; } } @@ -987,4 +995,52 @@ void TOV_Prepare_Fake_Evolution(CCTK_ARGUMENTS) "Could not prepare for fake evolution - steering failed\n"); } +inline static CCTK_REAL calc_coord_dist(CCTK_REAL prop_dist, CCTK_REAL M, CCTK_REAL R) +{ + CCTK_REAL xnp = prop_dist; + CCTK_REAL xn, f, fp; + do { + xn = xnp; + f = xn - prop_dist + M * log(xn/R - 1); + fp = 1 + M / (xn-R); + xnp = xn - f/fp; + } while (fabs((xnp-xn)/xn) > 10.e-14); + return xnp; +} + +/* Only works for equal-mass binary NS systems atm */ +void TOV_Set_ProperPositions(CCTK_ARGUMENTS) +{ + DECLARE_CCTK_ARGUMENTS + DECLARE_CCTK_PARAMETERS + + /* The specified parameters are given in proper distance */ + CCTK_REAL prop_dist = sqrt((TOV_Position_x[0]-TOV_Position_x[1])* + (TOV_Position_x[0]-TOV_Position_x[1])+ + (TOV_Position_y[0]-TOV_Position_y[1])* + (TOV_Position_y[0]-TOV_Position_y[1])+ + (TOV_Position_z[0]-TOV_Position_z[1])* + (TOV_Position_z[0]-TOV_Position_z[1])); + /* Now we need to calculate the coordinate distance from that */ + CCTK_REAL coord_dist = calc_coord_dist(prop_dist, + TOV_m_1d[TOV_Num_Radial-1], + TOV_RProp_Surface[0]); + /* Now steer the parameters before the stars are interpolated onto + * the Cactus grid using those very same variables */ + /* We do that by scaling every variable by the same factor */ + CCTK_REAL factor = coord_dist / prop_dist; + char tmp[1024]; + snprintf(tmp, 1023, "%g", coord_dist/2); + if (CCTK_ParameterSet("par_b", "TwoPunctures", tmp)) + CCTK_WARN(0, "Could not set par_b"); + snprintf(tmp, 1023, "%g", TOV_Position_x[0] * factor); + if (CCTK_ParameterSet("TOV_Position_x[0]", "Whisky_TOVSolverC", tmp)) + CCTK_WARN(0, "Could not set x[0]"); + snprintf(tmp, 1023, "%g", TOV_Position_x[1] * factor); + if (CCTK_ParameterSet("TOV_Position_x[1]", "Whisky_TOVSolverC", tmp)) + CCTK_WARN(0, "Could not set x[1]"); + + // printf("proper distance: %g, coordinate distance: %g\n", prop_dist, coord_dist); +} + #include "external.inc" diff --git a/src/utils.inc b/src/utils.inc index 646554a..62bbf2a 100644 --- a/src/utils.inc +++ b/src/utils.inc @@ -6,6 +6,7 @@ void TOV_C_AllocateMemory(CCTK_ARGUMENTS) assert(TOV_Surface==0); assert(TOV_R_Surface==0); + assert(TOV_RProp_Surface==0); assert(TOV_Atmosphere==0); assert(TOV_r_1d==0); @@ -18,6 +19,7 @@ void TOV_C_AllocateMemory(CCTK_ARGUMENTS) TOV_Surface = malloc(TOV_Num_TOVs * sizeof(CCTK_REAL)); TOV_R_Surface = malloc(TOV_Num_TOVs * sizeof(CCTK_REAL)); + TOV_RProp_Surface = malloc(TOV_Num_TOVs * sizeof(CCTK_REAL)); TOV_Atmosphere = malloc(TOV_Num_TOVs * sizeof(CCTK_REAL)); TOV_r_1d = malloc(TOV_Num_Radial * TOV_Num_TOVs * sizeof(CCTK_REAL)); @@ -36,6 +38,7 @@ void TOV_C_FreeMemory (CCTK_ARGUMENTS) assert(TOV_Surface!=0); assert(TOV_R_Surface!=0); + assert(TOV_RProp_Surface!=0); assert(TOV_Atmosphere!=0); assert(TOV_r_1d!=0); @@ -46,6 +49,7 @@ void TOV_C_FreeMemory (CCTK_ARGUMENTS) free(TOV_Surface); free(TOV_R_Surface); + free(TOV_RProp_Surface); free(TOV_Atmosphere); free(TOV_r_1d); @@ -58,6 +62,7 @@ void TOV_C_FreeMemory (CCTK_ARGUMENTS) TOV_Surface=0; TOV_R_Surface=0; + TOV_RProp_Surface=0; TOV_Atmosphere=0; TOV_r_1d=0; -- cgit v1.2.3