Codejam 2008 qualification round question C

其实每次codejam竞赛，qualification round的最后大题才是精华所在。因为qualification足有24个小时时间，最后一道题自然也必须要有一定复杂度。08年的题算是还好。

 

如果是数学牛人，也许一个积分式子就能解决问题，不过于在下这种烂工科生，这个分段积分有些太复杂了（真对不起数分老师他老人家啊）……无聊的时候想过用要求精度模拟无穷小量来分割运算，不过太慢 =_=

 

所以目前基本的思路是一象限内任给一个正方形，把它在圆心原点，半径r的圆内的面积计算出来，用圆的面积减去所有这些正方形的面积，再乘以四，就可以得到整个的面积。而正方形在圆内的面积总可以分为三种情况：1.全在圆内；2.全不在圆内；3.部分在圆内，此时面积可分割为多边形的面积和弓形的面积，于是可以计算了。

 

运行无问题，缺陷在于有些太慢，运行google提供的大数据输入花费6分钟左右，虽然在竞赛中还是能通过，不过应该还有很大优化的余地。除了优化算法，应该也可以改用多线程实现。Maybe全转换成C代码看看速度如何？嗯嗯……

更新：因为在每个象限内完全轴对称，所以只需算1/8的面积（半个象限）就可以了。改进后花费2分半左右，基本没问题了。哇哈哈哈

再次更新：吧calc_r函数从使用calc_segment_len()函数：

def calc_r(p):
    return calc_segment_len([p, [0, 0]])

改为直接计算：

def calc_r(p):
    return math.sqrt(p[0]**2+p[1]**2)

整个计算过程缩短到1分40秒左右……看来对一些使用频繁的操作的优化更有意义。

 

再次更新：尝试用map+lambda操作：

        inner_left_down_point = map(lambda x:x+frame_width, left_down_point)
        inner_right_up_point = map(lambda x:x+inner_square_side, inner_left_down_point)

替换原来的新建list然后直接相加的操作：

        inner_left_down_point = list(left_down_point)
        inner_left_down_point[0] += frame_width
        inner_left_down_point[1] += frame_width
        
        inner_right_up_point = list(inner_left_down_point)
        inner_right_up_point[0] += inner_square_side
        inner_right_up_point[1] += inner_square_side

计算时间延长了8秒左右。看来map+lambda在这种情况下效率略低。

 

'''
CodeJam Practice 2008 qualification round question B
Created on 2012-12-11

@author: festony
'''

import math
from cj_lib import *
from properties import *

#curr_file_name = 'C-small-practice'
curr_file_name = 'C-large-practice'
#curr_file_name = 'test'

# !!!!!!
# Here we only consider the situation in First Quadrant
# thus all points x >= 0, y >= 0

# point: [x, y]

# first we need to check whether a point is in circle or not
# circle center [0, 0], radius frame_width
# > 0: outside; < 0: inside; ==0: just on the circle
def point_in_circle(p, r):
    if p[0] < 0 or p[1] < 0:
        raise Exception('Point [%f, %f] not in first quadrant!' % tuple(p))
    d2 = p[0]**2 + p[1]**2
    d2 -= r**2
    if d2 == 0:
        return 0
    return d2 / math.fabs(d2)

# segment: [point1, point2], e.inner_side_len. [[0, 1], [1, 1]]

# then we need to check whether a segment is in/out circle, or they insect
# > 0: outside; < 0: inside; ==0: insect
def segment_in_circle(seg, r):
    d2 = point_in_circle(seg[0], r) * point_in_circle(seg[1], r)
    if d2 <= 0:
        return 0
    return d2 * point_in_circle(seg[0], r)

# then we need to calculate the insect point if segment insect with circle
# luckily we only need to consider horizontal and vertical segments here.
def segment_circle_insect_point(seg, r):
    if segment_in_circle(seg, r) != 0:
        raise Exception('Segment not insect with circle!')
    if seg[0][0] == seg[1][0]:
        # horizontal segment
        return [seg[0][0], math.sqrt(r**2 - seg[0][0]**2)]
    elif seg[0][1] == seg[1][1]:
        # vertical segment
        return [math.sqrt(r**2 - seg[0][1]**2), seg[0][1]]
    else:
        raise Exception('Segment not horizontal / vertical!')

# calculate length of seg
def calc_segment_len(seg):
    a = map(lambda x, y: float(x-y), seg[1], seg[0])
    return math.sqrt(a[0]**2+a[1]**2)

# calculate distance to origin point
def calc_r(p):
    return calc_segment_len([p, [0, 0]])

# calculate distance between mid-point of a seg and the origin point
def calc_distance_segment_mid(seg):
    mid = map(lambda x, y: float(x+y)/2, seg[1], seg[0])
    return math.sqrt(mid[0]**2+mid[1]**2)

# then we need to calculate the area of the circular segment cut off the circle
# by the segment.
def calc_circular_segment_area(seg, r):
    d = calc_segment_len(seg) / 2
    h = calc_distance_segment_mid(seg)
    fan_area = r**2 * math.acos(h/r)
    triangle_area = h * d
    return fan_area - triangle_area

# square : [point0, point1, point2, point3]
# first point is left-up vertex, then left down, right down, right up.
# anti-clockwise.
# could be created by left-up vertex and the side length
# also could access 4 segments: [seg0, seg1, seg2, seg3]
# first segment is left side, then down, right, up.
# also anti-clockwise.
class square:
    def __init__(self, left_down_point, side_len):
        if left_down_point[0] < 0 or left_down_point[1] < 0:
            str = 'The sequare should be in first quadrant! [%f, %f]' % tuple(left_down_point)
            raise Exception(str)
        self.inner_side_len = side_len
        self.points = []
        self.points.append([left_down_point[0], left_down_point[1] + side_len])
        self.points.append(left_down_point);
        self.points.append([left_down_point[0] + side_len, left_down_point[1]])
        self.points.append([left_down_point[0] + side_len, left_down_point[1] + side_len])
        self.sides = []
        self.sides.append([self.points[0], self.points[1]])
        self.sides.append([self.points[1], self.points[2]])
        self.sides.append([self.points[2], self.points[3]])
        self.sides.append([self.points[3], self.points[0]])
        
    # given 2 points on the sides, the segment of these 2 ponts cut square into 2 parts.
    # calculate the area left-down side as a polygon.
    def calc_half_area(self, seg):
        if seg[0][0] > seg[1][0]:
            seg.reverse()
        len0 = float(seg[0][0] - self.points[0][0])
        len1 = float(seg[1][0] - seg[0][0])
        len2 = float(seg[1][1] - self.points[1][1])
        len3 = float(seg[0][1] - seg[1][1])
        area = len0 * len2 + len1 * len2 + len0 * len3 + len1 * len3 / 2
        return area
        
    # then we need to check whether the square insects with circle or not.
    def check_insect_circle(self, r):
        d = sum(map(point_in_circle, self.points, [r]*4))
        if d >= 3:
            return 1
        elif d <= -3:
            return -1
        else:
            return 0
    # then we need to calculate 2 insection points if insects with circle.
    def insect_circle_area(self, r):
        if self.check_insect_circle(r) > 0:
            return 0
        elif self.check_insect_circle(r) < 0:
            return self.inner_side_len ** 2
        else:
            seg = []
            if segment_in_circle(self.sides[0], r) == 0:
                seg.append(segment_circle_insect_point(self.sides[0], r))
            else:
                seg.append(segment_circle_insect_point(self.sides[3], r))
                
            if segment_in_circle(self.sides[1], r) == 0:
                seg.append(segment_circle_insect_point(self.sides[1], r))
            else:
                seg.append(segment_circle_insect_point(self.sides[2], r))
                    
            return self.calc_half_area(seg) + calc_circular_segment_area(seg, r)
        
class wrapper_square:
    def __init__(self, left_down_point, inner_square_side, frame_width):
        self.inner_side_len = inner_square_side
        self.frame_width = frame_width
        inner_left_down_point = list(left_down_point)
        inner_left_down_point[0] += frame_width
        inner_left_down_point[1] += frame_width
        inner_right_up_point = list(inner_left_down_point)
        inner_right_up_point[0] += inner_square_side
        inner_right_up_point[1] += inner_square_side
        
        self.inner_square = square(inner_left_down_point, inner_square_side)
        self.inner_r = calc_r(inner_left_down_point)
        self.outer_r = calc_r(inner_right_up_point)
        
    def insect_circle_area(self, r):
        if r <= self.inner_r:
            return 0.0
        elif r >= self.outer_r:
            return self.inner_side_len**2
        else:
            return self.inner_square.insect_circle_area(r)

def input_dividing_func(input_lines):
    total_case = int(input_lines.pop(0))
    case_inputs = []
    for i in range(total_case):
        case_inputs += [map(float, input_lines.pop(0).split(' '))]
    return case_inputs

def process_func(func_input):
    f, R, t, r, g = func_input
    if f >= g/2:
        p = 1.0
    else:
#        inner_R = R - t - f
#        square_side = g + 2 * r
#        repeat_time = int(math.ceil(float(inner_R) / square_side))
#        blank_area = 0.0
#        for i in range(repeat_time):
#            for j in range(repeat_time):
#                point = [i * square_side + r + f, (j + 1) * square_side - r - f]
#                blank_area += square(point, g - 2 * f).insect_circle_area(inner_R)
#        blank_area *= 4

        # Because it's symmetric, we only need to do calc for half Quadrant (1/8 of the whole space).
        inner_R = R - t - f
        wrapper_square_side = g + 2 * r
        repeat_time = int(math.ceil(float(inner_R) / wrapper_square_side))
        blank_area = 0.0
        for i in range(repeat_time):
            for j in range(i + 1):
                point = [i * wrapper_square_side, j * wrapper_square_side]
                sub_blank_area = wrapper_square(point, g - 2 * f, r + f).insect_circle_area(inner_R)
                if i == j:
                    sub_blank_area /= 2
                blank_area += sub_blank_area
        blank_area *= 8
        
        p = (math.pi * (R**2) - blank_area) / (math.pi * (R**2))
    return '%.6f' % (p,)
        

#run_proc_m(process_func, input_dividing_func, curr_working_folder, curr_file_name)
run_proc(process_func, input_dividing_func, curr_working_folder, curr_file_name)