悬链线方程和C语言实现

来源:互联网 发布:php 日程管理系统 编辑:程序博客网 时间:2024/05/16 12:36

标准悬链线方程和C语言实现

cheungmine

 

(本文可以供大学本科学生和工程技术人员学习电线力学参考,也可供游戏物理引擎开发人员参考)


标准悬链线是指悬挂在2个固定端A,B的均质柔性曲线。根据理论力学原理,其形态仅与悬挂点的高差h,链条线密度w,最低点的水平拉力T0和水平档距l有关,即 f (T0, w, l, h)。


悬链线的几何方程决定了线上每一点在悬链线坐标系中的坐标:

                           y = g(T0, w, l, h,x).


悬链线的物理方程决定了线上每一点的轴向拉力,这个拉力可以分解为水平拉力T0和竖直拉力Tvx:

                          Tvx = T(T0, w, l, h,x).


悬链线坐标系是定义为以低悬点A为原点,指向高悬点在A的水平平面的投影B'为X轴,竖直向上与荷载平行为Y轴的平面直角坐标系。低悬点必须是A。


下面是标准悬链线方程的公式的函数化表示。参考书为:


《架空送电线路的电线力学计算(第二版)》 中国电力出版社,邵天晓

 

说明:

1)算法中的 P61 (3-6)指书中61页的公式(3-6)


2)本文仅仅是书中公式的C语言翻译。并未针对实际应用进行优化,因此存在大量充分计算的参数,如Lh0,s,cosB等。


3)在工程中应该采用双曲线版本,如果对精度要求不严格,可以采用斜抛物线版本。本文未提供平抛物线版本,平抛物线误差太大,在工程中没有实际意义

 


/***************************************************************************************

* 悬链线计算用坐标系为以A(低悬点)为原点,AB'B在与A等高的水平面上的投影点)为X轴,

* AB之间的高差 h = HB-HAh>=0。以下为变量命名规则:

*

* T 线上应力,单位截面上作用的张力:N/mm2

* L 线长:m

* H 高度:m

* w 比载,电线单位长度、单位截面上承受的荷载称为比载: N/mm2.m

* l 档距,相邻两悬挂点间垂直于荷载方向的投影距离称为档距:m

* h 高差,相邻两悬挂点沿荷载方向的高差,当右悬挂点高于左悬挂点时h为正值,反之为负值:m

*

* f 电线的弧垂,指自两悬挂点连线沿荷载方向到电线轴线间的距离:m

* fM 档距中央弧垂,档内电线的最大弧垂

* fO 悬链线最低点O的弧垂

*

* fAO 低悬点AO的高差

* fBO 高悬点BO的高差

*

* A 左(低)悬挂点,应力 TA

* B 右(高)悬挂点,应力 TB

* O 最低点,仅有水平应力 T0

*

* s 斜距|AB|,两悬挂点连线的直线距离:sqrt(l*l+h*h)

* a 平距|AO|,电线最低点到左悬挂点A垂直于荷载方向的投影距离

* b 平距|BO|,电线最低点到右悬挂点B垂直于荷载方向的投影距离

* cosB ls的比:1/s

*

* 其中:

* elec 电力

* catenary 悬链线

* hyperbolic 双曲线公式,精确计算公式

* parabolic  抛物线公式,高阶近似计算公式

* sag 弧垂

* suspension 悬挂点

****************************************************************************************/

 

/******************************************************************************

*                          Catenary CommonFunctions                       *

*                                                                           *

* asinh                                                                    *

* elec_catenary_inclined_dist                                              *

*                                                                           *

*****************************************************************************/

/**

* 反双曲正弦函数

* asinh(x) = ln[x + sqrt(x^2 + 1)]

*/

double asinh(double x)

{

    returnlog(x+sqrt(x*x+1));

}

 

/**

* 斜距s

* 已知A,B悬挂点的档距l和高差hA,B间的直线距离:斜距s

*/

double elec_catenary_inclined_dist ( double l,

                                    double h

                                   )

{

    returnsqrt(l*l + h*h);

}

 

 

/******************************************************************************

*                          Catenary HyperbolicFunctions                   *

*                        An Accurate Expression ofCatenary                *

*                                                                           *

* elec_catenary_hyperbolic_equation                                        *

*                                                                           *

*                                                                           *

*****************************************************************************/

/**

* 悬挂点等高时悬链线线长双曲线公式. P61 (3-6)

*/

double elec_catenary_hyperbolic_length_equal_high (

                             double T0,

                             double w,

                             double l

                             )

{

    return(2*T0/w)*sinh(w*l/2/T0);

}

 

/**

* 悬链线线长双曲线公式. P73 (3-37)

*/

double elec_catenary_hyperbolic_length (

                             double T0,

                             double w,

                             double l,

                             double h

                             )

{

    doubleLh0 = elec_catenary_hyperbolic_length_equal_high(T0, w, l);

    returnsqrt(h*h+Lh0*Lh0);

}

 

/**

* 悬链线最低点O到悬挂点A间垂直于荷载方向的投影 a (可以<0)

* 使用双曲线方程(精确式)计算 P61 (3-6)

* 说明 b = l - a

*/

double elec_catenary_hyperbolic_lowest_proj (

                             double T0,

                             double w,

                             double l,

                             double h

                             )

{

    doubleLh0 = elec_catenary_hyperbolic_length_equal_high (T0, w, l);

    doublea = l/2 - T0/w*asinh(h/Lh0);

    doubleb = l/2 + T0/w*asinh(h/Lh0);

    returna;

}

 

/**

* 悬链线双曲线坐标方程 P62 (3-9)

* 悬链线上任何一点对A点的高差

* y = f(x), x = [0, l]

* 验证:f(0)=0, f(l)=h

*/

double elec_catenary_hyperbolic_equation (

                             double T0,

                             double w,

                             double l,

                             double h,

                             double x

                             )

{

    doublea = elec_catenary_hyperbolic_lowest_proj (T0, w, l, h);

    doubley = T0*(cosh((x-a)*w/T0) - cosh(w*a/T0))/w;

    returny;

}

 

/**

* 悬链线双曲线力学方程 P66 (3-20)

* 悬链线上任一点轴向应力

* T = t(x), x = [0, l]

* 验证 TA = t(0), TB =t(l)

*/

double elec_catenary_hyperbolic_tension (

                             double T0,

                             double w,

                             double l,

                             double h,

                             double x

                             )

{

    doublea = elec_catenary_hyperbolic_lowest_proj (T0, w, l, h);

    doubleTx = T0*cosh(w*(x-a)/T0);

    returnTx;

}


/**

* 悬挂点应力TA, TB

* 使用双曲线方程(精确式)计算 P67 (3-24)

*/

void elec_catenary_hyperbolic_suspension_tension (

                             double T0,

                             double w,

                             double l,

                             double h,

                             double *TA,

                             double *TB

                             )

{

    doubleLh0 = elec_catenary_hyperbolic_length_equal_high (T0, w, l);

   

    doubleu = w*l/T0/2;

    doublev = asinh(h/Lh0);

   

    *TA = T0*cosh(u - v);

    *TB = T0*cosh(u + v);

}

 

/**

* 悬链线悬挂点处的竖向应力分量. P72 (3-36)

* Tva, Tvb

*/

void elec_catenary_hyperbolic_vert_suspension_tension (

                             double T0,

                             double w,

                             double l,

                             double h,

                             double *Tva,

                             double *Tvb

                             )

{

    doubleu = w*l/2/T0;

    doublev = sinh(u);

    doublek = h*w/2/T0;

 

    *Tva = T0*sinh(u - asinh(k/v));

    *Tvb = T0*sinh(u + asinh(k/v));

}

 

/**

* 悬链线双曲线弧垂计算方程 P63 (3-12)

* f = d(x), x = [0, l]

* 验证:fA = d(0) = 0,fB = d(l) = 0

*/

double elec_catenary_hyperbolic_sag (

                             double T0,

                             double w,

                             double l,

                             double h,

                             double x

                             )

{

    doubleLh0 = elec_catenary_hyperbolic_length_equal_high (T0, w, l);

    doublefx = h*x/l - T0*h/w/Lh0*(sinh(w*l/2/T0) + sinh(w*(2*x-l)/2/T0)) +

                (2*T0/w*sinh(w*x/2/T0)*sinh(w*(l-x)/2/T0))*sqrt(1+h*h/(Lh0*Lh0));

    returnfx;

}

 

/**

* 悬链线双曲线电线平均高度精确式 P78 (3-44a)

* 电线平均高度是指电线各点与电线弧垂最低点间的高差沿档距的积分面积被档距除得的商,

* 即全档电线高出最低点的平均高度

*/

double elec_catenary_hyperbolic_avg_high (

                             double T0,

                             double w,

                             double l,

                             double h

                             )

{

    doubleu = sinh(w*l/2/T0);

    doublev = cosh(w*l/2/T0);

    doublehav;

 

    /**calc hav from a

     * double a =elec_catenary_hyperbolic_lowest_proj(T0, w, l, h);

     * hav =T0/w/l*( T0/w * ( sinh(w*l/T0)*cosh(w*a/T0) - cosh(w*l/T0)*sinh(w*a/T0) +sinh(w*a/T0) ) - l );

     */

   

    /**calc hav from Lh0 */

    doubleLh0 = elec_catenary_hyperbolic_length_equal_high (T0, w, l);

    hav = 2*T0*T0*u/(w*w*l)*(v*v -u*u)*sqrt(1+h*h/Lh0/Lh0) - T0/w;

 

    returnhav;

}

 

/**

* 悬链线双曲线电线平均应力精确式 P80 (3-45)

* 全档电线各点不同应力下产生的全部弹性伸长通常用一个产生全部弹性伸长的等效应力

* 或平均应力Tav来计算。该应力的大小为各点应力沿线长积分被线长除后所得的平均值

* Tav = 1/L*∫(0~L)Tx.dL∫(0~L)表示积分符号

*/

double elec_catenary_hyperbolic_avg_tension (

                             double T0,

                             double w,

                             doublel,

                             double h

                             )

{

    doubleL = elec_catenary_hyperbolic_length(T0, w, l, h);

    doubleTav = T0/2/L*(l+ (L*L+h*h)*cosh(w*l/2/T0)/sqrt(L*L-h*h));

   

    /**An Approximate Expression of Tav

     * double Tav =T0/2/L*((l*l+h*h+L*L)/l);

     */

 

    returnTav;

}

 

/**

* 悬链线双曲线任一点轴向斜率 P65 (3-19)

* 电线切线与x轴正方向间的夹角的正切

*/

double elec_catenary_hyperbolic_tangent(

                             double T0,

                             double w,

                             double l,

                             double h,

                             double x

                             )

{

    doublea = elec_catenary_hyperbolic_lowest_proj(T0, w, l, h);

    returnsinh(w*(x-a)/T0);

}

 

 

/******************************************************************************

*                          Catenary ParabolicFunctions                    *

*                      An Approximate Expression ofCatenary               *

*                                                                           *

* elec_catenary_parabolic_equation                                         *

*                                                                           *

*                                                                           *

*****************************************************************************/

/**

* 悬链线线长抛物线公式. P74 (3-39)

* 工程上采用可以得到极为精确的近似值

*/

double elec_catenary_parabolic_length (

                             double T0,

                             double w,

                             double l,

                             double h

                             )

{

    doubles = elec_catenary_inclined_dist (l, h);

    returns + (w*w*l*l*l*l)/(24*s*T0*T0);

}

 

/**

* 悬链线抛物线坐标方程 P82 (3-48)

* 悬链线上任何一点对A点的高差

* y = f(x), x = [0, l]

* 验证:f(0)=0, f(l)=h

*/

double elec_catenary_parabolic_equation (

                             double T0,

                             double w,

                             double l,

                             double h,

                             double x

                             )

{

    doublet = w/(2*T0*l/sqrt(l*l+h*h));

    doubley = (x*t + (h/l - l*t)) * x;

    returny;

}

 

/**

* 悬链线抛物线力学方程. P85 (3-60a)

* 悬链线上任一点轴向应力

* T = t(x), x = [0, l]

*/

double elec_catenary_parabolic_tension (

                             double T0,

                             double w,

                             double l,

                             double h,

                             double x

                             )

{

    doublecosB = l / sqrt(l*l+h*h);

    doubleu = h/l - w*(l-2*x)/(2*T0*cosB);

    doubleTx = T0 * sqrt(1 + u*u);

    returnTx;

}

 

/**

* 悬链线按抛物线方程近似计算悬挂点应力TA, TB. P85(3-63, 3-64)

*

 */

void elec_catenary_parabolic_suspension_tension (

                             double T0,

                             double w,

                             double l,

                             double h,

                             double *TA,

                             double *TB

                             )

{

    doublecosB = l / sqrt(l*l+h*h);

    doubleu = T0/cosB;

    doublev = w*l*l/(8*T0*cosB);

   

    *TA = u + (v-h/2)*w;

    *TB = u + (v+h/2)*w;

}

 

/**

* 悬链线按抛物线方程近似计算悬挂点应力竖向分量 P84 (3-59)(3-60)

* Tva, Tvb

*/

void elec_catenary_parabolic_vert_suspension_tension (

                             double T0,

                             double w,

                             double l,

                             double h,

                             double *Tva,

                             double *Tvb

                             )

{

    doublecosB = l / sqrt(l*l+h*h);

    *Tva = w*l/2/cosB - T0*h/l;

    *Tvb = w*l/cosB -(*Tva);   

}

 

/**

* 给定高悬点应力TB求最低点容许应力T0 P86 (3-65)

* 按抛物线方程取近似值

*/

double elec_catenary_parabolic_lowest_tension (

                             double TB,

                             double w,

                             double l,

                             double h

                             )

{

    doubles = elec_catenary_inclined_dist ( l, h );

    doublecosB = l / s;

 

    doubleu = TB-w*h/2;

 

    doubleT0 = cosB*u/2 + sqrt(u*u*cosB*cosB - w*w*l*l/2)/2;

 

    returnT0;

}

 

/**

* 悬链线最低点O到悬挂点A间垂直于荷载方向的投影 a (可以<0)

* 使用抛物线方程(近似式)计算 P62 (3-7)

* 说明 b = l - a

*/

double elec_catenary_parabolic_lowest_proj (

                             double T0,

                             double w,

                             double l,

                             double h

                             )

{

    doublecosB = l / elec_catenary_inclined_dist ( l, h );

    doublea = l/2 - T0*h*cosB/(w*l);

    returna;

}

 

/**

* 悬链线电线平均高度

* 抛物线方程 P86 (3-66)

*/

double elec_catenary_parabolic_avg_high (

                             double T0,

                             double w,

                             double l,

                             double h

                             )

{

    doublecosB = l / elec_catenary_inclined_dist ( l, h );

    doublehav = w*l*l/(24*T0*cosB) + T0*cosB*h*h/(2*w*l*l);

    returnhav;   

}

 

/**

* 悬链线电线平均应力

* 抛物线方程 P86 (3-67)

*/

double elec_catenary_parabolic_avg_tension (

                             double T0,

                             double w,

                             double l,

                             double h

                             )

{

    doublecosB = l / elec_catenary_inclined_dist ( l, h );

    doubleTav = T0/cosB + w*w*l*l/(24*T0*cosB);

    returnTav;

}

 

/**

* 悬链线抛物线任一点轴向斜率 P83 (3-51)

* 电线切线与x轴正方向间的夹角的正切

*/

double elec_catenary_parabolic_tangent(

                             double T0,

                             double w,

                             double l,

                             double h,

                             double x

                             )

{

    doublecosB = l / elec_catenary_inclined_dist ( l, h );

    returnh/l - w*(l-2*x)/(2*T0*cosB);

}

 

 

int main(int argc,CHAR* argv[])

{

    doublel = 400;     /* m */

    doubleh = 100;     /* m */

    doublew = 0.06134; /* N/(m.mm2) */

    doubleT0 = 98.1;   /* N/mm2 */

 

    // 双曲线长精确式

    doubleL = elec_catenary_hyperbolic_length(T0, w, l, h);

   

    //抛物线长近似式

    doubleLp = elec_catenary_parabolic_length(T0, w, l, h);

   

    //双曲线悬点拉力TA, TB

    doubleTA, TB;

   elec_catenary_hyperbolic_suspension_tension(T0, w, l, h, &TA, &TB);

 

    // 抛物线悬点拉力TAp, TBp

    doubleTAp, TBp;

   elec_catenary_parabolic_suspension_tension(T0, w, l, h, &TAp, &TBp);

 

    // 已知高悬点拉力反算最低点应力T0p

    doubleT0p = elec_catenary_parabolic_lowest_tension(TB, w, l, h);

 

    // 计算弧垂 sag

    doublefA = elec_catenary_hyperbolic_sag(T0, w, l, h, 0);

    doublefB = elec_catenary_hyperbolic_sag(T0, w, l, h, l);

 

    // 档距中央弧垂

    doublefM = elec_catenary_hyperbolic_sag(T0, w, l, h, l/2);

 

    // 最低点弧垂

    doublea;

    a =elec_catenary_hyperbolic_lowest_proj(T0, w, l, h);

 

    doubleap;

    ap =elec_catenary_parabolic_lowest_proj(T0, w, l, h);

    doublefO = elec_catenary_hyperbolic_sag(T0, w, l, h, a);

 

    doubley;

    y =elec_catenary_hyperbolic_equation(T0, w, l, h, 0);

    y =elec_catenary_hyperbolic_equation(T0, w, l, h, a);

    y =elec_catenary_hyperbolic_equation(T0, w, l, h, l/2);

    y =elec_catenary_hyperbolic_equation(T0, w, l, h, l);

 

    // 电线平均高度

    doublehav = elec_catenary_hyperbolic_avg_high(T0, w, l, h);

    doubleTav = elec_catenary_hyperbolic_avg_tension(T0, w, l, h);

 

    doublehavp = elec_catenary_parabolic_avg_high(T0, w, l, h);

    doubleTavp = elec_catenary_parabolic_avg_tension(T0, w, l, h);

 

 

    // 悬挂点竖向应力

    doubleTva = elec_catenary_hyperbolic_vert_tension(T0, w, l, h, 0);

    doubleTvb = elec_catenary_hyperbolic_vert_tension(T0, w, l, h, l);

    doubleTva2;

    doubleTvb2;

   elec_catenary_hyperbolic_vert_suspension_tension(T0, w, l, h, &Tva2,&Tvb2);

 

   

    doubleTvap;

    doubleTvbp;

   elec_catenary_parabolic_vert_suspension_tension(T0, w, l, h, &Tvap,&Tvbp);

 

    // 计算验证 w×L = Tva+Tvb

 

    return0;

}

 

原创粉丝点击