Roland, thank you for your precious suggestions. I managed to solve the issue. After setting AHFinderDirect verbose level to "algorithm debug" I found out that some points of the horizon surface fell on an excised region:
AHFinderDirect::AHFinderDirect_find_horizons > in thorn AHFinderDirect, file > /marconi/home/userexternal/fcattori/EinsteinToolkit/ET_2018_09/Cactus/arrangements/EinsteinAnalysis/AHFinderDirect/src/gr/expansion.cc:949: > -> interpolate_geometry(): > one or more points on the trial horizon surface point > is/are in an excised region (or too close to the excision boundary) > Then I added to the .par file the string AHFinderDirect::mask_is_noshrink = "false" > and now everything is working properly: the surface is valid (sf_valid[2] = 1 from t_merger onwards) and outflow correctly computes the flow across it. Yours, Federico Il giorno lun 24 feb 2020 alle ore 16:20 Roland Haas <[email protected]> ha scritto: > Hello Federico, > > > I already track the SphericalSurface::sf_valid variable (line 570 of the > > parfile). The sf_valid[2] ascii yields a -1 for each timestep, even after > > t_merger = 1835.33. Yet, the other sf scalar outputs for the [2] surface > > yield real numbers after t_merger (e.g., sf_max_radius[2] yields r = > > 0.934217920021859 which is constant from the 13th iteration after the > > merger onwards, and is *-nan* before t_merger). > A value of "-1" means that (see the code snippet I had provided) that > the horizon finder tried to look for a horizon and did not find it in > that iteration. > > > As far as I can tell AHFinderDirect finds the 3rd horizon after the > merger. > > In my scalar output folder I have a BH_diagnostic.ah3.gp file which > yields > > real values starting from t_merger. > If you see "sensible" values in the other fields but sf_valid is -1 > then this would mean that the horizon finder found the horizon at least > once but is not currently finding it. Thus data on the surface if not > guaranteed to actually be a good representation of the horizon since > the "found it" could have been very long in the past. > > If you run eg > > git grep sf_max_radius > > in repos/einsteinanalysis/AHFinderDirect/src then you can see that > sf_max_radius is only accessed in one place, namely in > driver/BH_diagnostics.cc: > > --8<-- > 681 // only try to copy AH info if we've found AHs at this time level > 682 if (! AH_data.search_flag) { > 683 sf_valid[surface_number] = 0; > 684 return; > 685 } > 686 > 687 // did we actually *find* this horizon? > 688 if (! AH_data.found_flag) { > 689 sf_valid[surface_number] = -1; > 690 return; > 691 } > 692 > 693 sf_active [surface_number] = 1; > 694 sf_valid [surface_number] = 1; > 695 sf_origin_x [surface_number] = this->origin_x; > 696 sf_origin_y [surface_number] = this->origin_y; > 697 sf_origin_z [surface_number] = this->origin_z; > 698 sf_mean_radius [surface_number] = mean_radius; > 699 sf_min_radius [surface_number] = min_radius; > 700 sf_max_radius [surface_number] = max_radius; > --8<-- > > and you can see that sf_max_radius is left at its old value if sf_valid > is 0 (not trying to find) or -1 (not found). > > Thus if the horizon was ever found but is no longer then there will be > sane values in sf_XXX but sf_valid will be -1. > > If you would like Outflow to continue using the (now invalid) horizon > for a while until after it became invalid (for how long is up to you to > decide) then you will have to add some grid scalar > "horizon_last_valid_at" or so to Outflow (it must be checkpointed) and > track when you last saw a "valid" horizon. > > Yours, > Roland > > -- > My email is as private as my paper mail. I therefore support encrypting > and signing email messages. Get my PGP key from http://pgp.mit.edu . >
_______________________________________________ Users mailing list [email protected] http://lists.einsteintoolkit.org/mailman/listinfo/users
