claresinger / PySDM-examples

PySDM usage examples (mostly reproducing results from literature) depicting how to use PySDM in Python, in particular from Jupyter notebooks
GNU General Public License v3.0
0 stars 1 forks source link

Error with Ruehl comp film implementation #2

Closed claresinger closed 2 years ago

claresinger commented 2 years ago

Error in PySDM Ruehl surface tension implementation in the step where f_surf is found iteratively. @slayoo I'm not familiar with this numba error. Any help?

Cannot capture the non-constant value associated with variable 'A_iso' in a function that will escape.

File "../../../PySDM/PySDM/physics/surface_tension/compressed_film_Ruehl.py", line 32:
    def sigma(T, v_wet, v_dry, f_org):
        <source elided>
        # solve implicitly for fraction of organic at surface
        f = lambda f_surf: Cb_iso*(1-f_surf) - C0*np.exp(((A0**2 - (A_iso/f_surf)**2)*m_sigma*sci.N_A)/(2*sci.R*T))
        ^

https://github.com/claresinger/PySDM-examples/blob/cf_ruehl/PySDM_examples/Singer/fig_1.ipynb

slayoo commented 2 years ago

@claresinger please try using sigma.py_func() instead of sigma() to actually use the underlying Python function HTH

BTW, we've updated the arXiv e-print of the PySDM JOSS paper: https://arxiv.org/abs/2103.17238.pdf

github-actions[bot] commented 2 years ago

Stale issue message

claresinger commented 2 years ago

@slayoo I'm getting another funny error when I try to use the Ruehl surface tension. Error is that it doesn't recognize the scipy module.

I'm trying to run this notebook with the 4 different versions of surface tension to make 4 kohler curves. https://github.com/claresinger/PySDM-examples/blob/cf_ruehl/PySDM_examples/Singer/kohler.ipynb

slayoo commented 2 years ago

please try importing scipy within the function

claresinger commented 2 years ago

Thanks @slayoo. Importing within the function works, but only if I call sigma.py_func(...). What's the problem with importing outside the class like in this file? e.g. https://github.com/atmos-cloud-sim-uj/PySDM/blob/c2a5615484242b7f7a1d3560557f51fe2adb4e2d/PySDM/physics/coalescence_kernels/golovin.py#L5

Or how can I avoid needing the py_func()? I guess this turns off the numba compilation?

slayoo commented 2 years ago

In general, whatever is within the formulae is supposed to have a form of a one-liner that can be interpreted as C code on GPU. Whatever needs more than one-line logic would go to the backend API and have a separate implementation for CPU and GPU.

Having the above as a starting point, and aiming at tackling the Numba global caching problem, there is now an eval statement used within Formulae to implement constants substitution. This eval statement makes the outer imports "disappear". For use with parcel model, anyhow Numba is slowing things down rather than speeding up, and we do not use here the GPU backend, so I'd say that using py_func is OK?

slayoo commented 2 years ago

Sorry, it's exec not eval, here: https://github.com/atmos-cloud-sim-uj/PySDM/blob/master/PySDM/formulae.py#L130

claresinger commented 2 years ago

Just to be clear, does that mean the golovin kernel also will error *because scipy is not imported? **edit And numpy still works because it's explicitly included in that exec call, right?

For now, I can definitely just add the py_func, but it might be preferable to have a way of importing modules that still allows numba to be used, no?

slayoo commented 2 years ago

The kernels are not part of the Formulae "infrastructure", hence no code there is automagically njitted or converted to C. Numpy works indeed because it is explicitly included in that exec. The point is that all what is handled by the Formulae class is supposed to be njittable and best if C-convertible. Nothing besides numpy is available in Numba njitted code, same for translation to C.

slayoo commented 2 years ago

I've just noted here: https://github.com/atmos-cloud-sim-uj/PySDM/issues/726 that it would be best to move the kernels code out of "physics" not to mislead that it shares the logic with things handled from within Formulae class

slayoo commented 2 years ago

In other words (I'm aware this is puzzling):

claresinger commented 2 years ago

Is the solution to this that

  1. We write out the root solve w/o scipy, or
  2. Move the Ruehl surface tension outside of Formulae, or

Is there something else smarter/easier that can be done?

slayoo commented 2 years ago

Let's first make it work with py_func, have it test covered and merged, then we can think of refactoring.

claresinger commented 2 years ago

I opened 2 PRs to try to add tests. It's working well with py_func in my notebooks.

https://github.com/atmos-cloud-sim-uj/PySDM/pull/727 https://github.com/atmos-cloud-sim-uj/PySDM-examples/pull/118