[Commits] [svn:einsteintoolkit] TwoPunctures/trunk/src/ (Rev. 118)
knarf at cct.lsu.edu
knarf at cct.lsu.edu
Mon Jan 10 13:20:15 CST 2011
User: knarf
Date: 2011/01/10 01:20 PM
Modified:
/trunk/src/
TwoPunctures.c
Log:
patch from Bernard Kelly:
Output the individual puncture ADM mass post-solve in TwoPunctures
File Changes:
Directory: /trunk/src/
======================
File [modified]: TwoPunctures.c
Delta lines: +12 -1
===================================================================
--- trunk/src/TwoPunctures.c 2010-11-21 03:29:55 UTC (rev 117)
+++ trunk/src/TwoPunctures.c 2011-01-10 19:20:14 UTC (rev 118)
@@ -219,6 +219,7 @@
CCTK_REAL admMass;
if (! F) {
+ CCTK_REAL Mp_adm, Mm_adm, up, um;
/* Solve only when called for the first time */
F = dvector (0, ntotal - 1);
allocate_derivs (&u, ntotal);
@@ -255,7 +256,7 @@
target ADM masses target_M_plus and target_M_minus and with initial
guesses given by par_m_plus and par_m_minus. */
if(!(give_bare_mass)) {
- CCTK_REAL tmp, Mp_adm, Mm_adm, Mp_adm_err, Mm_adm_err, up, um;
+ CCTK_REAL tmp, Mp_adm_err, Mm_adm_err;
char valbuf[100];
CCTK_REAL M_p = target_M_plus;
@@ -316,6 +317,16 @@
"The two puncture masses are mp=%.17g and mm=%.17g",
(double) *mp, (double) *mm);
+ up = PunctIntPolAtArbitPosition(0, nvar, n1, n2, n3, v, par_b, 0., 0.);
+ um = PunctIntPolAtArbitPosition(0, nvar, n1, n2, n3, v,-par_b, 0., 0.);
+
+ /* Calculate the ADM masses from the current bare mass guess */
+ Mp_adm = (1 + up) * *mp + *mp * *mm / (4. * par_b);
+ Mm_adm = (1 + um) * *mm + *mp * *mm / (4. * par_b);
+
+ CCTK_VInfo (CCTK_THORNSTRING, "Puncture 1 ADM mass is %g", (double)Mp_adm);
+ CCTK_VInfo (CCTK_THORNSTRING, "Puncture 2 ADM mass is %g", (double)Mm_adm);
+
/* print out ADM mass, eq.: \Delta M_ADM=2*r*u=4*b*V for A=1,B=0,phi=0 */
admMass = (*mp + *mm
- 4*par_b*PunctEvalAtArbitPosition(v.d0, 0, 1, 0, 0, nvar, n1, n2, n3));
More information about the Commits
mailing list