aboutsummaryrefslogtreecommitdiff
path: root/src/gr/expansion.cc
diff options
context:
space:
mode:
Diffstat (limited to 'src/gr/expansion.cc')
-rw-r--r--src/gr/expansion.cc109
1 files changed, 70 insertions, 39 deletions
diff --git a/src/gr/expansion.cc b/src/gr/expansion.cc
index 490ff99..1946cfa 100644
--- a/src/gr/expansion.cc
+++ b/src/gr/expansion.cc
@@ -59,20 +59,25 @@ void setup_xyz_posns(patch_system& ps, bool print_msg_flag);
bool interpolate_geometry(patch_system* ps_ptr,
const struct cactus_grid_info& cgi,
const struct geometry_info& gi,
+ const struct error_info& error_info, bool initial_flag,
bool print_msg_flag);
-void convert_conformal_to_physical(patch_system& ps, bool print_msg_flag);
+void convert_conformal_to_physical(patch_system& ps,
+ bool print_msg_flag);
void Schwarzschild_EF_geometry(patch_system& ps,
const struct cactus_grid_info& cgi,
const struct geometry_info& gi,
bool print_msg_flag);
-bool h_is_finite(patch_system& ps, bool print_msg_flag);
+bool h_is_finite(patch_system& ps,
+ const struct error_info& error_info, bool initial_flag,
+ bool print_msg_flag);
bool geometry_is_finite(patch_system& ps,
- bool print_msg_flag);
+ const struct error_info& error_info, bool initial_flag,
+ bool print_msg_flag);
-bool compute_Theta(patch_system& ps,
- bool Jacobian_flag,
+bool compute_Theta(patch_system& ps, bool Jacobian_flag,
jtutil::norm<fp>* Theta_norms_ptr,
+ const struct error_info& error_info, bool initial_flag,
bool print_msg_flag);
}
@@ -122,6 +127,11 @@ bool compute_Theta(patch_system& ps,
// interpolator call are done, the latter with the number of
// interpolation points is set to 0 and all the output array
// pointers set to NULL.
+// initial_flag = true if this is the first evaluation of expansion()
+// for this horizon,
+// false otherwise;
+// this is used (only) to select which elements of error_info
+// are relevant
// Jacobian_flag = true to compute the Jacobian coefficients,
// false to skip this.
// print_msg_flag = true to print status messages,
@@ -145,6 +155,7 @@ bool compute_Theta(patch_system& ps,
bool expansion(patch_system* ps_ptr,
const struct cactus_grid_info& cgi,
const struct geometry_info& gi,
+ const struct error_info& error_info, bool initial_flag,
bool Jacobian_flag /* = false */,
bool print_msg_flag /* = false */,
jtutil::norm<fp>* Theta_norms_ptr /* = NULL */)
@@ -164,8 +175,9 @@ if (active_flag)
// fill in values of all ghosted gridfns in ghost zones
ps_ptr->synchronize();
- if (gi.check_that_h_is_finite
- && !h_is_finite(*ps_ptr, print_msg_flag))
+ if (gi.check_that_h_is_finite && !h_is_finite(*ps_ptr,
+ error_info, initial_flag,
+ print_msg_flag))
then return false; // *** ERROR RETURN ***
// set up xyz positions of grid points
@@ -181,15 +193,19 @@ case geometry__local_interp_from_Cactus_grid:
// ps_ptr (non-NULL vs NULL) to choose a normal vs dummy computation
if (!interpolate_geometry(ps_ptr,
cgi, gi,
+ error_info, initial_flag,
print_msg_flag))
then return false; // *** ERROR RETURN ***
if (active_flag && cgi.Cactus_conformal_metric)
- then convert_conformal_to_physical(*ps_ptr, print_msg_flag);
+ then convert_conformal_to_physical(*ps_ptr,
+ print_msg_flag);
break;
case geometry__Schwarzschild_EF:
if (active_flag)
- then Schwarzschild_EF_geometry(*ps_ptr, gi, print_msg_flag);
+ then Schwarzschild_EF_geometry(*ps_ptr,
+ gi,
+ print_msg_flag);
break;
default:
@@ -203,7 +219,9 @@ default:
if (active_flag)
then {
if (gi.check_that_geometry_is_finite
- && !geometry_is_finite(*ps_ptr, print_msg_flag))
+ && !geometry_is_finite(*ps_ptr,
+ error_info, initial_flag,
+ print_msg_flag))
then return false; // *** ERROR RETURN ***
// compute remaining gridfns --> $\Theta$
@@ -211,6 +229,7 @@ if (active_flag)
// by algebraic ops and angular finite differencing
if (!compute_Theta(*ps_ptr,
Jacobian_flag, Theta_norms_ptr,
+ error_info, initial_flag,
print_msg_flag))
then return false; // *** ERROR RETURN ***
}
@@ -322,6 +341,7 @@ namespace {
bool interpolate_geometry(patch_system* ps_ptr,
const struct cactus_grid_info& cgi,
const struct geometry_info& gi,
+ const struct error_info& error_info, bool initial_flag,
bool print_msg_flag)
{
const bool active_flag = (ps_ptr != NULL);
@@ -651,19 +671,26 @@ if (status == CCTK_ERROR_INTERP_POINT_OUTSIDE)
const char axis = "xyz"[error_axis];
assert((error_direction == -1) || (error_direction == +1));
- const char end = (error_direction == -1) ? '-' : '+';
+ const char direction = (error_direction == -1) ? '-' : '+';
if (print_msg_flag)
then {
- CCTK_VWarn(expansion_warning_levels::failure,
- __LINE__, __FILE__, CCTK_THORNSTRING,
+ const int warn_level
+ = initial_flag
+ ? error_info.warn_level__point_outside__initial
+ : error_info.warn_level__point_outside__subsequent;
+ CCTK_VWarn(warn_level, __LINE__, __FILE__, CCTK_THORNSTRING,
"\n"
-" the trial-horizon-surface point"
-" global_(x,y,z)=(%g,%g,%g)"
-" is outside the grid (or too close to the boundary) in the %c%c direction!"
+"interpolate_geometry():\n"
+" the trial-horizon-surface point\n"
+" global_(x,y,z)=(%g,%g,%g)\n"
+" is outside the grid (or too close to the boundary)"
+ " in the %c%c direction!\n"
+" (this may mean that driver::ghost_size is too small,\n"
+" especially if any horizon crosses a symmetry boundary)"
,
global_x, global_y, global_z,
- end, axis);
+ direction, axis);
}
return false; // *** ERROR RETURN ***
}
@@ -811,7 +838,7 @@ if (print_msg_flag)
// or more NaNs or infinities.
//
// Results:
-#ifdef Theta_AVE_FINITE
+#ifdef HAVE_FINITE
// This function returns true if all the h values are finite, false
// otherwise (i.e. if it contains any NaNs or infinities).
#else
@@ -820,12 +847,14 @@ if (print_msg_flag)
#endif
//
namespace {
-bool h_is_finite(patch_system& ps, bool print_msg_flag)
+bool h_is_finite(patch_system& ps,
+ const struct error_info& error_info, bool initial_flag,
+ bool print_msg_flag)
{
if (print_msg_flag)
then CCTK_VInfo(CCTK_THORNSTRING, " checking that h is finite");
-#ifdef Theta_AVE_FINITE
+#ifdef HAVE_FINITE
for (int pn = 0 ; pn < ps.N_patches() ; ++pn)
{
patch& p = ps.ith_patch(pn);
@@ -843,7 +872,7 @@ if (print_msg_flag)
const fp sigma = p.sigma_of_isigma(isigma);
const fp drho = jtutil::degrees_of_radians(rho);
const fp dsigma = jtutil::degrees_of_radians(sigma);
- CCTK_VWarn(expansion_warning_levels::nonfinite,
+ CCTK_VWarn(error_info.warn_level__nonfinite_geometry,
__LINE__, __FILE__, CCTK_THORNSTRING,
"\n"
" h=%g isn't finite!\n"
@@ -859,7 +888,7 @@ if (print_msg_flag)
}
return true; // *** all values finite ***
#else
-CCTK_VWarn(expansion_warning_levels::skip_nonfinite,
+CCTK_VWarn(error_info.warn_level__skipping_finite_check,
__LINE__, __FILE__, CCTK_THORNSTRING,
" no finite() fn ==> skipping is-h-finite check");
return true; // *** no check possible ***
@@ -878,7 +907,7 @@ return true; // *** no check possible ***
// or more NaNs or infinities.
//
// Results:
-#ifdef Theta_AVE_FINITE
+#ifdef HAVE_FINITE
// This function returns true if all the geometry variables are finite,
// false otherwise (i.e. if they contain any NaNs or infinities).
#else
@@ -887,12 +916,14 @@ return true; // *** no check possible ***
#endif
//
namespace {
-bool geometry_is_finite(patch_system& ps, bool print_msg_flag)
+bool geometry_is_finite(patch_system& ps,
+ const struct error_info& error_info, bool initial_flag,
+ bool print_msg_flag)
{
if (print_msg_flag)
then CCTK_VInfo(CCTK_THORNSTRING, " checking that geometry is finite");
-#ifdef Theta_AVE_FINITE
+#ifdef HAVE_FINITE
for (int pn = 0 ; pn < ps.N_patches() ; ++pn)
{
patch& p = ps.ith_patch(pn);
@@ -989,7 +1020,7 @@ if (print_msg_flag)
const fp global_x = ps.global_x_of_local_x(local_x);
const fp global_y = ps.global_y_of_local_y(local_y);
const fp global_z = ps.global_z_of_local_z(local_z);
- CCTK_VWarn(expansion_warning_levels::nonfinite,
+ CCTK_VWarn(error_info.warn_level__nonfinite_geometry,
__LINE__, __FILE__, CCTK_THORNSTRING,
"\n"
" geometry isn't finite at %s patch\n"
@@ -1041,7 +1072,7 @@ if (print_msg_flag)
}
return true; // *** no NaNs found ***
#else
-CCTK_VWarn(expansion_warning_levels::skip_nonfinite,
+CCTK_VWarn(error_info.warn_level__skipping_finite_check,
__LINE__, __FILE__, CCTK_THORNSTRING,
" no finite() ==> skipping is-geometry-finite check");
return true; // *** no check possible ***
@@ -1073,9 +1104,9 @@ return true; // *** no check possible ***
// g_ij isn't positive definite).
//
namespace {
-bool compute_Theta(patch_system& ps,
- bool Jacobian_flag,
+bool compute_Theta(patch_system& ps, bool Jacobian_flag,
jtutil::norm<fp>* Theta_norms_ptr,
+ const struct error_info& error_info, bool initial_flag,
bool print_msg_flag)
{
if (print_msg_flag)
@@ -1159,16 +1190,16 @@ if (print_msg_flag)
if (Theta_D <= 0)
then {
- CCTK_VWarn(expansion_warning_levels::failure,
- __LINE__, __FILE__, CCTK_THORNSTRING,
- "Theta_D <= 0 at %s patch rho=%g sigma=%g!",
- p.name(), double(rho), double(sigma));
- CCTK_VWarn(expansion_warning_levels::failure,
- __LINE__, __FILE__, CCTK_THORNSTRING,
- "(this probably means the interpolated g_ij");
- CCTK_VWarn(expansion_warning_levels::failure,
- __LINE__, __FILE__, CCTK_THORNSTRING,
- " isn't positive definite)");
+const int warn_level
+ = initial_flag ? error_info.warn_level__gij_not_positive_definite__initial
+ : error_info.warn_level__gij_not_positive_definite__subsequent;
+CCTK_VWarn(warn_level, __LINE__, __FILE__, CCTK_THORNSTRING,
+"\n"
+" compute_Theta(): Theta_D = $g^{ij} s_i s_j$ = %g <= 0\n"
+" at %s patch rho=%g sigma=%g!\n"
+" (i.e. the interpolated g_ij isn't positive definite)",
+ double(Theta_D),
+ p.name(), double(rho), double(sigma));
return false; // *** ERROR RETURN ***
}