My system
MacOS Version 10.15.7 (19H1323)
Blender Version 2.93.4 (2.93.4 2021-09-01)
Background
I am working on an idea I have for a scientific illustration within the field of theoretical chemistry. I will not go into the details, but give a short overview for some context. Figure 1 shows two potential energy surfaces that mostly are separated by a gap except for in a single point where they meet and eventually intersect. Interesting stuff regarding chemical dynamics take place near so-called conical intersections.
Figure 1. Conical intersection of two potential energy surfaces.
The energies of the two surfaces close to the intersection point can be approximated by a relatively simple model (in polar coordinates):
$$ E_{\mathrm{lower}}(r,\theta) = E^\times + r\left[ \sigma \cos(\theta_s - \theta) - \sqrt{1 + \Delta \cos 2\theta} \right] \\ E_{\mathrm{upper}}(r,\theta) = E^\times + r\left[ \sigma \cos(\theta_s - \theta) + \sqrt{1 + \Delta \cos 2\theta} \right] $$
with
$$ \sigma \in [] \\ \theta_s \in [0, 2\pi] \\ \Delta \in [0, 1] $$
where
$r$ and $\theta$ are the radial and angular coordinates, respectively. The model contains some parameters that define the topography of the conical intersection: $\sigma$ describes the amount of tilt, $\theta_s$ defines the direction of tilt, and $\Delta$ is an asymmetry parameter. The functions above can be plotted in order to visualize the conical intersections.
Python code for plotting:
import matplotlib.pyplot as plt
import numpy as np
class ConicalIntersection:
def __init__(self, size=1, nfaces=5, origo=0, asym=0, tilt=0, angle=0):
self.nfaces = nfaces
self.nverts = nfaces + 1
self.asym = asym
self.tilt = tilt
self.angle = angle
self.o = origo
# Set up mesh
_t = np.linspace(0, 2*np.pi, self.nverts)
_r = np.linspace(0, size, 2)
_r, _t = np.meshgrid(_r, _t)
self.x = (_r * np.sin(_t))
self.y = (_r * np.cos(_t))
self.zlower = (self.o +_r * (self.tilt * np.cos(_t - self.angle) - np.sqrt(1 + self.asym * np.cos(2*_t))))
self.zupper = (self.o + _r * (self.tilt * np.cos(_t - self.angle) + np.sqrt(1 + self.asym * np.cos(2*_t))))
ci = ConicalIntersection(nfaces=50, asym=0.6, angle=np.pi/3, tilt=0.5)
fig, ax = plt.subplots(dpi=200, subplot_kw={'projection': '3d'})
ax.set_axis_off()
ax.elev = 25
ax.azim = 25
ax.plot_surface(ci.x, ci.y, ci.zlower, color='crimson', edgecolor='black', lw=0.2)
ax.plot_surface(ci.x, ci.y, ci.zupper, color='teal', edgecolor='black', lw=0.2)
plt.savefig('ConicalIntersection.png')
An example of what this model looks like (slightly tilted and asymmetric):
My goal
I would like to create in Blender something that resembles the top-most figure of the crossings surfaces, using the models above to build the areas closest to the intersections. In addition, I envision having about 3-4 such conical intersections in my energy landscape, where all with different topographies as obtained by using different parameters to generate the data. The transition from the conical intersection should be smooth.
My attempted solution
I have adapted the code from above to add meshes to Blender (code below). My idea was to add a plane and apply and tweak an ocean modifier to get a smooth surface with a few hills and valleys, duplicate it and place the second one above the first one and slightly shifted in either x or y. I would then add a couple of conical intersections and place them where I wanted them.
However, I am not sure how to smoothly connect the conical intersections to the two surfaces, and would love some input on how to achieve this. The code below should work as is if you paste it into the script editor and run it, if you would like to play around with these conical intersections. It would be awesome to make a plugin where you can tweak the parameters in real-time when making these conical intersections, but that is beyond my Blender API capabilities.
Part of the difficulties is to match the topography of the surfaces to the topography of the conical intersections. This will make sense if you look at my preliminary render at the bottom. Should I sculpt the surfaces myself, or rely on some automatic mesh modification? Also, if I increase the number of faces on the cones (to make them smoother), I imagine it will be more work to actually get a smooth connection.
I would greatly appreciate some tips and guidance! If anything is confusing or unclear, let me know and I will try to make it clearer.
import os
import bpy
import numpy as np
import itertools as it
NFACES = 5
SIZE = 5
TILT = 0.6
ASYM = 0.4
ANGLE = 0
class ConicalIntersection:
"""
Class holding parameters defining the CI,
plus the vertices and faces that for showing
it in Blender.
"""
def __init__(self, size=1, nfaces=5, origo=0, asym=0, tilt=0, angle=0):
"""
Parameters
size : int, radius of cone
nfaces : int, number of faces around the cone (resolution)
origo : float, z-coordinate at the intersection point
asym : float, asymmetry parameter (\Delta)
tilt : float, how much tilt (\sigma)
angle : angle of tilt in radians (\theta_s)
"""
self.nfaces = nfaces
self.nverts = nfaces + 1
self.asym = asym
self.tilt = tilt
self.angle = angle
self.o = origo
# Set up mesh
_t = np.linspace(0, 2*np.pi, self.nverts)
_r = np.linspace(0, size, 2)
_r, _t = np.meshgrid(_r, _t)
# Get coordinates in Cartesian space
self.x = (_r * np.sin(_t))[:,-1]
self.y = (_r * np.cos(_t))[:,-1]
self.zlower = (self.o +_r * (self.tilt * np.cos(_t - self.angle) - np.sqrt(1 + self.asym * np.cos(2*_t))))[:,-1]
self.zupper = (self.o + _r * (self.tilt * np.cos(_t - self.angle) + np.sqrt(1 + self.asym * np.cos(2*_t))))[:,-1]
# Compute the vertices and faces
self.verts_lower = [(0.0, 0.0, 0.0)] + list(zip(self.x, self.y, self.zlower))
self.verts_upper = [(0.0, 0.0, 0.0)] + list(zip(self.x, self.y, self.zupper))
self.faces = self.filter_faces()
def filter_faces(self):
"""
I generated the faces by all possible combinations of
the vertices. But we just want to have the 'sides' of
the cones, and not the 'lids' or any other weird face.
The script below uses simple logic to only keep the
combinations that are relevant.
"""
combos = list(it.combinations(np.arange(self.nverts), 3))
maxval = np.max(np.asarray(combos))
result = []
for combo in combos:
c1, c2, c3 = combo
if 0 not in combo:
continue
elif abs(c2 - c3) != 1:
if c1 == 0 and c2 == 1 and c3 == maxval:
result.append(combo)
else:
continue
else:
result.append(combo)
return result
if __name__ == '__main__':
# Initialize conical intersection
ci = ConicalIntersection(size=SIZE, asym=ASYM, nfaces=NFACES)
# Add lower and upper meshes
for name in ['Lower', 'Upper']:
mesh = bpy.data.meshes.new(name)
obj = bpy.data.objects.new(name, mesh)
bpy.context.collection.objects.link(obj)
mesh.from_pydata(ci.verts_lower if name == 'Lower' else ci.verts_upper,
[],
ci.faces)
Quick render after I have added the conical intersections, and a couple of planes:
Progress update
After some reading and Youtubing, I managed to make a relatively smooth connection. Here is my workflow:
- Add conical intersection using the script (
NFACES = 10), and join the two parts together - Add a plane mesh, and subdivide a couple of times
- Sculpt plane to roughly match the lower part of the conical intersection, and move the conical intersection so that it fully overlaps with the plane
- Add a sculpting mask to the entire conical intersection
- Add boolean modifier to plane mesh, and union join with the conical intersection
- Following the ideas here, I remove the vertices that previously belonged to the conical intersection but that now is located below the plane
- Clean up the resulting mesh by deleting dangling edges, recalculating normals, and merge vertices within a certain distance.
- Then smooth the mesh in sculpt mode to make the transition smooth (here the mask becomes important, because I want to retain a certain "physical correctness" by not distorting the conical intersection topography from the model).
Then I did a quick render, and here is the result:
I think I could use this method to achieve what I want. I just merge with a much larger plane that contains some topographical features of its own. Then repeat with the upper surface.
However, there are some artefacts from the merging. Also, notice that the intersection point is almost black, even though the light source should put light on those areas. Perhaps some flipped faces?
Again, any feedback would be much appreciated! I really don't know what I am doing here :)
Progress update 2
I have sculpted the energy landscapes, since I have a very specific shape in mind (to resemble published research). Screenshot:
Here is my blender file: 
