I need to compute ~2 million coordinates (64-bit float tuples) They're the intersection points of lines that form a parallelogram-shaped grid (See visual below) The shape of the parallelogram is different for every iteration. I have all the boundary coordinates (in green) at my disposition. The lines running east-west have constant latitudes. I am sure this has been solved many times but alas, I don't know what keywords I should be searching for to find the solution. Logical approaches I've considered so far: Calculate north-south line slopes (i.e.: (y2 - y1) / (x2 - x1)). Solve for x (i.e.: east-west line longitude at intersection). Calculate east-west line lengths. Subtract total constant widths (i.e. 6x narrow strips). Divide remainder by 6. Apply result as longitudinal offset for each point on east-west line. Tried using shapely's split() function to split the parallelogram-shaped polygon with line strings but it quickly became clear that it would be a nightmare to process the GeometryCollection that are returned and having to sort through them. I'm sure you guys know of much better ways though. Questions: What is the best mathematical computational approach to do this? and Are there any python libraries I could use to speed things up? Bonus: I will probably need to run this as a multi-process when the time comes to generate the dataset. Can you confirm that these types of computations are best done on the CPU cores and not the GPU ones? Continue reading...