# 计算直线塔不平衡张力
# 新版输电线路大手册 P328

import math
import data
import numpy as np

# h_i 悬点高差
# l_i 悬点档距
# _alpha 导线膨胀系数  1/°C
# _elastic 弹性系数 N/mm2
# _t_e 架线时考虑初伸长的降温，取正值。单位°C
# lambda_i 计算不平衡张力时导线比载 N/(m.mm)
# sigma_i 计算不平衡张力时最低点水平应力 单位N/mm2
# t_i 计算不平衡张力时导线温度 单位°C
# _lambda_m 导线架线时时导线比载 N/(m.mm)
# _sigma_m 导线架线时时最低点水平应力 单位N/mm2
# _t_m 导线架线时时导线温度 单位°C


def delta_li(
    h_i: float,
    l_i,
    lambda_i: float,
    _alpha: float,
    _elastic: float,
    _t_e: float,
    t_i: float,
    sigma_i,
    _lambda_m: float,
    _t_m: float,
    _sigma_m: float,
) -> float:
    beta_i = math.atan(h_i / l_i)  # 高差角
    _delta_li = (
        l_i
        / ((math.cos(beta_i) ** 2) * (1 + (lambda_i * l_i / sigma_i) ** 2 / 8))
        * (
            (l_i * math.cos(beta_i)) ** 2
            / 24
            * ((_lambda_m / _sigma_m) ** 2 - (lambda_i / sigma_i) ** 2)
            + ((sigma_i - _sigma_m) / _elastic / math.cos(beta_i))
            + _alpha * (t_i + _t_e - _t_m)
        )
    )
    return _delta_li


# area 导线截面 单位mm2
# sigma_i 第i档内水平应力 单位N/mm2
# b_i 悬垂串沿线路方向水平偏移距离，沿大号方向为正，反之为负。 单位m
# length_i 第i基直线塔串长 单位m
# g_i 第i基直线塔串重 单位N
# h_i 悬垂串处千中垂位置时，,第 i 基对第 i-1 杆塔上导线悬挂点间的高差大号比小号杆塔悬挂点高者h本身为正值，反之为负值。
# lambda_i 第i档导线比载 N/(m.mm)
# h_i1 悬垂串处千中垂位置时，第 i+1 基对第 i 杆塔上导线悬挂点间的高差大号比小号杆塔悬挂点高者h本身为正值，反之为负值。
# lambda_i1
def fun_sigma_i1(
    area: float,
    sigma_i: float,
    b_i: float,
    length_i: float,
    g_i: float,
    h_i: float,
    l_i: float,
    lambda_i: float,
    h_i1: float,
    l_i1: float,
    lambda_i1: float,
) -> float:
    beta_i = math.atan(h_i / l_i)
    beta_i1 = math.atan(h_i1 / l_i1)
    _sigma_i1 = (
        (
            g_i / 2 / area  # g_i传入时已考虑导线分裂数
            + lambda_i * l_i / 2 / math.cos(beta_i)
            + lambda_i1 * l_i1 / 2 / math.cos(beta_i1)
            + sigma_i * h_i / l_i
        )
        + sigma_i / b_i * math.sqrt(length_i ** 2 - b_i ** 2)
    ) / (math.sqrt(length_i ** 2 - b_i ** 2) / b_i + h_i1 / l_i1)
    return _sigma_i1


# 求解循环。
def cal_loop():
    loop_end = data.loop_end  # 最大循环次数
    # 架线时的状态
    # 取外过无风
    string_length = data.string_length  # 串长 单位m
    string_g = data.string_g  # 串重 单位N
    t_m = data.t_m  # 导线架设时的气温。单位°C
    t_e = data.t_e  # 架线时考虑初伸长的降温，取正值。单位°C
    alpha = data.alpha  # 导线膨胀系数  1/°C
    elastic = data.elastic  # 弹性系数 N/mm2
    area = data.area  # 导线面积 mm2
    lambda_m = data.lambda_m  # 导线比载 N/(m.mm)
    sigma_m = data.sigma_m  # 架线时，初伸长未释放前的最低点水平应力。单位N/mm2
    span_count = data.span_count  # 几个档距
    # n个档距,n-1个直线塔
    h_array = data.h_array
    l_array = data.l_array
    t_i = data.t
    lambda_i_array = data.lambda_i_array
    loop_count = 1
    sigma_0 = sigma_m * 0.2
    while True:
        b_i = 0
        # 一次增加0.1N
        sigma_0 = sigma_0 + 0.01 / data.area
        sigma_array = [sigma_0 for _ in range(span_count)]
        delta_l_i_array = []
        for i in range(span_count):
            h_i = h_array[i]
            l_i = l_array[i]
            lambda_i = lambda_i_array[i]
            t_i = t_i
            sigma_i = sigma_array[i]
            _delta_l_i = delta_li(
                h_i,
                l_i,
                lambda_i,
                alpha,
                elastic,
                t_e,
                t_i,
                sigma_i,
                lambda_m,
                t_m,
                sigma_m,
            )
            delta_l_i_array.append(_delta_l_i)
            b_i += _delta_l_i
            length_i = string_length
            g_i = string_g / data.conductor_n
            if i < span_count - 1:
                lambda_i1 = lambda_i_array[i + 1]
                h_i1 = h_array[i + 1]
                l_i1 = l_array[i + 1]
                try:
                    sigma_i1 = fun_sigma_i1(
                        area,
                        sigma_i,
                        b_i,
                        length_i,
                        g_i,
                        h_i,
                        l_i,
                        lambda_i,
                        h_i1,
                        l_i1,
                        lambda_i1,
                    )
                except ValueError:
                    break
                    pass
                sigma_array[i + 1] = sigma_i1
        if math.fabs(b_i) < data.epsilon:
            print("迭代{loop_count}次找到解。".format(loop_count=loop_count))
            print("悬垂串偏移累加bi为{b_i}".format(b_i=b_i))
            for i in range(span_count):
                print("第{i}档导线应力为{tension}".format(i=i, tension=sigma_array[i]))
            for i in range(span_count - 1):
                print(
                    "第{i}串偏移值为{bias}".format(i=i, bias=math.fsum(delta_l_i_array[0:i]))
                )
            verify(
                area,
                h_array,
                l_array,
                string_length,
                string_g / data.conductor_n,
                sigma_array,
                delta_l_i_array,
                lambda_i_array,
                t_i,
                alpha,
                elastic,
                t_e,
                lambda_m,
                t_m,
                sigma_m,
            )
            break
        loop_count += 1
        if loop_count >= loop_end:
            print("!!!未找到解。")
            print(sigma_array)
            print(b_i)
            break


# 检验等式。
def verify(
    area: float,
    h_array: [float],
    l_array: [float],
    string_length: [float],
    string_g: [float],
    sigma_array: [float],
    delta_l_i_array: [float],
    lambda_array: [float],
    t_i: float,
    alpha: float,
    elastic: float,
    t_e: float,
    lambda_m: float,
    t_m: float,
    sigma_m: float,
    epsilon: float = 1e-4,
):
    # 用新版大手册p329页(5-61)第一个公式校验
    b_i = 0
    if math.fabs((math.fsum(delta_l_i_array))) > 1e-5:
        print("偏移累加不等于0")
        return
    for i in range(len(delta_l_i_array)):
        sigma_i = sigma_array[i]
        _delta_l_i = delta_l_i_array[i]
        t_i = t_i
        # 此处用新版大手册p329页(5-58)校验偏移值。
        lambda_i = lambda_array[i]
        h_i = h_array[i]
        l_i = l_array[i]
        cal_delta_l_i = delta_li(
            h_i,
            l_i,
            lambda_i,
            alpha,
            elastic,
            t_e,
            t_i,
            sigma_i,
            lambda_m,
            t_m,
            sigma_m,
        )
        if math.fabs(cal_delta_l_i - _delta_l_i) > 1e-4:
            print("!!!偏移等式不满足。")
        if i < len(delta_l_i_array) - 1:
            sigma_i1 = sigma_array[i + 1]
            left_equ = sigma_array[i + 1]
            lambda_i1 = lambda_array[i + 1]
            h_i1 = h_array[i + 1]
            l_i1 = l_array[i + 1]
            beta_i = math.atan(h_i / l_i)
            beta_i1 = math.atan(h_i1 / l_i1)
            w_i = (
                lambda_i * l_i / 2 / math.cos(beta_i)
                + sigma_i * h_i / l_i
                + (lambda_i1 * l_i1 / 2 / math.cos(beta_i1) - sigma_i1 * h_i1 / l_i1)
            )
            b_i += delta_l_i_array[i]
            # 新版大手册p329 (5-61) 最上方公式
            right_equ = sigma_i + b_i / math.sqrt(string_length ** 2 - b_i ** 2) * (
                string_g / 2 / area + w_i  # string_g已在传入时考虑了导线分裂数
            )
            # TODO 等式允许误差是否可调？
            if math.fabs(right_equ - left_equ) > epsilon:
                print(math.fabs(right_equ - left_equ))
                print("!!!应力等式不满足")
                return
    print("等式满足。")
    print("sigma")
    print(np.array(sigma_array) * area)
    print("delta_li")
    print(delta_l_i_array)
    print("串偏移")
    b_i = 0
    b_i_list = []
    for l_i in delta_l_i_array[:-1]:
        b_i += l_i
        b_i_list.append(b_i)
    print(b_i_list)


if __name__ == "__main__":
    cal_loop()