From f6290b709fba33db2fc3b4b3e2a73cc7ec9e8432 Mon Sep 17 00:00:00 2001 From: jthorn Date: Fri, 8 Jul 2005 10:41:34 +0000 Subject: Add an option to control whether we generate a conformal metric (StaticConformal::conformal_state = 3) or a physical metric ( = 0), as per http://www.cactuscode.org/pipermail/developers/2005-July/001178.html and subsequent discussion git-svn-id: http://svn.einsteintoolkit.org/cactus/EinsteinInitialData/IDAxiBrillBH/trunk@79 0a4070d5-58f5-498f-b6c0-2693e757fa0f --- doc/documentation.tex | 30 ++++++++++++++++++++++++++---- interface.ccl | 15 +++++++++++++++ param.ccl | 12 ++++++++++++ src/IDAxiBrillBH.F | 25 +++++++++++++++++++++---- 4 files changed, 74 insertions(+), 8 deletions(-) diff --git a/doc/documentation.tex b/doc/documentation.tex index 0331d01..c14f6a6 100644 --- a/doc/documentation.tex +++ b/doc/documentation.tex @@ -17,14 +17,21 @@ % Do not delete next line % START CACTUS THORNGUIDE +% Add all definitions used in this documentation here +% \def\mydef etc + +\def\thorn#1{\textbf{#1}} +\def\arrangement#1{\textbf{#1}} + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% \begin{abstract} - Thorn IDAxiBrillBH provides analytic initial data for a vacuum + Thorn \thorn{IDAxiBrillBH} provides analytic initial data for a vacuum black hole spacetime: a single Schwarzschild black hole in isotropic coordinates plus Brill wave. This initial data is - provided for the 3-conformal metric, it's spatial derivatives, and - extrinsic curvature. + provided for \thorn{ADMBase} 3-metric and extrinsic curvature, + and optionally also for the \thorn{StaticConformal} conformal + factor and its 1st and 2nd~spatial derivatives. \end{abstract} %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% @@ -216,7 +223,7 @@ together with either 4th~order Lagrange or 3rd~order Hermite interpolation (provided by thorn \textbf{AEIThorns/AEILocalInterp}) to get sufficient accuracy. -One problem with such high resolutions is that \textbf{IDAxiBrillBH} +One problem with such high resolutions is that \thorn{IDAxiBrillBH} uses an internal multigrid solver which allocates local arrays on the stack, whose size depends on the $\eta$ and $\theta$ resolutions. For high resolutions these arrays may exceed system- and/or shell-imposed @@ -233,6 +240,21 @@ and \verb|ulimit -s unlimited|. %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +\section{Physical or Conformal Metric} + +By default, \thorn{IDAxiBrillBH} generates initial data which uses +a nontrivial static conformal factor (as defined by thorn +\thorn{StaticConformal}). This initial data includes both the +conformal factor and its 1st and 2nd~spatial derivatives, +so \thorn{IDAxiBrillBH} sets \verb|conformal_state| to 3. + +However, if the Boolean parameter \verb|generate_StaticConformal_metric| +is set to \verb|false|, then \thorn{IDAxiBrillBH} generates a pure +physical 3-metric (and sets \verb|conformal_state| to 0). This is +useful if you have other thorns which don't grok a conformal metric. + +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + \section{Debugging Parameters} This thorn has options to print very detailed debugging information diff --git a/interface.ccl b/interface.ccl index 1ed1806..75ff091 100644 --- a/interface.ccl +++ b/interface.ccl @@ -5,3 +5,18 @@ implements: idaxibrillbh inherits: ADMBase grid StaticConformal public: + +# +# functions from StaticConformal +# +void FUNCTION ConfToPhysInPlace (CCTK_INT IN nx, \ + CCTK_INT IN ny, \ + CCTK_INT IN nz, \ + CCTK_REAL ARRAY IN psi, \ + CCTK_REAL ARRAY INOUT gxx, \ + CCTK_REAL ARRAY INOUT gxy, \ + CCTK_REAL ARRAY INOUT gxz, \ + CCTK_REAL ARRAY INOUT gyy, \ + CCTK_REAL ARRAY INOUT gyz, \ + CCTK_REAL ARRAY INOUT gzz) +uses function ConfToPhysInPlace diff --git a/param.ccl b/param.ccl index b20a0df..db7206b 100644 --- a/param.ccl +++ b/param.ccl @@ -167,3 +167,15 @@ INT interpolation_order "Order for interpolation" STEERABLE = ALWAYS { 0:9 :: "any integer accepted by the interpolator" } 1 + +############################################################ + +# +# ***** miscellaneous parameters ***** +# + +Boolean generate_StaticConformal_metric \ + "should we generate a StaticConformal conformal metric (true), \ + or a pure ADMBase physical metric (false)" +{ +} true diff --git a/src/IDAxiBrillBH.F b/src/IDAxiBrillBH.F index f6553cb..93f6db8 100644 --- a/src/IDAxiBrillBH.F +++ b/src/IDAxiBrillBH.F @@ -521,6 +521,21 @@ c dxgxx = 1/2 \partial gxx / \partial x enddo enddo +c convert to physical metric if StaticConformal isn't wanted + if (generate_StaticConformal_metric .eq. 0) then + call CCTK_INFO("converting to physical metric") + call ConfToPhysInPlace(nx, ny, nz, + $ psi, + $ gxx, gxy, gxz, + $ gyy, gyz, + $ gzz) +c record that we now have a physical metric + conformal_state = 0 + else +c record that we computed psi and its 1st and 2nd derivatives + conformal_state = 3 + end if + c Extrinsic Curvature is identically zero kxx = 0.0d0 kxy = 0.0d0 @@ -531,7 +546,12 @@ c Extrinsic Curvature is identically zero if (debug .ge. 6) then print *, '### final results (again at this 3-D grid point) ###' - print *, ' psi =', psi(debug_i,debug_j,debug_k) + if (conformal_state .gt. 0) then + print *, '### ... conformal metric with' + print *, ' psi =', psi(debug_i,debug_j,debug_k) + else + print *, '### ... physical metric with' + end if print *, ' gxx =', gxx(debug_i,debug_j,debug_k) print *, ' gxy =', gxy(debug_i,debug_j,debug_k) print *, ' gxz =', gxz(debug_i,debug_j,debug_k) @@ -562,9 +582,6 @@ c ADM mass alp = (2.0d0*r - adm)/(2.0d0*r+adm) endif -c conformal_state = 'all' derivatives were calculated - conformal_state = 3 - deallocate(cc,ce,cw,cn,cs,rhs,psi2d,detapsi2d,dqpsi2d, $ detaetapsi2d,detaqpsi2d,dqqpsi2d, $ etagrd,qgrd, -- cgit v1.2.3