How to slice a 3D volume at a curved 2D surface using Python?

Question

Consider a closed surface mesh in 3D (mesh1), consisting of two merged blocks as shown in the image.

The mesh is given as STL file and is subdivided into different patches. Each patch is saved as a separate body in the STL file.

Further, I have a curved 2D surface mesh (mesh2) given by an STL file, which completely intersects the closed surface mesh and therefore divides it into an upper and lower part as shown in the next image:

How can I slice the closed surface at the given 2D surface while keeping the information about the patches?

For mesh manipulation I use the trimesh package in Python in combination with the Blender API, but I'm open for suggestions with other packages/programs. What I tried so far is the following process:

1. Load mesh1 and mesh2 in trimesh
2. Scale up mesh2 by 1% to guarantee that it completely intersects with mesh1 and save it
3. Applying the trimesh boolean difference modifier to mesh1 with mesh2 as tool

Using this method I get the trimesh error that mesh1 and mesh2 are no volumes. The next method I tried was the following: 3. Thicken mesh2 by a small amount using the Blender solidify modifier. 4. Merge vertices of mesh1 using trimesh.grouping.merge_vertices(mesh1). The problem here was, that trimesh doesn't consider mesh1 to be watertight since the faces between different patches are not connected. 5. Applying the trimesh boolean difference modifier to mesh1 with mesh2 as tool.

During step 4 of this method I lose the information about patch affiliation.

To conclude: Is there a way to split up the bigger body at a curved 2D surface while keeping the information about the patches?

Edit: I tried the method in @Andras Deak's comment with clip_surface() in pyvista which was successful for the simple clipping surface in the original post. Trying it with the "real" surface that I want to use it for was not successful until now. The "real" clipping surface is much more complex and looks like this:

It can be assumed that the clipping surface is always located in the upper block of the outer geometry. For simplicity and stability reason I split up the STL of the outer geometry and only calculated the clipped mesh between the upper block and the clipping surface. This looks like this:

At this point the clipping surface is scaled up by 1% to ensure intersection and then the clip_surface() command can be applied. The result is shown in the next image:

How can the gaps between the clipping surface and the clipped mesh be filled or how can the clipping accuracy be increased?