Texture generated from a script

Is there a way in Python to generate texture data procedurally from a script? Like, a custom mathematical formula?
In particular, I’m interested in a way of generating voxel data for a volumetric texture. Any resources about that?

Edit: OK after some research I found some possible directions I can investigate.
I created a Volume Absorbtion shader in the material node editor and attached a Voronoi pattern to it – it seems to work. Beside the density, I can also modify the colour: turns out that I can attach another procedural texture generator to the Color input. I tried with a colorful checkerboard and again, it seems to work.

Now it’s where things get messy… Because now I’d like to have my own custom function for generating the colour and density data from a quite complicated mathematical function (it’s for a scientific visualization purposes). Simply assembling the function from Math nodes won’t quite work, because I need some more advanced math to do the job, and the formula would depend on some parameters I choose. But all that stuff I can do in a Python script, as long as I can use Python’s external math libraries. So no worries about that.

The problem is, however: supposing I already have my fancy function written in Python, that can take the XYZ coordinates as the input (as well as a bunch of other parameters) and spit out the colour and density values for that point, how can I use it for the volumetric rendering in Blender? :q

After some digging through through the documentation I found this:
https://docs.blender.org/api/2.79/bpy.types.VoxelDataTexture.html
which seems to be one possible way of doing this. Usually the voxel data are being read from a binary file. I could generate a binary file, but I bet this would take a lot of space, and probably it will only be static, while I’d rather like to animate it in the animation controller by tweaking the parameters of my function. So I’d rather like to generate the voxel data inside of Blender from a Python script with my fancy function. Is it possible to use this class for that purpose?

Another problem is that none of the supported formats for the voxel data seems to support colours for the voxels :q How can I supply the colour information then? I’d like to do something of the likes of those built-in procedural textures that can generate volumetric textures somehow, but with my own custom function. Is there some sort of a “volumetric texture” supported in Blender? Is there a way to generate its colour data from a script?

Annother possibility I tried is to use an OSL shader with Cycles. I already tested a simple shader that converts the coordinates of a point into a colour and returns it on its output, which I can then connect with the Volume Absorbtion’s “Color” input and it seems to work. But I’m not sure how far I can go with this approach, since the OSL language only supports some standard math functions, but none of the more advanced “special functions”, like Bessel functions or associated Legendre polynomials, or the Gamma function, etc. – it’s not exactly Python :q So I’d rather go with Python if possible, since I can already to the required calculations there. But I need some way of feeding their results through Blender now.

Any ideas?