diff --git a/QA/grad-vector-rectangle-rotate.py b/QA/grad-vector-rectangle-rotate.py
index 488f65b..9447064 100644
--- a/QA/grad-vector-rectangle-rotate.py
+++ b/QA/grad-vector-rectangle-rotate.py
@@ -41,7 +41,7 @@ def forward(length1, length2, angle, c1):
     length1 = length1.type(torch.complex128)
     length2 = length2.type(torch.complex128)
 
-    obj_list = mee.rectangle_rotate(*c1, length1, length2, 5, 5, n_index, angle)
+    obj_list = mee.rectangle(*c1, length1, length2, 5, 5, n_index, angle)
 
     layer_info_list = [[layer_base, obj_list]]
 
diff --git a/QA/grad-vector.py b/QA/grad-vector.py
index d3fd35a..f36bb45 100644
--- a/QA/grad-vector.py
+++ b/QA/grad-vector.py
@@ -42,8 +42,8 @@ def forward(input_length1, input_length2, input_length3, input_length4, period,
     length3 = input_length3.type(torch.complex128)
     length4 = input_length4.type(torch.complex128)
 
-    obj1_list = ModelingTorch.rectangle(*c1, length1, length2, n_index)
-    obj2_list = ModelingTorch.rectangle(*c2, length3, length4, n_index + 2)
+    obj1_list = ModelingTorch.rectangle_no_approximation(*c1, length1, length2, n_index)
+    obj2_list = ModelingTorch.rectangle_no_approximation(*c2, length3, length4, n_index + 2)
 
     obj_list = obj1_list + obj2_list
 
diff --git a/meent/on_torch/modeler/modeling.py b/meent/on_torch/modeler/modeling.py
index 071c830..20c0314 100644
--- a/meent/on_torch/modeler/modeling.py
+++ b/meent/on_torch/modeler/modeling.py
@@ -40,59 +40,58 @@ def __init__(self, period=None, *args, **kwargs):
         self.period = period
 
     @staticmethod
-    def rectangle(cx, cy, lx, ly, base):
-
-        # a = torch.hstack([cy - ly / 2, cx - lx / 2])  # row, col
-        # b = torch.hstack([cy + ly / 2, cx + lx / 2])  # row, col
-        #
-        # res = [[[a, b, base]]]  # top_left, bottom_right
+    def rectangle_no_approximation(cx, cy, lx, ly, base):
 
         a = [cy - ly / 2, cx - lx / 2]  # row, col
         b = [cy + ly / 2, cx + lx / 2]  # row, col
 
         res = [[a, b, base]]  # top_left, bottom_right
 
-
         return res
 
-    def rectangle_rotate(self, cx, cy, lx, ly, n_split_triangle, n_split_parallelogram, n_index, angle=None, angle_margin=1E-5):
-        # if type(lx) in (int, float):
-        #     lx = torch.tensor(lx)
-        # if type(ly) in (int, float):
-        #     ly = torch.tensor(ly)
-        # if type(angle) in (int, float):
-        #     angle = torch.tensor(angle)
+    def rectangle(self, cx, cy, lx, ly, n_index, angle=0, n_split_triangle=2, n_split_parallelogram=2, angle_margin=1E-5):
+
         if type(lx) in (int, float):
-            lx = torch.tensor([lx])
+            lx = torch.tensor(lx).reshape(1)
+        elif type(lx) is torch.Tensor:
+            lx = lx.reshape(1)
+
         if type(ly) in (int, float):
-            ly = torch.tensor([ly])
+            ly = torch.tensor(ly).reshape(1)
+        elif type(ly) is torch.Tensor:
+            ly = ly.reshape(1)
+
         if type(angle) in (int, float):
-            angle = torch.tensor([angle])
+            angle = torch.tensor(angle).reshape(1)
+        elif type(angle) is torch.Tensor:
+            angle = angle.reshape(1)
 
         if lx.type not in (torch.complex64, torch.complex128):
             lx = lx.type(self.type_complex)  # TODO
         if ly.type not in (torch.complex64, torch.complex128):
             ly = ly.type(self.type_complex)
 
-        n_split_triangle, n_split_parallelogram = n_split_triangle + 2, n_split_parallelogram + 2
-
-        if angle is None:
-            angle = torch.tensor(0 * torch.pi / 180)
+        # n_split_triangle, n_split_parallelogram = n_split_triangle + 2, n_split_parallelogram + 2
 
-        # if 0 * torch.pi / 2 - angle_margin <= abs(angle) % (2 * torch.pi) <= 0 * torch.pi / 2 + angle_margin:
-        #     return self.rectangle(cx, cy, lx, ly, n_index)
-        # elif 1 * torch.pi / 2 - angle_margin <= abs(angle) % (2 * torch.pi) <= 1 * torch.pi / 2 + angle_margin:
-        #     return self.rectangle(cx, cy, ly, lx, n_index)
-        # elif 2 * torch.pi / 2 - angle_margin <= abs(angle) % (2 * torch.pi) <= 2 * torch.pi / 2 + angle_margin:
-        #     return self.rectangle(cx, cy, lx, ly, n_index)
-        # elif 3 * torch.pi / 2 - angle_margin <= abs(angle) % (2 * torch.pi) <= 3 * torch.pi / 2 + angle_margin:
-        #     return self.rectangle(cx, cy, ly, lx, n_index)
-        # else:
-        #     pass
+        # if angle is None:
+        #     angle = torch.tensor(0 * torch.pi / 180)
 
         angle = angle % (2 * torch.pi)
 
-        rotate = torch.ones((2, 2), dtype=torch.complex128)
+        # No rotation
+        if 0 * torch.pi / 2 - angle_margin <= abs(angle) % (2 * torch.pi) <= 0 * torch.pi / 2 + angle_margin:
+            return self.rectangle_no_approximation(cx, cy, lx, ly, n_index)
+        elif 1 * torch.pi / 2 - angle_margin <= abs(angle) % (2 * torch.pi) <= 1 * torch.pi / 2 + angle_margin:
+            return self.rectangle_no_approximation(cx, cy, ly, lx, n_index)
+        elif 2 * torch.pi / 2 - angle_margin <= abs(angle) % (2 * torch.pi) <= 2 * torch.pi / 2 + angle_margin:
+            return self.rectangle_no_approximation(cx, cy, lx, ly, n_index)
+        elif 3 * torch.pi / 2 - angle_margin <= abs(angle) % (2 * torch.pi) <= 3 * torch.pi / 2 + angle_margin:
+            return self.rectangle_no_approximation(cx, cy, ly, lx, n_index)
+        else:
+            pass
+
+        # Yes rotation
+        rotate = torch.ones((2, 2), dtype=self.type_complex)
         rotate[0, 0] = torch.cos(angle)
         rotate[0, 1] = -torch.sin(angle)
         rotate[1, 0] = torch.sin(angle)
@@ -171,129 +170,425 @@ def rectangle_rotate(self, cx, cy, lx, ly, n_split_triangle, n_split_parallelogr
 
         # point in region 1(top1~top2), 2(top2~top3) and 3(top3~top4)
 
-        xxx, yyy = [], []
-        xxx_cp, yyy_cp = [], []
+        # xxx, yyy = [], []
+        # xxx_cp, yyy_cp = [], []
         if top2_left:
 
             length = length_top12 / torch.sin(angle_inside)
             top3_cp = [top3[0] - length, top3[1]]
 
-            for i in range(n_split_triangle + 1):
-                x = top1[0] - (top1[0] - top2[0]) / n_split_triangle * i
-                y = top1[1] - (top1[1] - top2[1]) / n_split_parallelogram * i
-                xxx.append(x)
-                yyy.append(y)
+            # for i in range(n_split_triangle + 1):
+            #     x = top1[0] - (top1[0] - top2[0]) / n_split_triangle * i
+            #     y = top1[1] - (top1[1] - top2[1]) / n_split_parallelogram * i
+            #     xxx.append(x)
+            #     yyy.append(y)
+            #
+            #     xxx_cp.append(x + length / n_split_triangle * i)
+            #     yyy_cp.append(y)
+            #
+            # for i in range(n_split_parallelogram + 1):
+            #
+            #     x = top2[0] + (top3_cp[0] - top2[0]) / n_split_triangle * i
+            #     y = top2[1] - (top2[1] - top3_cp[1]) / n_split_parallelogram * i
+            #     xxx.append(x)
+            #     yyy.append(y)
+            #
+            #     xxx_cp.append(x + length)
+            #     yyy_cp.append(y)
+            #
+            # for i in range(n_split_triangle + 1):
+            #     x = top3_cp[0] + (top4[0] - top3_cp[0]) / n_split_triangle * i
+            #     y = top3_cp[1] - (top3_cp[1] - top4[1]) / n_split_parallelogram * i
+            #     xxx.append(x)
+            #     yyy.append(y)
+            #
+            #     xxx_cp.append(x + length / n_split_triangle * (n_split_triangle - i))
+            #     yyy_cp.append(y)
+
+            # 1: Upper triangle
+            xxx1 = top1[0] - (top1[0] - top2[0]) / n_split_triangle * torch.arange(n_split_triangle+1).reshape((-1, 1))
+            yyy1 = top1[1] - (top1[1] - top2[1]) / n_split_parallelogram * torch.arange(n_split_triangle+1).reshape((-1, 1))
+            xxx_cp1 = xxx1 + length / n_split_triangle * torch.arange(n_split_triangle+1).reshape((-1, 1))
+            yyy_cp1 = yyy1 * torch.ones(n_split_triangle+1).reshape((-1, 1))
+
+            # 2: Mid parallelogram
+            xxx2 = top2[0] + (top3_cp[0] - top2[0]) / n_split_triangle * torch.arange(n_split_parallelogram+1).reshape((-1, 1))
+            yyy2 = top2[1] - (top2[1] - top3_cp[1]) / n_split_parallelogram * torch.arange(n_split_parallelogram+1).reshape((-1, 1))
+            xxx_cp2 = (xxx2 + length) * torch.ones(n_split_parallelogram+1).reshape((-1, 1))
+            yyy_cp2 = yyy2 * torch.ones(n_split_parallelogram+1).reshape((-1, 1))
+
+            # 3: Lower triangle
+            xxx3 = top3_cp[0] + (top4[0] - top3_cp[0]) / n_split_triangle * torch.arange(n_split_triangle + 1).reshape(
+                (-1, 1))
+            yyy3 = top3_cp[1] - (top3_cp[1] - top4[1]) / n_split_parallelogram * torch.arange(n_split_triangle + 1).reshape(
+                (-1, 1))
+
+            xxx_cp3 = xxx3 + length / n_split_triangle * torch.arange(n_split_triangle, -1, -1).reshape((-1, 1))
+            yyy_cp3 = yyy3 * torch.ones(n_split_triangle + 1).reshape((-1, 1))
+
+            xxx = torch.concat((xxx1, xxx2, xxx3))
+            yyy = torch.concat((yyy1, yyy2, yyy3))
+
+            xxx_cp = torch.concat((xxx_cp1, xxx_cp2, xxx_cp3))
+            yyy_cp = torch.concat((yyy_cp1, yyy_cp2, yyy_cp3))
+
+            # # #####
+            #
+            # t00 = time.time()
+            # obj_list1 = []
+            #
+            # for i in range(len(xxx)):
+            #     if i == len(xxx) - 1:
+            #         break
+            #     x, y = xxx[i], yyy[i]
+            #     x_cp, y_cp = xxx_cp[i], yyy_cp[i]
+            #
+            #     x_next, y_next = xxx[i + 1], yyy[i + 1]
+            #     x_cp_next, y_cp_next = xxx_cp[i + 1], yyy_cp[i + 1]
+            #
+            #     x_mean = (x + x_next) / 2
+            #     x_cp_mean = (x_cp + x_cp_next) / 2
+            #     obj_list1.append([[y_cp_next, x_mean], [y, x_cp_mean], n_index])
+            # t01 = time.time()
+            #
+            #
+            # t0=time.time()
+            # obj_list1 = []
+            # x_mean_arr = (xxx + torch.roll(xxx, -1)) / 2
+            # x_cp_mean_arr = (xxx_cp + torch.roll(xxx_cp, -1)) / 2
+            # y_cp_next_arr = torch.roll(yyy_cp, -1)
+            #
+            # for i in range(len(xxx)-1):
+            #     obj_list1.append([[y_cp_next_arr[i], x_mean_arr[i]], [yyy[i], x_cp_mean_arr[i]], n_index])
+            #
+            # t1 =time.time()
+
+            x_mean_arr = (xxx + torch.roll(xxx, -1)) / 2
+            x_cp_mean_arr = (xxx_cp + torch.roll(xxx_cp, -1)) / 2
+            y_cp_next_arr = torch.roll(yyy_cp, -1)
+
+            obj_list1 = [[[y_cp_next_arr[i], x_mean_arr[i]], [yyy[i], x_cp_mean_arr[i]], n_index] for i in range(len(xxx)-1)]
+
+            # t2 =time.time()
+            # print(t01-t00, t1-t0, t2-t1)
+
+            # return obj_list1
 
-                xxx_cp.append(x + length / n_split_triangle * i)
-                yyy_cp.append(y)
+        else:
+            length = length_top12 / torch.cos(angle_inside)
+            top3_cp = [top3[0] + length, top3[1]]
 
-            for i in range(n_split_parallelogram + 1):
+            # 1: Top triangle
+            xxx1 = top1[0] + (top2[0] - top1[0]) / n_split_triangle * torch.arange(n_split_triangle + 1).reshape(
+                (-1, 1))
+            yyy1 = top1[1] - (top1[1] - top2[1]) / n_split_parallelogram * torch.arange(n_split_triangle + 1).reshape(
+                (-1, 1))
+            xxx_cp1 = xxx1 - length / n_split_triangle * torch.arange(n_split_triangle + 1).reshape((-1, 1))
+            yyy_cp1 = yyy1 * torch.ones(n_split_triangle + 1).reshape((-1, 1))
+
+            # for i in range(n_split_triangle + 1):
+            #     x = top1[0] + (top2[0] - top1[0]) / n_split_triangle * i
+            #     y = top1[1] - (top1[1] - top2[1]) / n_split_parallelogram * i
+            #     xxx.append(x)
+            #     yyy.append(y)
+            #
+            #     xxx_cp.append(x - length / n_split_triangle * i)
+            #     yyy_cp.append(y)
+
+            # 2: Mid parallelogram
+            xxx2 = top2[0] - (top2[0] - top3_cp[0]) / n_split_triangle * torch.arange(
+                n_split_parallelogram + 1).reshape((-1, 1))
+            yyy2 = top2[1] - (top2[1] - top3_cp[1]) / n_split_parallelogram * torch.arange(
+                n_split_parallelogram + 1).reshape((-1, 1))
+            xxx_cp2 = xxx2 - length * torch.ones(n_split_parallelogram + 1).reshape((-1, 1))
+            yyy_cp2 = yyy2 * torch.ones(n_split_parallelogram + 1).reshape((-1, 1))
+
+            # for i in range(n_split_parallelogram + 1):
+            #
+            #     x = top2[0] - (top2[0] - top3_cp[0]) / n_split_triangle * i
+            #     y = top2[1] - (top2[1] - top3_cp[1]) / n_split_parallelogram * i
+            #     xxx.append(x)
+            #     yyy.append(y)
+            #
+            #     xxx_cp.append(x - length)
+            #     yyy_cp.append(y)
+
+            # 3: Lower triangle
+            xxx3 = top3_cp[0] - (top3_cp[0] - top4[0]) / n_split_triangle * torch.arange(n_split_triangle + 1).reshape(
+                (-1, 1))
+            yyy3 = top3_cp[1] - (top3_cp[1] - top4[1]) / n_split_parallelogram * torch.arange(
+                n_split_triangle + 1).reshape(
+                (-1, 1))
+
+            xxx_cp3 = xxx3 - length / n_split_triangle * torch.arange(n_split_triangle, -1, -1).reshape((-1, 1))
+            yyy_cp3 = yyy3 * torch.ones(n_split_triangle + 1).reshape((-1, 1))
+
+            xxx = torch.concat((xxx1, xxx2, xxx3))
+            yyy = torch.concat((yyy1, yyy2, yyy3))
+
+            xxx_cp = torch.concat((xxx_cp1, xxx_cp2, xxx_cp3))
+            yyy_cp = torch.concat((yyy_cp1, yyy_cp2, yyy_cp3))
+
+            # for i in range(n_split_triangle + 1):
+            #     x = top3_cp[0] - (top3_cp[0] - top4[0]) / n_split_triangle * i
+            #     y = top3_cp[1] - (top3_cp[1] - top4[1]) / n_split_parallelogram * i
+            #     xxx.append(x)
+            #     yyy.append(y)
+            #
+            #     xxx_cp.append(x - length / n_split_triangle * (n_split_triangle - i))
+            #     yyy_cp.append(y)
+
+            x_mean_arr = (xxx + torch.roll(xxx, -1)) / 2
+            x_cp_mean_arr = (xxx_cp + torch.roll(xxx_cp, -1)) / 2
+            y_cp_next_arr = torch.roll(yyy_cp, -1)
+
+            obj_list1 = [[[y_cp_next_arr[i], x_cp_mean_arr[i]], [yyy[i], x_mean_arr[i]], n_index] for i in
+                         range(len(xxx) - 1)]
+
+        # obj_list1 = []
+        #
+        # for i in range(len(xxx)):
+        #     if i == len(xxx) - 1:
+        #         break
+        #     x, y = xxx[i], yyy[i]
+        #     x_cp, y_cp = xxx_cp[i], yyy_cp[i]
+        #
+        #     x_next, y_next = xxx[i + 1], yyy[i + 1]
+        #     x_cp_next, y_cp_next = xxx_cp[i + 1], yyy_cp[i + 1]
+        #
+        #     x_mean = (x + x_next) / 2
+        #     x_cp_mean = (x_cp + x_cp_next) / 2
+        #     obj_list1.append([[y_cp_next, x_cp_mean], [y, x_mean], n_index])
+
+        return obj_list1
+
+    # def ellipse(self, cx, cy, lx, ly, n_index, angle=0, n_split_w=2, n_split_h=2, angle_margin=1E-5):
+    #
+    #     if type(lx) in (int, float):
+    #         lx = torch.tensor(lx).reshape(1)
+    #     elif type(lx) is torch.Tensor:
+    #         lx = lx.reshape(1)
+    #
+    #     if type(ly) in (int, float):
+    #         ly = torch.tensor(ly).reshape(1)
+    #     elif type(ly) is torch.Tensor:
+    #         ly = ly.reshape(1)
+    #
+    #     if type(angle) in (int, float):
+    #         angle = torch.tensor(angle).reshape(1)
+    #     elif type(angle) is torch.Tensor:
+    #         angle = angle.reshape(1)
+    #
+    #     if lx.type not in (torch.complex64, torch.complex128):
+    #         lx = lx.type(self.type_complex)  # TODO
+    #     if ly.type not in (torch.complex64, torch.complex128):
+    #         ly = ly.type(self.type_complex)
+    #
+    #     angle = angle % (2 * torch.pi)
+    #
+    #
+    #     points_x_origin = lx/2 * torch.cos(torch.linspace(torch.pi/2, 0, n_split_w))
+    #     points_y_origin = ly/2 * torch.sin(torch.linspace(-torch.pi/2, torch.pi/2, n_split_h))
+    #
+    #     vv = torch.sin(torch.linspace(torch.pi/2, 0, n_split_w))
+    #     # hh = torch.cos(torch.linspace(torch.pi, 0, n_split_w))
+    #
+    #     # horizontal_length = cx + lx/2 * hh
+    #     vertical_length = ly/2 * vv
+    #
+    #     axis_x_origin = torch.vstack([points_x_origin, torch.ones(len(points_x_origin))])
+    #     axis_y_origin = torch.vstack([torch.ones(len(points_y_origin)), points_y_origin])
+    #
+    #     rotate = torch.ones((2, 2), dtype=points_x_origin.dtype)
+    #     rotate[0, 0] = torch.cos(angle)
+    #     rotate[0, 1] = -torch.sin(angle)
+    #     rotate[1, 0] = torch.sin(angle)
+    #     rotate[1, 1] = torch.cos(angle)
+    #
+    #     axis_x_origin_rot = rotate @ axis_x_origin
+    #     axis_y_origin_rot = rotate @ axis_y_origin
+    #
+    #
+    #     axis_x_rot = axis_x_origin_rot[:,:,None]
+    #     axis_x_rot[0] += cx
+    #     axis_x_rot[1] += cy
+    #
+    #     axis_y_rot = axis_y_origin_rot[:,:,None]
+    #     axis_y_rot[0] += cx
+    #     axis_y_rot[1] += cy
+    #
+    #     # points_origin_contour_rot = rotate @ points_origin_contour
+    #     # points_contour_rot = points_origin_contour_rot[:, :, None]
+    #     # points_contour_rot[0] += cx
+    #     # points_contour_rot[1] += cy
+    #
+    #     points = rotate @ torch.vstack((points_x_origin, points_y_origin))
+    #     points[0] += cx
+    #     points[1] += cy
+    #
+    #     # import matplotlib.pyplot as plt
+    #     #
+    #     # plt.scatter(*axis_x_rot.detach().numpy())
+    #     # plt.scatter(*axis_y_rot.detach().numpy())
+    #
+    #     # plt.scatter(*points_contour_rot.detach().numpy())
+    #
+    #     points = points[:, :, None]
+    #
+    #     # res = [[[points[1][i], points[0][i]], [points[1][i+1], points[0][i+1]], n_index] for i in range(len(points[0])-1)]
+    #
+    #     res = [[[axis_x_rot[1][i] - vertical_length[i]/torch.cos(angle),
+    #              axis_x_rot[0][i]],
+    #             [axis_x_rot[1][i] + vertical_length[i]/torch.cos(angle),
+    #              axis_x_rot[0][i+1]],
+    #             n_index] for i in range(len(points[0])-1)]
+    #
+    #     res = [[[axis_x_rot[1][i] - vertical_length[i]/torch.cos(angle),
+    #              axis_x_rot[0][i]],
+    #             [axis_x_rot[1][i] + vertical_length[i]/torch.cos(angle),
+    #              axis_x_rot[0][i+1]],
+    #             n_index] for i in range(len(points[0])-1)]
+    #
+    #     ress = []
+    #     for i in range(len(axis_x_rot[0])-1):
+    #         LL = [axis_x_rot[1][i] - vertical_length[i], axis_x_rot[0][i]]
+    #         UR = [axis_x_rot[1][i] + vertical_length[i], axis_x_rot[0][i+1]]
+    #
+    #         center_x = (LL[1] + UR[1])/2
+    #         center_y = (LL[0] + UR[0])/2
+    #
+    #         alpha = center_x - cx
+    #         beta = center_y - cy
+    #
+    #         # reflection over the origin
+    #         LL_pair = [LL[0] - 2*beta, LL[1] - 2*alpha]
+    #         UR_pair = [UR[0] - 2*beta, UR[1] - 2*alpha]
+    #
+    #         ress.append([LL, UR, n_index])
+    #         ress.append([LL_pair, UR_pair, n_index])
+    #
+    #
+    #     return ress, (axis_x_rot, axis_y_rot)
+
+    def ellipse(self, cx, cy, lx, ly, n_index, angle=0, n_split_w=2, n_split_h=2, angle_margin=1E-5, debug=False):
+
+        if type(lx) in (int, float):
+            lx = torch.tensor(lx).reshape(1)
+        elif type(lx) is torch.Tensor:
+            lx = lx.reshape(1)
 
-                x = top2[0] + (top3_cp[0] - top2[0]) / n_split_triangle * i
-                y = top2[1] - (top2[1] - top3_cp[1]) / n_split_parallelogram * i
-                xxx.append(x)
-                yyy.append(y)
+        if type(ly) in (int, float):
+            ly = torch.tensor(ly).reshape(1)
+        elif type(ly) is torch.Tensor:
+            ly = ly.reshape(1)
 
-                xxx_cp.append(x + length)
-                yyy_cp.append(y)
+        if type(angle) in (int, float):
+            angle = torch.tensor(angle).reshape(1)
+        elif type(angle) is torch.Tensor:
+            angle = angle.reshape(1)
 
-            for i in range(n_split_triangle + 1):
-                x = top3_cp[0] + (top4[0] - top3_cp[0]) / n_split_triangle * i
-                y = top3_cp[1] - (top3_cp[1] - top4[1]) / n_split_parallelogram * i
-                xxx.append(x)
-                yyy.append(y)
+        if lx.type not in (torch.complex64, torch.complex128):
+            lx = lx.type(self.type_complex)  # TODO
+        if ly.type not in (torch.complex64, torch.complex128):
+            ly = ly.type(self.type_complex)
 
-                xxx_cp.append(x + length / n_split_triangle * (n_split_triangle - i))
-                yyy_cp.append(y)
+        angle = angle % (2 * torch.pi)
 
-            obj_list1 = []
+        points_x_origin = lx/2 * torch.cos(torch.linspace(torch.pi/2, 0, n_split_w))
+        points_y_origin = ly/2 * torch.sin(torch.linspace(-torch.pi/2, torch.pi/2, n_split_h))
 
-            for i in range(len(xxx)):
-                if i == len(xxx) - 1:
-                    break
-                x, y = xxx[i], yyy[i]
-                x_cp, y_cp = xxx_cp[i], yyy_cp[i]
+        points_x_origin_contour = lx/2 * torch.cos(torch.linspace(-torch.pi, torch.pi, n_split_w))[:-1]
+        points_y_origin_contour = ly/2 * torch.sin(torch.linspace(-torch.pi, torch.pi, n_split_h))[:-1]
+        points_origin_contour = torch.vstack([points_x_origin_contour, points_y_origin_contour])
 
-                x_next, y_next = xxx[i + 1], yyy[i + 1]
-                x_cp_next, y_cp_next = xxx_cp[i + 1], yyy_cp[i + 1]
+        axis_x_origin = torch.vstack([points_x_origin, torch.ones(len(points_x_origin))])
+        axis_y_origin = torch.vstack([torch.ones(len(points_y_origin)), points_y_origin])
 
-                x_mean = (x + x_next) / 2
-                x_cp_mean = (x_cp + x_cp_next) / 2
-                obj_list1.append([[y_cp_next, x_mean], [y, x_cp_mean], n_index])
+        rotate = torch.ones((2, 2), dtype=points_x_origin.dtype)
+        rotate[0, 0] = torch.cos(angle)
+        rotate[0, 1] = -torch.sin(angle)
+        rotate[1, 0] = torch.sin(angle)
+        rotate[1, 1] = torch.cos(angle)
 
-            return obj_list1
+        axis_x_origin_rot = rotate @ axis_x_origin
+        axis_y_origin_rot = rotate @ axis_y_origin
 
-        else:
+        axis_x_rot = axis_x_origin_rot[:, :, None]
+        axis_x_rot[0] += cx
+        axis_x_rot[1] += cy
 
-            length = length_top12 / torch.cos(angle_inside)
-            top3_cp = [top3[0] + length, top3[1]]
+        axis_y_rot = axis_y_origin_rot[:, :, None]
+        axis_y_rot[0] += cx
+        axis_y_rot[1] += cy
 
-            for i in range(n_split_triangle + 1):
-                x = top1[0] + (top2[0] - top1[0]) / n_split_triangle * i
-                y = top1[1] - (top1[1] - top2[1]) / n_split_parallelogram * i
-                xxx.append(x)
-                yyy.append(y)
+        points_origin_contour_rot = rotate @ points_origin_contour
+        points_contour_rot = points_origin_contour_rot[:, :, None]
+        points_contour_rot[0] += cx
+        points_contour_rot[1] += cy
 
-                xxx_cp.append(x - length / n_split_triangle * i)
-                yyy_cp.append(y)
+        y_highest_index = torch.argmax(points_contour_rot.real, dim=1)[1, 0]
 
-            for i in range(n_split_parallelogram + 1):
+        points_contour_rot = torch.roll(points_contour_rot, (points_contour_rot.shape[1] // 2 - y_highest_index).item(), dims=1)
+        y_highest_index = torch.argmax(points_contour_rot.real, dim=1)[1, 0]
 
-                x = top2[0] - (top2[0] - top3_cp[0]) / n_split_triangle * i
-                y = top2[1] - (top2[1] - top3_cp[1]) / n_split_parallelogram * i
-                xxx.append(x)
-                yyy.append(y)
+        right = points_contour_rot[:, y_highest_index-1]
+        left = points_contour_rot[:, y_highest_index+1]
 
-                xxx_cp.append(x - length)
-                yyy_cp.append(y)
+        right_y = right[1].real
+        left_y = left[1].real
 
-            for i in range(n_split_triangle + 1):
-                x = top3_cp[0] - (top3_cp[0] - top4[0]) / n_split_triangle * i
-                y = top3_cp[1] - (top3_cp[1] - top4[1]) / n_split_parallelogram * i
-                xxx.append(x)
-                yyy.append(y)
+        left_array = []
+        right_array = []
 
-                xxx_cp.append(x - length / n_split_triangle * (n_split_triangle - i))
-                yyy_cp.append(y)
+        res = []
 
-            obj_list1 = []
+        if left_y > right_y:
+            right_array.append(points_contour_rot[:, y_highest_index])
+        elif left_y < right_y:
+            left_array.append(points_contour_rot[:, y_highest_index])
+        else:
+            raise ValueError  # TODO
 
-            for i in range(len(xxx)):
-                if i == len(xxx) - 1:
-                    break
-                x, y = xxx[i], yyy[i]
-                x_cp, y_cp = xxx_cp[i], yyy_cp[i]
+        for i in range(points_contour_rot.shape[1]//2):
+            left_array.append(points_contour_rot[:, (y_highest_index+i+1) % points_contour_rot.shape[1]])
+            right_array.append(points_contour_rot[:, (y_highest_index-i-1) % points_contour_rot.shape[1]])
 
-                x_next, y_next = xxx[i + 1], yyy[i + 1]
-                x_cp_next, y_cp_next = xxx_cp[i + 1], yyy_cp[i + 1]
+        arr = torch.zeros((2, len(right_array) + len(left_array), 1), dtype=points_contour_rot.dtype)
 
-                x_mean = (x + x_next) / 2
-                x_cp_mean = (x_cp + x_cp_next) / 2
-                obj_list1.append([[y_cp_next, x_cp_mean], [y, x_mean], n_index])
+        if left_y > right_y:
+            arr[:, ::2] = torch.stack(right_array, dim=1)
+            arr[:, 1::2] = torch.stack(left_array, dim=1)
+        elif left_y < right_y:
+            arr[:, ::2] = torch.stack(left_array, dim=1)
+            arr[:, 1::2] = torch.stack(right_array, dim=1)
 
-            return obj_list1
+        arr_roll = torch.roll(arr, -1, 1)
 
-    def ellipse(self, cx, cy, lx, ly, dx, dy, base, rotation_angle=0 * torch.pi / 180):
-        points_x = torch.arange(cx-lx, cx+lx, dx+2)
-        points_y = torch.arange(cy-ly, cy+ly, dy+2)
+        for i in range(arr.shape[1]):
+            ax, ay = arr[:, i]
+            bx, by = arr_roll[:, i]
 
-        rotate = torch.tensor([torch.cos(rotation_angle), -torch.sin(rotation_angle)],
-                              [torch.sin(rotation_angle),  torch.cos(rotation_angle)],
-                              )
+            LL = [min(ay.real, by.real)+0j, min(ax.real, bx.real)+0j]
+            UR = [max(ay.real, by.real)+0j, max(ax.real, bx.real)+0j]
 
-        points = rotate @ torch.vstack((points_x, points_y))
-        res = [points[0], points[1], base]
-        return res
+            res.append([LL, UR, n_index])
+
+        if debug:
+            return res[:-1], (axis_x_rot, axis_y_rot, points_contour_rot)
+        else:
+            return res[:-1]
 
     def vector(self, layer_info, x64=True):
 
         # TODO: activate and apply 'x64' option thru this function and connect to meent class.
         if x64:
             datatype = torch.complex128
-            perturbation = 1E-14
+            perturbation = 0
+            perturbation_unit = 1E-14
         else:
             datatype = torch.complex64
-            perturbation = 1E-6
+            perturbation = 0
+            perturbation_unit = 1E-6
 
         pmtvy_base, obj_list = layer_info
 
@@ -301,6 +596,7 @@ def vector(self, layer_info, x64=True):
         row_list = []
         col_list = []
 
+        # overlap check and apply perturbation
         for obj in obj_list:
             top_left, bottom_right, _ = obj
 
@@ -308,11 +604,15 @@ def vector(self, layer_info, x64=True):
             for _ in range(100):
                 index = bisect_left(row_list, top_left[0].real, key=lambda x: x.real)
                 if len(row_list) > index and top_left[0] == row_list[index]:
+                    perturbation += perturbation_unit
                     if top_left[0] == 0:
                         top_left[0] = top_left[0] + perturbation
+
                     else:
                         # top_left[0] = top_left[0] - (top_left[0] * perturbation)  # TODO: plus or minus?
-                        top_left[0] = top_left[0] + (top_left[0] * perturbation)
+                        top_left[0] = top_left[0] + (top_left[0] * perturbation)  # TODO: change; save how many perturbations were applied in a variable
+                    row_list.insert(index, top_left[0])
+                    break
                 else:
                     row_list.insert(index, top_left[0])
                     break
@@ -325,6 +625,7 @@ def vector(self, layer_info, x64=True):
             for _ in range(100):
                 index = bisect_left(row_list, bottom_right[0].real, key=lambda x: x.real)
                 if len(row_list) > index and bottom_right[0] == row_list[index]:
+                    perturbation += perturbation_unit
                     # if bottom_right[0] == 0:
                     #     bottom_right[0] = bottom_right[0] + perturbation
                     # else:
@@ -333,6 +634,9 @@ def vector(self, layer_info, x64=True):
 
                     # bottom_right[0] = bottom_right[0] + (bottom_right[0] * perturbation)
                     bottom_right[0] = bottom_right[0] - (bottom_right[0] * perturbation)
+                    row_list.insert(index, bottom_right[0])
+                    break
+
                 else:
                     row_list.insert(index, bottom_right[0])
                     break
@@ -345,11 +649,15 @@ def vector(self, layer_info, x64=True):
             for _ in range(100):
                 index = bisect_left(col_list, top_left[1].real, key=lambda x: x.real)
                 if len(col_list) > index and top_left[1] == col_list[index]:
+                    perturbation += perturbation_unit
+
                     if top_left[1] == 0:
                         top_left[1] = top_left[1] + perturbation
                     else:
                         # top_left[1] = top_left[1] - (top_left[1] * perturbation)
                         top_left[1] = top_left[1] + (top_left[1] * perturbation)
+                    col_list.insert(index, top_left[1])
+                    break
                 else:
                     col_list.insert(index, top_left[1])
                     break
@@ -362,6 +670,7 @@ def vector(self, layer_info, x64=True):
             for _ in range(100):
                 index = bisect_left(col_list, bottom_right[1].real, key=lambda x: x.real)
                 if len(col_list) > index and bottom_right[1] == col_list[index]:
+                    perturbation += perturbation_unit
                     # if bottom_right[1] == 0:
                     #     bottom_right[1] = bottom_right[1] + perturbation
                     # else:
@@ -370,6 +679,8 @@ def vector(self, layer_info, x64=True):
 
                     # bottom_right[1] = bottom_right[1] + (bottom_right[1] * perturbation)
                     bottom_right[1] = bottom_right[1] - (bottom_right[1] * perturbation)
+                    col_list.insert(index, bottom_right[1])
+                    break
                 else:
                     col_list.insert(index, bottom_right[1])
                     break
@@ -379,9 +690,9 @@ def vector(self, layer_info, x64=True):
                 col_list.insert(index, bottom_right[1])
 
         if not row_list or row_list[-1] != self.period[1]:
-            row_list.append(self.period[1])
+            row_list.append(self.period[1].reshape(1))
         if not col_list or col_list[-1] != self.period[0]:
-            col_list.append(self.period[0])
+            col_list.append(self.period[0].reshape(1))
 
         if row_list and row_list[0] == 0:
             row_list = row_list[1:]
@@ -409,14 +720,8 @@ def vector(self, layer_info, x64=True):
 
             ucell_layer[row_begin:row_end, col_begin:col_end] = pmty
 
-        x_list = torch.zeros((len(col_list), 1), dtype=datatype)
-        y_list = torch.zeros((len(row_list), 1), dtype=datatype)
-
-        for i in range(len(col_list)):
-            x_list[i] = col_list[i]
-
-        for i in range(len(row_list)):
-            y_list[i] = row_list[i]
+        x_list = torch.cat(col_list).reshape((-1, 1))
+        y_list = torch.cat(row_list).reshape((-1, 1))
 
         return ucell_layer, x_list, y_list
 
diff --git a/setup.py b/setup.py
index a021190..732292c 100644
--- a/setup.py
+++ b/setup.py
@@ -12,7 +12,7 @@
 }
 setup(
     name='meent',
-    version='0.9.8',
+    version='0.9.9',
     url='https://github.com/kc-ml2/meent',
     author='KC ML2',
     author_email='yongha@kc-ml2.com',