Segfault when using 'distortion_method' = 'none'

[matthiasfabry]

Phoebe segfaults in the mesh_radiosity_problem_nbody_convex when selecting no distortion.

MWE:

import phoebe as p

lggr = p.logger()

b = p.default_binary()
b.add_dataset('lc', compute_phases=p.linspace(-0.5, 0.5, 50))
b.set_value_all('distortion_method', 'none')
b.run_compute()

stack trace:

0   libphoebe.cpython-311-darwin.so	       0x15a0ccdcb void PyArray_To3DPointVector<double>(tagPyArrayObject*, std::__1::vector<T3Dpoint<double>, std::__1::allocator<T3Dpoint<double>>>&) + 27
1   libphoebe.cpython-311-darwin.so	       0x15a088439 mesh_radiosity_problem_nbody_convex(_object*, _object*, _object*) + 1097
2   python                        	       0x1047c04e8 cfunction_call + 56
3   python                        	       0x1048b5683 _PyEval_EvalFrameDefault + 235091
4   python                        	       0x10475d898 _PyFunction_Vectorcall + 3064
5   python                        	       0x104762b4b method_vectorcall + 363
6   python                        	       0x1048bb546 _PyEval_EvalFrameDefault + 259350
7   python                        	       0x10475d898 _PyFunction_Vectorcall + 3064
8   python                        	       0x104762b4b method_vectorcall + 363
9   python                        	       0x1048bb546 _PyEval_EvalFrameDefault + 259350
10  python                        	       0x10475d898 _PyFunction_Vectorcall + 3064
11  python                        	       0x104762b4b method_vectorcall + 363
12  python                        	       0x1048bb546 _PyEval_EvalFrameDefault + 259350
13  python                        	       0x10475d898 _PyFunction_Vectorcall + 3064
14  python                        	       0x1048bb162 _PyEval_EvalFrameDefault + 258354
15  python                        	       0x10487a0d4 _PyEval_Vector + 2996
16  python                        	       0x1048794a9 PyEval_EvalCode + 249
17  python                        	       0x1048743ad builtin_exec + 381
18  python                        	       0x1048b982a _PyEval_EvalFrameDefault + 251898
19  python                        	       0x104778aea gen_send_ex2 + 202
20  python                        	       0x10488a322 _PyEval_EvalFrameDefault + 58098
21  python                        	       0x104778aea gen_send_ex2 + 202
22  python                        	       0x10488a322 _PyEval_EvalFrameDefault + 58098
23  python                        	       0x104778aea gen_send_ex2 + 202
24  python                        	       0x1047792ba gen_send + 42
25  python                        	       0x10476eaf5 method_vectorcall_O + 117
26  python                        	       0x1048b4a44 _PyEval_EvalFrameDefault + 231956
27  python                        	       0x10487a0d4 _PyEval_Vector + 2996
28  python                        	       0x1048794a9 PyEval_EvalCode + 249
29  python                        	       0x1049285b8 run_mod + 168
30  python                        	       0x1049283c5 pyrun_file + 133
31  python                        	       0x104927ea3 _PyRun_SimpleFileObject + 275
32  python                        	       0x10492783f _PyRun_AnyFileObject + 143
33  python                        	       0x10494c00b pymain_run_file_obj + 267
34  python                        	       0x10494ba85 pymain_run_file + 85
35  python                        	       0x10494b3fe Py_RunMain + 1966
36  python                        	       0x1046e96a8 main + 56
37  dyld                          	    0x7ff805d32530 start + 3056

[matthiasfabry]

maybe related:
when also disabling irradiation (as a strategy to avoid the segfault), phoebe computes nan's for the fluxes:

import phoebe as p

lggr = p.logger()

b = p.default_binary()
b.add_dataset('lc', compute_phases=p.linspace(-0.5, 0.5, 50))
b.set_value_all('distortion_method', 'none')
b.set_value_all('irrad_method', 'none')

b.run_compute()

b['value@fluxes@model']

outputs:

array([nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan,
       nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan,
       nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan,
       nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan])

[kecnry]

What is the use case for disabling the meshes for both stars? Does this work if you only set to none for the secondary as in the example in the docs?

https://phoebe-project.org/docs/2.4/examples/distortion_method_none

We could add a check to skip computations entirely and just return 0 for fluxes (still allowing orbits) or we could just include this in compute checks and prevent calling run compute at all.

[matthiasfabry]

I guess this was me not understanding the options of deformation_method, this is not documented in the physics section BTW. I guess sphere is what I would consider "no distortion". None apparently skips building a mesh at all. Instead I'd expect something like enable_mesh=False without a none option for the distortion_method. It's less confusing but requires some renaming and adding a ChoiceParameter. I can try and implement this if you want, including the Docs.

[kecnry]

ok, I can see how that would be confusing. In order to try to introduce the least amount of breaking changes, how about renaming distortion_method='none' to distortion_method='point-mass'?

Even with a potential rename, there still is the issue above where setting both components should either be caught by checks or handled gracefully with zeros for fluxes.

[matthiasfabry]

Good suggestion, but I still don't quite like it, a point-mass is not something that can be distorted.

The checks I fully agree with.
I'll try to take care of these things:

-[ ] make a new BoolParameter enable_mesh and remove distortion_method='none' in favor of enable_mesh=False

• catch when running compute without meshes
• write a section in the physics docs on mesh_shapes

[kecnry]

From what I've seen from questions in the past, I worry that an additional parameter just adds too much complexity for this edge-case though and lowers confusion for this case in exchange for confusion in the much more common use-case. Let's discuss at a telecon and come to a consensus on a plan forward. Either way, we also need to write migration logic to maintain backwards compatibility.

[kecnry]

If rewriting this now, I'd probably do something like "mesh_shape" as the parameter, but I worry that that is too big of a change that isn't warranted.

[kecnry]

We can also easily modify the description on the parameter. It is currently: "Method to use for distorting stars", but could be modified to "Method to use for determining the shape of the mesh representing the star. roche: account for gravitational and rotation distortion according to roche equipotentials, rotstar: account for rotational distortion only, sphere: ignore rotation and gravitational distortion, point-source: do not generate a mesh for the star and treat as point-source mass without any light or eclipsing contributions."

[matthiasfabry]

This last one I can get behind. I'll write a section in the physics docs on this

[matthiasfabry]

Hmm. I've found in universe.py this passage:

            visibilities = meshes.get_column_flat('visibilities')

            if np.all(visibilities==0):
                # then no triangles are visible, so we should return nan -
                # probably shouldn't ever happen for lcs
                return {'flux': np.nan}

meaning nan is returned on purpose when "observing" no triangles. Do we wanna change this to zero? Or do the explicit check on the distortion_method of the meshableref. However I don't think one can access those within universe.observe. This has to be in backend._run_single_time then no?

[kecnry]

Right, this section of the code likely never gets triggered, so populating the lightcurve directly in _run_single_time, likely with nans, makes sense to me to be consistent.