HeteQuadrature/gjj/analyticSolution.py

import math


class Sphere:
    def __init__(self, radius, center):
        self.radius = radius
        self.center = center


# specify class for the input parameters
def caseSpheresBOnAFace(rA, rB):
    # intersection:
    rAsq = rA**2
    rBsq = rB**2
    crownHeightA = rBsq / (2 * rA)
    crownHeightB = rB - crownHeightA
    crownVolumeA = math.pi * crownHeightA**2 * (3 * rA - crownHeightA) / 3
    crownVolumeB = math.pi * crownHeightB**2 * (3 * rB - crownHeightB) / 3
    volumeInter = crownVolumeA + crownVolumeB

    # union:
    volumeA = 4 / 3 * math.pi * rAsq * rA
    volumeB = 4 / 3 * math.pi * rBsq * rB
    volumeUnion = volumeA + volumeB - volumeInter

    # difference:
    volumeDiff = volumeA - volumeInter

    return volumeInter, volumeUnion, volumeDiff


def caseExtrude():
    x1 = -32450.000001122418
    y1 = -3.8391231093737305e07
    x2 = 32449.999998877582
    y2 = -3.8392495305561452e07
    x3 = 33538.999998877582
    y3 = 3.8392534654100421e07
    x4 = -33539.000001122418
    y4 = 3.8391228016184551e07
    area = 0.5 * abs(
        x1 * y2 + x2 * y3 + x3 * y4 + x4 * y1 - (y1 * x2 + y2 * x3 + y3 * x4 + y4 * x1)
    )
    return area * 3300


if __name__ == "__main__":
    # volumeInter, volumeUnion, volumeDiff = caseSpheresBOnAFace(800000, 1000)
    # print("volumeInter = ", volumeInter)
    # print("volumeUnion = ", volumeUnion)
    # print("volumeDiff = ", volumeDiff)
    print(caseExtrude())