On 21 Oct 2010, at 20:58, Bernard Kelly wrote:
Hi all.
I'm trying some GW extraction from a Carpet/McLachlan black-hole binary evolution. I thought the combination of thorns WeylScal4 and Multipole would give me something reasonable, but I'm having trouble with it.
I'm attaching a parameter file for my run. As you can see, the extraction is supposed to be over a set of coordinate spheres at radii 45, 50 55, 60 65, & 70.
What I get is NaN in the WFs at all timesteps after t=0. I'm also outputting 2D data (though not as frequently), and Psi4 seems well- behaved. So what's going wrong?
I have noted that though the routine "Multipole_Calc" should happen after "calc_np" according to its schedule.ccl, it appears -before- GROUP WeylScal4_Calculate in the run's initial printing of scheduled routines.
Could you post the schedule output? According to the copy of the code that I have from the release version, Multipole_Calc is called in ANALYSIS and WeylScal4_Calculate is called in MoL_PseudoEvolution, which is called in EVOL (and some other places). The calc_np is a relic from when both were called in ANALYSIS.
Suggestions appreciated,
I tried to run your parameter file to test, but unfortunately I have some work-in-progress hacks to my local copy of WeylScal4, so it didn't run. Nevertheless, I think I know what the problem is. The spatial triad used by WeylScal4 is a coordinate triad (d/dph, d/dr, d/ dth), suitably orthonormalised. This triad is singular on the z axis, and the way it is computed there leads to a NaN. When you compute spherical harmonic modes, points on the z axis (th = 0) are multiplied by sin(th) = 0, so the value on the axis is never used, but numerically the NaN still causes problems.
There are several approaches to this:
1. Replace the triad at that point with a different one. This way, Psi4 on th = 0 is still "Psi4", just with respect to a different tetrad at that point. This is probably the best approach, and I have it implemented in my local copy. I can commit this to the development version.
2. Use an open integration formula on (0,Pi) on the 2-sphere when doing the mode decomposition. This means that the point at th = 0 is not used. Multipole does not do this, but it wouldn't be too hard to modify it.
3. Do a numerical trick to make sure that the point at th = 0 doesn't have a NaN. The particular trick I usually use is to offset the coordinates used to construct the spatial triad by 10^-15 in one direction.
WeylScal4::offset = 1e-15
This is what I did before I implemented the triad change on the axis.
What do other people do to get around this problem? Would some more physical tetrad (e.g. quasi-Kinnersley?) help?
P.S. Is this appropriate to einsteintoolkit, or Cactus at large? I'm using the ET_2010_06 release ...
WeylScal4 and Multipole are both part of the Einstein Toolkit (in the EinsteinAnalysis arrangement), so discussion is probably best on the ET list. I have CC'd this message to that list.
Other comments about your parameter file:
* You have asked for the mode decomposition of Psi4r and Psi4i individually, and you have not set the spin weight parameter, so it will default to 0, which is not what you want. Try
Multipole::variables = "WeylScal4::Psi4r{sw=-2 cmplx='WeylScal4::Psi4i' name='Psi4'}"
instead. This tells Multipole to use spin weight -2 for Psi4, and to treat Psi4r and Psi4i as the real and imaginary parts of a complex number. This means you don't have to recombine the modes of the real and imaginary parts by hand. WeylScal4 and Multipole could be modified to work with Cactus complex grid functions, but I don't have any experience with using them, and Kranc does not currently support them.
* Your time_refinement_factors parameter is quite aggressive - you are using time refinement only for finest four levels. You can probably get away with time refinement for the last 6, or more depending on your finite differencing order and value of eta.
* I recommend using the parameters
TimerReport::out_every = 32 TimerReport::n_top_timers = 40
which will give you a quick overview of the parts of the code which are taking the most time.
Apart from those points, the parameter file looks fine.
Hi Ian.
Thanks for looking at this so carefully.
On 10/22/10 3:55 AM, Ian Hinder wrote:
Could you post the schedule output? According to the copy of the code that I have from the release version, Multipole_Calc is called in ANALYSIS and WeylScal4_Calculate is called in MoL_PseudoEvolution, which is called in EVOL (and some other places). The calc_np is a relic from when both were called in ANALYSIS.
Aha. Yes, I changed that back (i.e., put WeylScal4_Calculate back into ANALYSIS), because I wanted both in the same schedule bin. Why did that change? Are variables in MoL_PseudoEvolution available for ANALYSIS automatically (i.e. will they retain their stored values)?
I think part of my problem is that I haven't caught up with the modern set of Cactus bins -- it used to be much simpler. I'm attaching my local copy of WeylScal4's schedule.ccl, as well as the overall scheduling header for the run. The ANALYSIS-bin lines I didn't understand are 415-417 & 620-622.
I've restored the original schedule.ccl now, have recompiled, and am now rerunning. I'll let you know if this helps.
I tried to run your parameter file to test, but unfortunately I have some work-in-progress hacks to my local copy of WeylScal4, so it didn't run. Nevertheless, I think I know what the problem is. The spatial triad used by WeylScal4 is a coordinate triad (d/dph, d/dr, d/ dth), suitably orthonormalised. This triad is singular on the z axis, and the way it is computed there leads to a NaN. When you compute spherical harmonic modes, points on the z axis (th = 0) are multiplied by sin(th) = 0, so the value on the axis is never used, but numerically the NaN still causes problems.
I had a different Carpet grid structure before, and then the same extraction didn't give NaNs. I wonder why?
There are several approaches to this:
- Replace the triad at that point with a different one. This way,
Psi4 on th = 0 is still "Psi4", just with respect to a different tetrad at that point. This is probably the best approach, and I have it implemented in my local copy. I can commit this to the development version.
- Use an open integration formula on (0,Pi) on the 2-sphere when
doing the mode decomposition. This means that the point at th = 0 is not used. Multipole does not do this, but it wouldn't be too hard to modify it.
- Do a numerical trick to make sure that the point at th = 0 doesn't
have a NaN. The particular trick I usually use is to offset the coordinates used to construct the spatial triad by 10^-15 in one direction.
WeylScal4::offset = 1e-15
This is what I did before I implemented the triad change on the axis.
What do other people do to get around this problem? Would some more physical tetrad (e.g. quasi-Kinnersley?) help?
I don't think quasi-Kinnersley will be a panacea; theta = 0,pi is still a problem unless you take special precautions. Personally, I think the most "correct" thing to do is your option 2 -- use an open integration formula that avoids evaluation on the poles.
Other comments about your parameter file:
- You have asked for the mode decomposition of Psi4r and Psi4i
individually, and you have not set the spin weight parameter, so it will default to 0, which is not what you want. Try
Multipole::variables = "WeylScal4::Psi4r{sw=-2 cmplx='WeylScal4::Psi4i' name='Psi4'}"
instead. This tells Multipole to use spin weight -2 for Psi4, and to treat Psi4r and Psi4i as the real and imaginary parts of a complex number. This means you don't have to recombine the modes of the real and imaginary parts by hand. WeylScal4 and Multipole could be modified to work with Cactus complex grid functions, but I don't have any experience with using them, and Kranc does not currently support them.
Thanks; for some reason, I thought the grid function metadata had the appropriate spin-weight baked in.
Thanks for the time-refinement info also; I'll bear it in mind, but we have some (possibly misguided) reason for what I listed.
Beany
Replying to myself on this ...
I reverted the WeylScal4 schedule.ccl to its repository version didn't fix things -- still NaNs almost all the time.
I realised I wasn't accurate when I said I saw NaNs. What I actually see is that at t=0, components extracted at r = 45, 50, 55, and SOMETIMES 60 have finite values. At all later times, all modes are NaN. All r=65, 70 modes are always NaN, even at t=0.
I've tried your "offset" parameter, but that doesn't seem to help, I'm afraid (actually, it's on by default anyway). I've just enabled Multipole's "out_1d_every" option, and in fact the data has nan at the -equator- (that is, at points theta = 1.530519 & 1.611073 -- those closest to the x-y plane).
Bernard
Hi Ian.
Thanks for looking at this so carefully.
On 10/22/10 3:55 AM, Ian Hinder wrote:
Could you post the schedule output? According to the copy of the code that I have from the release version, Multipole_Calc is called in ANALYSIS and WeylScal4_Calculate is called in MoL_PseudoEvolution, which is called in EVOL (and some other places). The calc_np is a relic from when both were called in ANALYSIS.
Aha. Yes, I changed that back (i.e., put WeylScal4_Calculate back into ANALYSIS), because I wanted both in the same schedule bin. Why did that change? Are variables in MoL_PseudoEvolution available for ANALYSIS automatically (i.e. will they retain their stored values)?
I think part of my problem is that I haven't caught up with the modern set of Cactus bins -- it used to be much simpler. I'm attaching my local copy of WeylScal4's schedule.ccl, as well as the overall scheduling header for the run. The ANALYSIS-bin lines I didn't understand are 415-417& 620-622.
I've restored the original schedule.ccl now, have recompiled, and am now rerunning. I'll let you know if this helps.
I tried to run your parameter file to test, but unfortunately I have some work-in-progress hacks to my local copy of WeylScal4, so it didn't run. Nevertheless, I think I know what the problem is. The spatial triad used by WeylScal4 is a coordinate triad (d/dph, d/dr, d/ dth), suitably orthonormalised. This triad is singular on the z axis, and the way it is computed there leads to a NaN. When you compute spherical harmonic modes, points on the z axis (th = 0) are multiplied by sin(th) = 0, so the value on the axis is never used, but numerically the NaN still causes problems.
I had a different Carpet grid structure before, and then the same extraction didn't give NaNs. I wonder why?
There are several approaches to this:
- Replace the triad at that point with a different one. This way,
Psi4 on th = 0 is still "Psi4", just with respect to a different tetrad at that point. This is probably the best approach, and I have it implemented in my local copy. I can commit this to the development version.
- Use an open integration formula on (0,Pi) on the 2-sphere when
doing the mode decomposition. This means that the point at th = 0 is not used. Multipole does not do this, but it wouldn't be too hard to modify it.
- Do a numerical trick to make sure that the point at th = 0 doesn't
have a NaN. The particular trick I usually use is to offset the coordinates used to construct the spatial triad by 10^-15 in one direction.
WeylScal4::offset = 1e-15
This is what I did before I implemented the triad change on the axis.
What do other people do to get around this problem? Would some more physical tetrad (e.g. quasi-Kinnersley?) help?
I don't think quasi-Kinnersley will be a panacea; theta = 0,pi is still a problem unless you take special precautions. Personally, I think the most "correct" thing to do is your option 2 -- use an open integration formula that avoids evaluation on the poles.
Other comments about your parameter file:
- You have asked for the mode decomposition of Psi4r and Psi4i
individually, and you have not set the spin weight parameter, so it will default to 0, which is not what you want. Try
Multipole::variables = "WeylScal4::Psi4r{sw=-2 cmplx='WeylScal4::Psi4i' name='Psi4'}"
instead. This tells Multipole to use spin weight -2 for Psi4, and to treat Psi4r and Psi4i as the real and imaginary parts of a complex number. This means you don't have to recombine the modes of the real and imaginary parts by hand. WeylScal4 and Multipole could be modified to work with Cactus complex grid functions, but I don't have any experience with using them, and Kranc does not currently support them.
Thanks; for some reason, I thought the grid function metadata had the appropriate spin-weight baked in.
Thanks for the time-refinement info also; I'll bear it in mind, but we have some (possibly misguided) reason for what I listed.
Beany
On 22 Oct 2010, at 17:09, Bernard Kelly wrote:
Replying to myself on this ...
I reverted the WeylScal4 schedule.ccl to its repository version didn't fix things -- still NaNs almost all the time.
I realised I wasn't accurate when I said I saw NaNs. What I actually see is that at t=0, components extracted at r = 45, 50, 55, and SOMETIMES 60 have finite values. At all later times, all modes are NaN. All r=65, 70 modes are always NaN, even at t=0.
I've tried your "offset" parameter, but that doesn't seem to help, I'm afraid (actually, it's on by default anyway). I've just enabled Multipole's "out_1d_every" option, and in fact the data has nan at the -equator- (that is, at points theta = 1.530519 & 1.611073 -- those closest to the x-y plane).
Sounds like the symmetry boundaries are not being applied. Aha! You need to set the parameters
WeylScal4::Psi4r_group_bound = "flat" WeylScal4::Psi4i_group_bound = "flat"
Have a look at the test suite parameter files for WeylScal4. This is not necessary in the current development branch, as this is done automatically by Kranc now.
It occurs to me that we could do with an example for the use of WeylScal4 and Multipole, e.g. by using WeylScal4 in one of the McLachlan QC0 example parameter files.
Bernard
Hi Ian.
Thanks for looking at this so carefully.
On 10/22/10 3:55 AM, Ian Hinder wrote:
Could you post the schedule output? According to the copy of the code that I have from the release version, Multipole_Calc is called in ANALYSIS and WeylScal4_Calculate is called in MoL_PseudoEvolution, which is called in EVOL (and some other places). The calc_np is a relic from when both were called in ANALYSIS.
Aha. Yes, I changed that back (i.e., put WeylScal4_Calculate back into ANALYSIS), because I wanted both in the same schedule bin. Why did that change? Are variables in MoL_PseudoEvolution available for ANALYSIS automatically (i.e. will they retain their stored values)?
I think part of my problem is that I haven't caught up with the modern set of Cactus bins -- it used to be much simpler. I'm attaching my local copy of WeylScal4's schedule.ccl, as well as the overall scheduling header for the run. The ANALYSIS-bin lines I didn't understand are 415-417& 620-622.
I've restored the original schedule.ccl now, have recompiled, and am now rerunning. I'll let you know if this helps.
I tried to run your parameter file to test, but unfortunately I have some work-in-progress hacks to my local copy of WeylScal4, so it didn't run. Nevertheless, I think I know what the problem is. The spatial triad used by WeylScal4 is a coordinate triad (d/dph, d/ dr, d/ dth), suitably orthonormalised. This triad is singular on the z axis, and the way it is computed there leads to a NaN. When you compute spherical harmonic modes, points on the z axis (th = 0) are multiplied by sin(th) = 0, so the value on the axis is never used, but numerically the NaN still causes problems.
I had a different Carpet grid structure before, and then the same extraction didn't give NaNs. I wonder why?
There are several approaches to this:
- Replace the triad at that point with a different one. This way,
Psi4 on th = 0 is still "Psi4", just with respect to a different tetrad at that point. This is probably the best approach, and I have it implemented in my local copy. I can commit this to the development version.
- Use an open integration formula on (0,Pi) on the 2-sphere when
doing the mode decomposition. This means that the point at th = 0 is not used. Multipole does not do this, but it wouldn't be too hard to modify it.
- Do a numerical trick to make sure that the point at th = 0
doesn't have a NaN. The particular trick I usually use is to offset the coordinates used to construct the spatial triad by 10^-15 in one direction.
WeylScal4::offset = 1e-15
This is what I did before I implemented the triad change on the axis.
What do other people do to get around this problem? Would some more physical tetrad (e.g. quasi-Kinnersley?) help?
I don't think quasi-Kinnersley will be a panacea; theta = 0,pi is still a problem unless you take special precautions. Personally, I think the most "correct" thing to do is your option 2 -- use an open integration formula that avoids evaluation on the poles.
Other comments about your parameter file:
- You have asked for the mode decomposition of Psi4r and Psi4i
individually, and you have not set the spin weight parameter, so it will default to 0, which is not what you want. Try
Multipole::variables = "WeylScal4::Psi4r{sw=-2 cmplx='WeylScal4::Psi4i' name='Psi4'}"
instead. This tells Multipole to use spin weight -2 for Psi4, and to treat Psi4r and Psi4i as the real and imaginary parts of a complex number. This means you don't have to recombine the modes of the real and imaginary parts by hand. WeylScal4 and Multipole could be modified to work with Cactus complex grid functions, but I don't have any experience with using them, and Kranc does not currently support them.
Thanks; for some reason, I thought the grid function metadata had the appropriate spin-weight baked in.
Thanks for the time-refinement info also; I'll bear it in mind, but we have some (possibly misguided) reason for what I listed.
Beany
--
Bernard J. Kelly
NASA Goddard Space Flight Center, Code 663 8800 Greenbelt Road Greenbelt, MD 20771, U.S.A.
phone: +1 (301) 286-7243
That was it, Ian. Setting those symmetry boundaries fixed the NaNs.
Thanks again for spending time on this.
Bernard
P.S. I think a binary WF-extraction testsuite would be a good idea. If you like, I'll apply what you've shown me to the qc0_mclachlan.par file, and send in the results ...
(do I have commit privileges for EinsteinAnalysis/Multipole? I suspect I'm doing everything anonymously right now)
On 10/22/10 1:15 PM, "Ian Hinder" ian.hinder@aei.mpg.de wrote:
On 22 Oct 2010, at 17:09, Bernard Kelly wrote:
I've tried your "offset" parameter, but that doesn't seem to help, I'm afraid (actually, it's on by default anyway). I've just enabled Multipole's "out_1d_every" option, and in fact the data has nan at the -equator- (that is, at points theta = 1.530519 & 1.611073 -- those closest to the x-y plane).
Sounds like the symmetry boundaries are not being applied. Aha! You need to set the parameters
WeylScal4::Psi4r_group_bound = "flat" WeylScal4::Psi4i_group_bound = "flat"
Have a look at the test suite parameter files for WeylScal4. This is not necessary in the current development branch, as this is done automatically by Kranc now.
It occurs to me that we could do with an example for the use of WeylScal4 and Multipole, e.g. by using WeylScal4 in one of the McLachlan QC0 example parameter files.
-- Ian Hinder ian.hinder@aei.mpg.de
On 22 Oct 2010, at 20:06, Kelly, Bernard J. (GSFC-660.0)[UNIVERSITY OF MARYLAND BALTIMORE COUNTY] wrote:
That was it, Ian. Setting those symmetry boundaries fixed the NaNs.
Thanks again for spending time on this.
Bernard
P.S. I think a binary WF-extraction testsuite would be a good idea. If you like, I’ll apply what you’ve shown me to the qc0_mclachlan.par file, and send in the results ...
(do I have commit privileges for EinsteinAnalysis/Multipole? I suspect I’m doing everything anonymously right now)
There are testsuites already; what we need is a practical example. So I would just add the wave extraction parameters to the qc0- mclachlan.par parameter file in McLachlan/par (after checking with Erik).
On 10/22/10 1:15 PM, "Ian Hinder" ian.hinder@aei.mpg.de wrote:
On 22 Oct 2010, at 17:09, Bernard Kelly wrote:
I've tried your "offset" parameter, but that doesn't seem to help, I'm afraid (actually, it's on by default anyway). I've just enabled Multipole's "out_1d_every" option, and in fact the data has nan at the -equator- (that is, at points theta = 1.530519 & 1.611073 -- those closest to the x-y plane).
Sounds like the symmetry boundaries are not being applied. Aha! You need to set the parameters
WeylScal4::Psi4r_group_bound = "flat" WeylScal4::Psi4i_group_bound = "flat"Have a look at the test suite parameter files for WeylScal4. This is not necessary in the current development branch, as this is done automatically by Kranc now.
It occurs to me that we could do with an example for the use of WeylScal4 and Multipole, e.g. by using WeylScal4 in one of the McLachlan QC0 example parameter files.
-- Ian Hinder ian.hinder@aei.mpg.de
On Sat, Oct 23, 2010 at 5:28 AM, Ian Hinder ian.hinder@aei.mpg.de wrote:
On 22 Oct 2010, at 20:06, Kelly, Bernard J. (GSFC-660.0)[UNIVERSITY OF MARYLAND BALTIMORE COUNTY] wrote:
That was it, Ian. Setting those symmetry boundaries fixed the NaNs.
Thanks again for spending time on this.
Bernard
P.S. I think a binary WF-extraction testsuite would be a good idea. If you like, I’ll apply what you’ve shown me to the qc0_mclachlan.par file, and send in the results ...
(do I have commit privileges for EinsteinAnalysis/Multipole? I suspect I’m doing everything anonymously right now)
There are testsuites already; what we need is a practical example. So I would just add the wave extraction parameters to the qc0- mclachlan.par parameter file in McLachlan/par (after checking with Erik).
Yes, please do that. Wave extraction should definitely be in there.
-erik
On 22 Oct 2010, at 15:32, Bernard Kelly wrote:
Hi Ian.
Thanks for looking at this so carefully.
On 10/22/10 3:55 AM, Ian Hinder wrote:
Could you post the schedule output? According to the copy of the code that I have from the release version, Multipole_Calc is called in ANALYSIS and WeylScal4_Calculate is called in MoL_PseudoEvolution, which is called in EVOL (and some other places). The calc_np is a relic from when both were called in ANALYSIS.
Aha. Yes, I changed that back (i.e., put WeylScal4_Calculate back into ANALYSIS), because I wanted both in the same schedule bin.
Scheduling in Cactus, especially with mesh refinement, is very complicated. The schedule for WeylScal4 has been carefully chosen to work (by me, Roland and Erik) - if you change it, you're on your own :)
Why did that change?
The difference is related to prolongation of mesh refinement boundaries. Psi4 is computed using spatial derivatives, and hence needs the boundaries of the grid to be filled in via prolongation from the next coarser grid. Functions scheduled in ANALYSIS cannot have time interpolation, due to the time at which ANALYSIS is scheduled relative to the Berger-Oliger time stepping. In order that all points in the grid have sensible values, analysis quantities which take spatial derivatives should be scheduled in EVOL so that time interpolation onto their refinement boundaries is possible. It is also necessary to recompute these quantities after regridding and in several other places. To make this convenient, a new schedule group, MoL_PseudoEvolution, was introduced in EVOL in which such quantities should be computed. Variables set in this group should have three timelevels, just like evolved variables, since time interpolation requires this.
Erik: is this description accurate?
Are variables in MoL_PseudoEvolution available for ANALYSIS automatically (i.e. will they retain their stored values)?
This depends on how you schedule storage for the variables. If you schedule storage for the variables unconditionally, then yes. In fact, since output of the variables will happen after ANALYSIS, the storage had better be permanent otherwise you won't be able to output them.
I think part of my problem is that I haven't caught up with the modern set of Cactus bins -- it used to be much simpler.
We used to not have mesh refinement :)
I'm attaching my local copy of WeylScal4's schedule.ccl, as well as the overall scheduling header for the run. The ANALYSIS-bin lines I didn't understand are 415-417 & 620-622.
I have noted that though the routine "Multipole_Calc" should happen after "calc_np" according to its schedule.ccl, it appears -before- GROUP WeylScal4_Calculate in the run's initial printing of scheduled routines.
WeylScal4_Calculate is no longer scheduled "as calc_np" as it used to be, so there is no dependency information between it and Multipole. In the unmodified version, the dependency was automatically enforced because WeylScal4 was scheduled in MoL_PseudoEvolution and Multipole was scheduled in ANALYSIS.
I've restored the original schedule.ccl now, have recompiled, and am now rerunning. I'll let you know if this helps.
I tried to run your parameter file to test, but unfortunately I have some work-in-progress hacks to my local copy of WeylScal4, so it didn't run. Nevertheless, I think I know what the problem is. The spatial triad used by WeylScal4 is a coordinate triad (d/dph, d/dr, d/ dth), suitably orthonormalised. This triad is singular on the z axis, and the way it is computed there leads to a NaN. When you compute spherical harmonic modes, points on the z axis (th = 0) are multiplied by sin(th) = 0, so the value on the axis is never used, but numerically the NaN still causes problems.
I had a different Carpet grid structure before, and then the same extraction didn't give NaNs. I wonder why?
When Multipole calls the interpolator to find the values of Psi4 on the extraction spheres, Carpet uses the finest grid that the point exists on. This could mean the last points in a buffer zone, and if you haven't filled these correctly with prolongation, these will be poisoned (or undefined, if you're not using poisoning). So if your extraction spheres intersect your refinement boundaries (very difficult to avoid in 3D), you might pick up poison from them.
In my opinion, the interpolator should not give points from buffer zones, and instead use the coarser grid. I think this might be the case in the new (Mercurial) version of Carpet. Erik: is that the case?
I don't think quasi-Kinnersley will be a panacea; theta = 0,pi is still a problem unless you take special precautions. Personally, I think the most "correct" thing to do is your option 2 -- use an open integration formula that avoids evaluation on the poles.
Right, or compute sin(th) Psi4, taking the appropriate limit in the neighbourhood of the axis. But I was more interested in computing the gravitational waves themselves, not the mode decomposition. At Scri+, Psi4 is supposed to represent the waves, a physical observable, right?
Other comments about your parameter file:
- You have asked for the mode decomposition of Psi4r and Psi4i
individually, and you have not set the spin weight parameter, so it will default to 0, which is not what you want. Try
Multipole::variables = "WeylScal4::Psi4r{sw=-2 cmplx='WeylScal4::Psi4i' name='Psi4'}"
instead. This tells Multipole to use spin weight -2 for Psi4, and to treat Psi4r and Psi4i as the real and imaginary parts of a complex number. This means you don't have to recombine the modes of the real and imaginary parts by hand. WeylScal4 and Multipole could be modified to work with Cactus complex grid functions, but I don't have any experience with using them, and Kranc does not currently support them.
Thanks; for some reason, I thought the grid function metadata had the appropriate spin-weight baked in.
No - the spin weight is not stored in the grid function. But that's a good idea!
On Friday, October 22, 2010, Ian Hinder ian.hinder@aei.mpg.de wrote:
The difference is related to prolongation of mesh refinement boundaries. Psi4 is computed using spatial derivatives, and hence needs the boundaries of the grid to be filled in via prolongation from the next coarser grid. Functions scheduled in ANALYSIS cannot have time interpolation, due to the time at which ANALYSIS is scheduled relative to the Berger-Oliger time stepping. In order that all points in the grid have sensible values, analysis quantities which take spatial derivatives should be scheduled in EVOL so that time interpolation onto their refinement boundaries is possible. It is also necessary to recompute these quantities after regridding and in several other places. To make this convenient, a new schedule group, MoL_PseudoEvolution, was introduced in EVOL in which such quantities should be computed. Variables set in this group should have three timelevels, just like evolved variables, since time interpolation requires this.
Erik: is this description accurate?
Yes, this is accurate.
When Multipole calls the interpolator to find the values of Psi4 on the extraction spheres, Carpet uses the finest grid that the point exists on. This could mean the last points in a buffer zone, and if you haven't filled these correctly with prolongation, these will be poisoned (or undefined, if you're not using poisoning). So if your extraction spheres intersect your refinement boundaries (very difficult to avoid in 3D), you might pick up poison from them.
In my opinion, the interpolator should not give points from buffer zones, and instead use the coarser grid. I think this might be the case in the new (Mercurial) version of Carpet. Erik: is that the case?
This is partly correct. The centre of an interpolation stencil will not be in a buffer zone in the Mercurial version. However, the stencil itself may extend into buffer / ghost / symmetry / outer boundary zones. That means that one still needs to prolongate, synchronise and apply symmetry and boundary conditions before interpolating.
Thanks for the explanations, Ian!
-erik
users@lists.einsteintoolkit.org
-
Bernard Kelly -
Erik Schnetter -
Ian Hinder -
Kelly, Bernard J. (GSFC-660.0)[UNIVERSITY OF MARYLAND BALTIMORE COUNTY]