Dear all,
I am having issues with the thorn Outflow running a BHB simulation with WhiskyMHD (parfile attached). Outflow correctly computes the flows on the single orbiting BHs during the inspiral, but after the merger it fails to compute it on the newly formed BH. I get the following warning in the .out:
WARNING level 1 from host r149c18s02.marconi.cineca.it process 0 while executing schedule bin CCTK_ANALYSIS, routine Outflow::outflow in thorn Outflow, file /marconi/home/userexternal/fcattori/EinsteinToolkit/ET_2018_09/Cactus/arrangements/EinsteinAnalysis/Outflow/src/outflow.c:958: -> didn't find valid detector surface for sn=2, det=2
In my parfile I implement Outflow as
outflow::compute_every = 64 outflow::num_detectors = 3 outflow::surface_index[0] = 0 outflow::surface_index[1] = 1 outflow::surface_index[2] = 2 outflow::interpolator_name = "Lagrange polynomial interpolation" outflow::interpolator_pars = "order=4"
AHFinderDirect uses three spherical surfaces (n. 0,1,2). PunctureTracker tracks and write its values on n. 0 and 1.
Does anyone know how to address this issue?
Thanks a lot,
Federico
let me add the AHFinderDirect correctly finds an AH after merger (the one that should be stored in surface 2).
Moreover, Federico also run a sim in which AHFinderDirect was writing on surfaces 0, 1, and 2 while PunctureTracker was writing on surfaces 3 and 4. In this case Outflow was giving the same error on all three surfaces 0, 1, and 2.
My naive impression is that we are missing something and AHFinderDirect is not saving the surfaces properly and therefore Outflow is not able to read them.
I had a look at previous par files from my BNS simulations (where Outflow worked) and the only difference I saw was in the number of points in theta and phi for the surfaces.
Thanks, Bruno
Il giorno gio 20 feb 2020 alle ore 18:37 Federico Cattorini < f.cattorini@campus.unimib.it> ha scritto:
Dear all,
I am having issues with the thorn Outflow running a BHB simulation with WhiskyMHD (parfile attached). Outflow correctly computes the flows on the single orbiting BHs during the inspiral, but after the merger it fails to compute it on the newly formed BH. I get the following warning in the .out:
WARNING level 1 from host r149c18s02.marconi.cineca.it process 0 while executing schedule bin CCTK_ANALYSIS, routine Outflow::outflow in thorn Outflow, file /marconi/home/userexternal/fcattori/EinsteinToolkit/ET_2018_09/Cactus/arrangements/EinsteinAnalysis/Outflow/src/outflow.c:958: -> didn't find valid detector surface for sn=2, det=2
In my parfile I implement Outflow as
outflow::compute_every = 64 outflow::num_detectors = 3 outflow::surface_index[0] = 0 outflow::surface_index[1] = 1 outflow::surface_index[2] = 2 outflow::interpolator_name = "Lagrange polynomial interpolation" outflow::interpolator_pars = "order=4"
AHFinderDirect uses three spherical surfaces (n. 0,1,2). PunctureTracker tracks and write its values on n. 0 and 1.
Does anyone know how to address this issue?
Thanks a lot,
Federico
Users mailing list Users@einsteintoolkit.org http://lists.einsteintoolkit.org/mailman/listinfo/users
Hello Federico, Bruno,
let me add the AHFinderDirect correctly finds an AH after merger (the one that should be stored in surface 2).
Looking at the source for Outflow (the line number given by the warning):
--8<-- 957 if (sf_valid[sn]<=0) { 958 CCTK_VWarn(1, __LINE__, __FILE__, CCTK_THORNSTRING, 959 "didn't find valid detector surface for sn=%d, det=%d",(int)sn,(int)det); 960 continue; 961 } --8<--
this warning is output exactly when the surface is not marked as not "valid" which AHFinderDirect does if it cannot find a horizon or did not look for it, in driver/BH_diagnostics.cc (the lack of a "return" in line 678 looks like a bug to me):
--8<-- 674 if ((my_dont_find_after >= 0 and cctk_iteration > my_dont_find_after) or 675 (my_dont_find_after_time > my_find_after_time and cctk_time > my_dont_find_after_time)) 676 { 677 assert (! AH_data.search_flag); 678 sf_active[surface_number] = 0; 679 } 680 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; --8<--
Moreover, Federico also run a sim in which AHFinderDirect was writing on surfaces 0, 1, and 2 while PunctureTracker was writing on surfaces 3 and 4. In this case Outflow was giving the same error on all three surfaces 0, 1, and 2.
So as far as I can tell AHFinderDirect did indeed not find it or not try to find it.
My naive impression is that we are missing something and AHFinderDirect is not saving the surfaces properly and therefore Outflow is not able to read them.
Are you perhaps asking it to stop searching (the parfile does not seem to indicate this but still)? In that case valid=0 and Outflow will no longer use that surface. You could request (scalar) output for SphericalSurface::sf_valid to track the state of the variable.
I had a look at previous par files from my BNS simulations (where Outflow worked) and the only difference I saw was in the number of points in theta and phi for the surfaces.
That would not make a difference to Outflow.
Yours, Roland
Roland, thank you for your feedback. It is my understanding that AHFinderDirect was able to find the horizons since it produced the BH diagnostic files. Federico: can you double check it?
Adding SphericalSurface::sf_valid to the scalar output seems to be a very good suggestion.
Thank you very much, Bruno
Il giorno ven 21 feb 2020 alle ore 19:29 Roland Haas rhaas@illinois.edu ha scritto:
Hello Federico, Bruno,
let me add the AHFinderDirect correctly finds an AH after merger (the one that should be stored in surface 2).
Looking at the source for Outflow (the line number given by the warning):
--8<-- 957 if (sf_valid[sn]<=0) { 958 CCTK_VWarn(1, __LINE__, __FILE__, CCTK_THORNSTRING, 959 "didn't find valid detector surface for sn=%d, det=%d",(int)sn,(int)det); 960 continue; 961 } --8<--
this warning is output exactly when the surface is not marked as not "valid" which AHFinderDirect does if it cannot find a horizon or did not look for it, in driver/BH_diagnostics.cc (the lack of a "return" in line 678 looks like a bug to me):
--8<-- 674 if ((my_dont_find_after >= 0 and cctk_iteration > my_dont_find_after) or 675 (my_dont_find_after_time > my_find_after_time and cctk_time > my_dont_find_after_time)) 676 { 677 assert (! AH_data.search_flag); 678 sf_active[surface_number] = 0; 679 } 680 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; --8<--
Moreover, Federico also run a sim in which AHFinderDirect was writing on surfaces 0, 1, and 2 while PunctureTracker was writing on surfaces 3 and
In this case Outflow was giving the same error on all three surfaces 0,
1,
and 2.
So as far as I can tell AHFinderDirect did indeed not find it or not try to find it.
My naive impression is that we are missing something and AHFinderDirect
is
not saving the surfaces properly and therefore Outflow is not able to
read
them.
Are you perhaps asking it to stop searching (the parfile does not seem to indicate this but still)? In that case valid=0 and Outflow will no longer use that surface. You could request (scalar) output for SphericalSurface::sf_valid to track the state of the variable.
I had a look at previous par files from my BNS simulations (where Outflow worked) and the only difference I saw was in the number of points in
theta
and phi for the surfaces.
That would not make a difference to Outflow.
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 .
Roland, Bruno thank you very much.
You could request (scalar) output for
SphericalSurface::sf_valid to track the state of the variable.
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).
I also track the SphericalSurface::sf_active variable which is 0 up to t = 1834.67, and 1 from t_merger.
It is my understanding that AHFinderDirect was able to find the horizons
since it produced the BH diagnostic files. Federico: can you double check it?
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.
Federico
Il giorno lun 24 feb 2020 alle ore 09:58 Bruno Giacomazzo < bruno.giacomazzo@unimib.it> ha scritto:
Roland, thank you for your feedback. It is my understanding that AHFinderDirect was able to find the horizons since it produced the BH diagnostic files. Federico: can you double check it?
Adding SphericalSurface::sf_valid to the scalar output seems to be a very good suggestion.
Thank you very much, Bruno
Il giorno ven 21 feb 2020 alle ore 19:29 Roland Haas rhaas@illinois.edu ha scritto:
Hello Federico, Bruno,
let me add the AHFinderDirect correctly finds an AH after merger (the
one
that should be stored in surface 2).
Looking at the source for Outflow (the line number given by the warning):
--8<-- 957 if (sf_valid[sn]<=0) { 958 CCTK_VWarn(1, __LINE__, __FILE__, CCTK_THORNSTRING, 959 "didn't find valid detector surface for sn=%d, det=%d",(int)sn,(int)det); 960 continue; 961 } --8<--
this warning is output exactly when the surface is not marked as not "valid" which AHFinderDirect does if it cannot find a horizon or did not look for it, in driver/BH_diagnostics.cc (the lack of a "return" in line 678 looks like a bug to me):
--8<-- 674 if ((my_dont_find_after >= 0 and cctk_iteration > my_dont_find_after) or 675 (my_dont_find_after_time > my_find_after_time and cctk_time > my_dont_find_after_time)) 676 { 677 assert (! AH_data.search_flag); 678 sf_active[surface_number] = 0; 679 } 680 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; --8<--
Moreover, Federico also run a sim in which AHFinderDirect was writing on surfaces 0, 1, and 2 while PunctureTracker was writing on surfaces 3
and 4.
In this case Outflow was giving the same error on all three surfaces 0,
1,
and 2.
So as far as I can tell AHFinderDirect did indeed not find it or not try to find it.
My naive impression is that we are missing something and AHFinderDirect
is
not saving the surfaces properly and therefore Outflow is not able to
read
them.
Are you perhaps asking it to stop searching (the parfile does not seem to indicate this but still)? In that case valid=0 and Outflow will no longer use that surface. You could request (scalar) output for SphericalSurface::sf_valid to track the state of the variable.
I had a look at previous par files from my BNS simulations (where
Outflow
worked) and the only difference I saw was in the number of points in
theta
and phi for the surfaces.
That would not make a difference to Outflow.
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 .
--
Prof. Bruno Giacomazzo Department of Physics University of Milano-Bicocca Piazza della Scienza 3 20126 Milano Italy
email: bruno.giacomazzo@unimib.it phone: (+39) 02 6448 2321 web: http://www.brunogiacomazzo.org
There are only 10 types of people in the world: Those who understand binary, and those who don't
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
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 rhaas@illinois.edu 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 .
Hello Federico,
glad to hear that this worked out and that you could correct the cause of the issue.
Yours, Roland
----- Original Message ----- From: Federico Cattorini f.cattorini@campus.unimib.it Sent: 2020-03-05 - 02:50 To: rhaas@illinois.edu Subject: Re: [Users] issue with Outflow after BHB merger
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 rhaas@illinois.edu 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@lists.einsteintoolkit.org
-
Bruno Giacomazzo -
Federico Cattorini -
Roland Haas