龙格库塔法计算固体火箭发动机内弹道
- 格式:doc
- 大小:96.00 KB
- 文档页数:8
固体火箭发动机内弹道学方丁酉固体火箭发动机内弹道学方丁酉近年来,随着航天事业的蓬勃发展,固体火箭发动机内弹道学成为越来越受关注的领域。
方丁酉,中国工程院院士、火箭技术专家,曾在这一领域做出了卓越的贡献。
本文将从简单到复杂、由浅入深地探讨固体火箭发动机内弹道学,以帮助读者全面、深入地理解这一话题。
一、固体火箭发动机内弹道学的概念固体火箭发动机内弹道学是研究固体火箭发动机内燃烧过程、燃料燃烧特性、燃烧产物排放、工作环境等方面的学科。
方丁酉指出,固体火箭发动机内弹道学是固体火箭技术中的重要基础理论,对于提高固体火箭发动机的性能、可靠性和运载能力具有重要意义。
二、固体火箭发动机内弹道学的理论基础在固体火箭发动机内弹道学的研究中,燃烧动力学是重要的理论基础之一。
燃烧动力学研究了燃料在燃烧过程中的变化规律,以及燃烧反应对火箭发动机内部气体流动和压力变化的影响。
燃烧产物排放和燃烧室内部气体动力学也是固体火箭发动机内弹道学的重要内容之一。
研究这些理论基础可以帮助我们更好地理解固体火箭发动机内部的工作原理和特性。
三、固体火箭发动机内弹道学的关键技术在固体火箭发动机内弹道学研究中,燃烧稳定性和效率是两个关键技术。
燃烧稳定性是指在固体火箭发动机工作过程中保持燃烧的稳定性,避免出现燃烧不均匀或燃烧失稳等问题。
而燃烧效率则是指在燃烧过程中尽可能提高燃料的利用率,减少燃料的浪费。
方丁酉在固体火箭发动机内弹道学的研究中,提出了一系列有效的技术方案,使固体火箭发动机在燃烧稳定性和效率方面取得了显著的进展。
四、固体火箭发动机内弹道学的应用前景固体火箭发动机内弹道学的研究成果已经在我国的航天事业中得到了广泛的应用。
在长征系列火箭、嫦娥探月工程等多个航天工程中,固体火箭发动机内弹道学的研究成果为提高火箭的性能和可靠性做出了重要贡献。
未来,随着我国航天事业的不断发展,固体火箭发动机内弹道学将继续发挥重要作用,推动我国航天事业迈向新的高度。
龙格库塔法解轨道参数引言龙格库塔法(Runge-Kutta method)是一种常用的数值计算方法,用于求解常微分方程的数值解。
在天体力学中,我们经常需要通过数值方法来计算天体的轨道参数,如轨道椭圆的长短轴、离心率、倾角等。
本文将介绍龙格库塔法在解轨道参数中的应用,并详细探讨该方法的原理和实现过程。
基本原理龙格库塔法是一种迭代求解的方法,在每个时间步长内利用当前的状态来估计下一个状态。
具体而言,龙格库塔法将微分方程的求解问题转化为一个迭代的求解问题,通过逐步迭代来逼近精确解。
在解轨道参数的问题中,我们通常需要根据已知的初始条件以及天体的质量和力学模型来求解天体的轨道参数。
常用的力学模型有开普勒模型和牛顿模型。
龙格库塔法可以根据力学模型的不同进行相应的求解。
开普勒模型下的轨道参数求解步骤一:确定初始条件在使用龙格库塔法求解轨道参数之前,我们需要确定一些初始条件。
这些初始条件包括天体的质量、位置和速度。
步骤二:选择时间步长在求解过程中,我们需要选择一个合适的时间步长。
时间步长越小,计算的精度会越高,但计算的时间会增加。
步骤三:迭代求解利用龙格库塔法进行迭代求解的具体步骤如下:1.根据当前时刻的位置和速度,计算天体在该时刻的加速度。
2.根据当前时刻的位置、速度和加速度,计算下一个时刻的位置和速度。
3.更新当前时刻的位置和速度为新的位置和速度。
4.重复上述步骤,直到达到指定的终止条件。
步骤四:计算轨道参数通过迭代求解,我们可以得到天体在不同时刻的位置和速度。
根据这些位置和速度,我们可以计算出轨道参数,如离心率、倾角、长短轴等。
常用的轨道参数计算公式如下:1.离心率:e=√1+2El2μ(GM⊕)2)2.倾角:i=arccos(ℎzℎ3.长轴:a=−μT22E4.短轴:b=a√1−e2其中,E表示能量,l表示轨道角动量,μ表示标准引力参数,G表示引力常数,M⊕表示地球的质量,ℎ和ℎz分别表示轨道角动量和轨道角动量在z轴上的分量。
内弹道是指射程较短的导弹或火箭弹在飞行过程中受到大气阻力和重力等作用的飞行轨迹。
内弹道理论研究的是导弹或火箭弹在发射后到离开大气层再进入大气层末时的飞行过程。
内弹道包括导弹或火箭弹在发射后的加速、稳定、制导、飞行以及飞行过程中的动力学性能仿真等诸多内容。
内弹道有着复杂的飞行特性和动力学方程,在实际工程中需要进行准确的计算和仿真。
内弹道的计算中,龙格库塔(Runge-Kutta)法是一种常用的数值积分方法,在求解微分方程等领域有着广泛的应用。
龙格库塔法是由数学家奥特翁格(C. W. Runge)和马丁庫塔(M. W. J. Kutta)于1900年提出的,用于求解常微分方程初值问题,其优点是精度较高,适用范围广。
在内弹道计算中,可以利用龙格库塔法对导弹或火箭弹的飞行轨迹进行数值模拟和计算,得到较为准确的飞行轨迹数据。
在实际工程中,为了方便进行内弹道的计算,可以使用Matlab等数学建模和仿真软件。
Matlab是一种常用的科学计算软件,具有强大的数值计算和仿真功能,可以用于内弹道计算中的龙格库塔法数值模拟。
在Matlab中,可以编写相应的程序,利用龙格库塔法对导弹或火箭弹的飞行过程进行仿真和计算,得到准确的飞行轨迹和动力学性能数据。
内弹道计算是导弹或火箭弹研究设计中的重要内容,龙格库塔法是一种常用的数值积分方法,Matlab是一种常用的科学计算软件,它们的应用能够有效地进行内弹道的计算和仿真,为导弹或火箭弹的研制提供重要的技术支持。
随着技术的不断发展,内弹道计算已经成为导弹或火箭弹研究设计中不可或缺的一部分。
在内弹道计算中,龙格库塔法是一种常用的数值积分方法,可以对导弹或火箭弹的飞行轨迹进行数值模拟和计算,提供准确的飞行轨迹数据。
而Matlab作为一种强大的科学计算软件,对于内弹道的计算和仿真也有着重要的应用价值。
在实际工程中,使用Matlab编写程序,利用龙格库塔法对导弹或火箭弹的飞行轨迹进行数值模拟和计算,将为导弹或火箭弹的研制提供重要的技术支持。
摘要在日常生活和工作中,我们经常遇到很多问题。
例如机械、电讯、核能、火箭、人造卫星、生物、医学及若干社会学科(如人口理论、经济预测等)的各个领域的问题,尤其是弹道轨道的定位、大型机械振动的分析、自动控制的设计、气象数值预报、按龄人口增长宏观预测问题等等。
这些问题都归结于求解常微分方程。
常微分方程的研究与应用已经深入到自然科学和社会科学的众多领域,并且成功地揭示了许多自然和社会现象的内在规律。
数值方法是求解常微分方程的一种有效的方法。
数值方法有很多种,龙格库塔法是解决初值问题的一种有效的数值方法。
我们可以通过选取不同参数从而得到不同的龙格库塔格式。
对于一类特定的初值问题,我们可以用选取恰当参数的龙格库塔格式求解这一类特定的初值问题。
每一种方法都有其自身的优点和缺点。
通过与解决这一类特定的初值问题的欧拉法的分析比较,我们可得到这两种方法的优点和缺点。
关键词:常微分方程,龙格库塔法,初值问题,导弹动力学弹道方程组IIIAbstractIn daily life and work, we often encounter many problems. Such as machinery, telecommunications, nuclear energy, rockets, satellites, biotechnology, medicine, and number of social disciplines (such as population theory, economic forecasts, etc.) in all areas, Especially ballistic orbit location, large machinery vibration analysis, control design, numerical weather forecasting, macroeconomic forecasts by age problems like population growth. These problems are due to solution of ordinary differential equations. Research and application of ordinary differential equations have been deep into many natural and social sciences fields, and successfully reveals the many natural and social phenomena of the inherent laws. Numerical solution of ordinary differential equations is an effective method. There are many numerical methods, Runge-Kutta method to solve initial value problems is an effective numerical method. We can choose different parameters resulting in different Runge-Kutta format. For a class-specific initial value problem, we can select the appropriate parameters of Runge-Kutta Schemes for this type of specific initial value. Each method has its own advantages and disadvantages. With the initial value to solve this particular type of problem analysis and comparison of the Euler method, we can obtain the advantages and disadvantages of both methods.Keywords:Ordinary Differential Equations, Runge-Kutta method, initial value problem,Ballistic missile dynamics equationsIV目录第一章前言1.1 常微分方程的应用 (1)1.2 数值方法的基本思想和内容 (1)1.3 龙格库塔法研究导弹轨迹的目的和意义 (3)第二章研究内容与方法2.1 选参构建龙格库塔格式 (4)2.2 误差分析 (7)2.3 简化导弹动力学方程组 (11)2.4 欧拉法求解导弹动力学方程组数值 (12)2.5 龙格库塔法求解导弹动力学方程组数值 (13)第三章研究结果与讨论3.1 欧拉法求解结果 (17)3.2 龙格库塔法求解结果 (17)3.3 两种方法求解结果的讨论 (17)第四章结论与展望4.1 结论 (19)4.2 改进与展望 (19)参考文献 (20)致谢 (21)I第一章前言1.1 常微分方程的应用求解常微分方程是解决现实生活中很多问题的有效的方法。
2009年10月第10期电 子 测 试ELEC TRONIC T ES TOct.2009No.10固体火箭发动机内弹道性能的仿真研究刘宝华1,杨志菊2(1辽宁葫芦岛市海军飞行学院教研部 葫芦岛 125000,2辽宁葫芦岛92941部队 葫芦岛 125001) 摘 要:研究发动机燃烧室内压强随时间变化的规律,是固体火箭发动机工作过程分析的主要方面,一般采用实验法可以获得直观可靠的数据,但实验法耗资大、周期长,不易操作。
本文根据零维内弹道数学模型,运用龙格库塔法及M A T L A B语言对某型固体导弹发动机内弹道工作过程进行数值仿真,画出燃烧室内压强随时间变化的曲线,并进一步分析得出影响发动机内弹道性能的因素,仿真结果与发动机燃烧室内工作情况相符,影响因素分析与实验相一致,为固体火箭发动机内弹道性能研究提供了参考。
关键词:发动机;内弹道;数值仿真;燃烧室中图分类号:T N957 文献标识码:BSimulate research of solide rocket eng ine inside trajecto ry characteristicLiu Baohua1,Y ang Zhiju2(1T eaching and scientific resea rch ministry,N aval Flying A czdemy,Huludao125001,China;292941U nit,PL A,H uludao125001,China)A bstract:Researching the rule of the pressure with the time changing o f the eng ine's firebox is the m ost aspect of analy zing the process o f the solid missile engine.Co mmo nly,experimenta-tion can achieve the intuitionistic and credible data,but this method consume m uch tim e and money,and not easy to manipulate.According to zero-dimensional inside trajecto ry numerical simulate model,to o ur co untry some ty pe solid missile engine,applacating Lo nger-Kuta me thod and the lang uage of m atlab launches the inside trajecto ry numerical simulate of the w o rking pro cess,and draw ed the burner pressure-time curve,passing further analy sis o uted the influence facto rs of inside trajecto ry characteristic.The result of simulatio n is consistent to the fact in en-gine's firebox and analysis of the influence facto rs is consistent to experim entation.It offer s a new idea to researching the capability o f the ballistic trajecto ry.Keywords:engine;inside trajectory;numerical sim ulate;firebox0 引 言固体火箭发动机内弹道学的核心是研究发动机燃烧室内压强随时间变化的规律,是固体火箭发动机工作过程分析的主要方面。
龙格-库塔方法3.2 Runge-Kutta法3.2.1 显式Runge-Kutta法的一般形式上节已给出与初值问题(1.2.1)等价的积分形式(3.2.1)只要对右端积分用不同的数值求积公式近似就可得到不同的求解初值问题(1.2.1)的数值方法,若用显式单步法(3.2.2)当,即数值求积用左矩形公式,它就是Euler法(3.1.2),方法只有一阶精度,若取(3.2.3)就是改进Euler法,这时数值求积公式是梯形公式的一种近似,计算时要用二个右端函数f的值,但方法是二阶精度的.若要得到更高阶的公式,则求积分时必须用更多的f值,根据数值积分公式,可将(3.2.1)右端积分表示为注意,右端f中还不能直接得到,需要像改进Euler法(3.1.11)一样,用前面已算得的f值表示为(3.2.3),一般情况可将(3.2.2)的 表示为(3.2.4)其中这里均为待定常数,公式(3.2.2),(3.2.4)称为r级的显式Runge-Kutta法,简称R-K方法.它每步计算r个f值(即),而k由前面(i-1)个已算出的表示,故公式是显式的.例i如当r=2时,公式可表示为(3.2.5) 其中.改进Euler 法(3.1.11)就是一个二级显式R-K 方法.参数取不同的值,可得到不同公式.3.2.2 二、三级显式R-K 方法对r=2的显式R-K 方法(3.2.5),要求选择参数,使公式的精度阶p 尽量高,由局部截断误差定义11122211()()[(,())(,)]n n n n n n n T y x y x h c f x y x c f x a h y b hk ++=--+++ (3.2.6) 令,对(3.2.6)式在处按Taylor 公式展开,由于将上述结果代入(3.2.6)得要使公式(3.2.5)具有的阶p=2,即,必须(3.2.7)即由此三式求的解不唯一.因r=2,由(3.2.5)式可知,于是有解(3.2.8)它表明使(3.2.5)具有二阶的方法很多,只要都可得到二阶精度R-K方法.若取,则,则得改进Euler法(3.1.11),若取,则得,此时(3.2.5)为(3.2.9)其中称为中点公式.改进的Euler法(3.1.11)及中点公式(3.2.9)是两个常用的二级R-K方法,注意二级R-K方法只能达到二阶,而不可能达到三阶.因为r=2只有4个参数,要达到p=3则在(3.2.6)的展开式中要增加3项,即增加三个方程,加上(3.2.7)的三个方程,共计六个方程求4个待定参数,验证得出是无解的.当然r=2,p=2的R-K方法(3.2.5)当取其他数时,也可得到其他公式,但系数较复杂,一般不再给出.对r=3的情形,要计算三个k值,即其中将按二元函数在处按Taylor公式展开,然后代入局部截断误差表达式,可得可得三阶方法,其系数共有8个,所应满足的方程为(3.2.10)这是8个未知数6个方程的方程组,解也是不唯一的,通常.一种常见的三级三阶R-K方法是下面的三级Kutta方法:(3.2.11)附:R-K 的三级Kutta 方法程序如下function y = DELGKT3_kuta(f, h,a,b,y0,varvec) format long; N = (b-a)/h;y = zeros(N+1,1); y(1) = y0; x = a:h:b;var = findsym(f); for i=2:N+1K1 = Funval(f,varvec,[x(i-1) y(i-1)]);K2 = Funval(f,varvec,[x(i-1)+h/2 y(i-1)+K1*h/2]); K3 = Funval(f,varvec,[x(i-1)+h y(i-1)-h*K1+K2*2*h]);y(i) = y(i-1)+h*(K1+4*K2+K3)/6; %满足c1+c2+c3=1,(1/6 4/6 1/6)endformat short; 3.2.3 四阶R-K 方法及步长的自动选择利用二元函数Taylor 展开式可以确定(3.2.4)中r=4,p=4的R-K 方法,其迭代公式为111223344()n n y y h c k c k c k c k +=++++其中1(,)n n k f x y =,2221(,(,))n n n n k f x a h y b hf x y =++,而33311322(,)n n k f x a h y b hk b hk =+++ 44411422433(,)n n k f x a h y b hk b hk b hk =++++共计13个参数待定,Taylor 展开分析局部截断误差,使得精度达到四阶,即误差为5()O h 。
---------------------------------------------------------------范文最新推荐------------------------------------------------------ 双脉冲火箭发动机设计+文献综述摘要20世纪50年代,固体火箭在技术上取得突破,并广泛应用于弹道导弹和运载火箭。
20世纪60年代起,国外在可控推力固体火箭发动机的理论和实验研究方面作了大量的工作,探索出了很多技术途径和设计方案。
固体火箭发动机具有使用安全性好、可靠性高、储存性能好、密度比冲高及勤务处理方便等优点,使其在战略、战术导弹武器领域内成为主要的动力装臵。
脉冲固体火箭发动机是兴起于20世纪末期的新型固体火箭发动机,它的出现为推进系统的发展开辟了新的发展方向,在航天和军事领域都有着潜在的应用价值。
本论文主要研究双脉冲火箭发动机总体结构、装药结构设计、燃烧室设计、喷管设计、点火装臵设计以及内弹道计算。
通过对双脉冲固体火箭发动机的设计了解了双脉冲固体火箭发动机的结构以及关键技术,为研制高性能固体火箭发动机奠定坚实的基础。
60131 / 16关键词脉冲固体火箭发动机动力装臵总体结构装药结构毕业设计说明书(论文)外文摘要TitleDouble Pulse Solid Rocket Motor DesignAbstractIn the 1950's, the solid rocket breakthrough a lot in the technical field , and the solid rocket motor are widely used in ballistic missiles and launch vehicles. In the 1960's, the foreign countries made a great mass of works in theoretical and experimental studies of controllable thrust solid rocket motor to explore a lot of technical approaches and designs. Solid rocket motor with the advantage of security, reliability, and storage performance density than the ascribed and service processing advantages of convenience, it become a major power plant in the strategic field of tactical missiles weapons. Pulse solid rocket motor is a new type of solid rocket motor which risen in the late 20th century, its---------------------------------------------------------------范文最新推荐------------------------------------------------------propulsion system development has opened a new direction of development, have potential applications in the aerospace and military fields. This paper mainly studies the dual-pulse rocket engine general structure, propellant structure design, the combustion chamber design, nozzle design, the ignition system design and interior ballistic calculation. The purpose of designing dual-pulse solid rocket motor is to understand the structure and key technology of dual-pulse solid rocket motor, and lay a solid foundation for the development of high-performance solid rocket motor.4.7燃烧室的强度校核255发动机喷管设计265.1喷管的型面设计265.2喷管壁厚273 / 165.3喷管的热防护276点火装臵的设计296.1点火药的选择296.2点火药量的计算296.3点火药盒设计307内弹道的计算317.1内弹道计算的基本方程31 7.2四阶龙格-库塔法介绍31 7.3计算步骤327.4曲线绘制33结论35---------------------------------------------------------------范文最新推荐------------------------------------------------------ 致谢36参考文献371绪论1.1课题的研究背景20世纪50年代,固体火箭在技术上取得突破,并广泛应用于弹道导弹和运载火箭。
固体火箭发动机零维两相内弹道研究陈军【摘要】为方便应用两相内弹道流动模型对火箭发动机进行性能预示以及提高性能预示精度,利用火箭喷管内的两相流动性能计算公式,建立了零维两相内弹道模型,包括零维两相内弹道微分方程和平衡压强公式,给出了模型中涉及到的两相特性参数的计算方法.相比于一维两相内弹道模型,该零维内弹道模型简单且满足必要的精度,适于工程应用的快速估算.利用该模型对某远程火箭发动机进行了内弹道计算与分析,计算结果与实验数据吻合良好,表明该两相内弹道模型可以有效地降低纯气相模型引起的理论与实际之间的模型偏差,有利于快速计算固体推进剂火箭发动机的两相内弹道参数以及提高预示精度.%To predict internal ballistic properties with two-phase flow for convenience and to improve prediction precision of internal ballistics in solid rocket motor(SRM),a model of zero-dimensional two-phase internal ballistics was built.The model includes internal ballistic differential equation and corresponding equilibrium pressure formula,and the computational methods of twophase property parameters in the model were pared with one-dimensional twophase internal ballistic model,the model is very simple and more suitable to quick computation in project applications,and more accurate.Based on the model,an internal ballistic computation and analysis for a long-range solid rocket engine were carried out.The computed results accord well with experimental data.The zero-dimensional two-phase internal ballistic model can availably reduce theoretical error caused by single-gas-phase model.By the model,two-phase internal ballisticproperties can be quickly computed,and the prediction precision of internal ballistics in SRM can be improved.【期刊名称】《弹道学报》【年(卷),期】2013(025)002【总页数】5页(P39-43)【关键词】内弹道;两相流;固体火箭发动机;固体推进剂【作者】陈军【作者单位】南京理工大学机械工程学院,南京210094【正文语种】中文【中图分类】V435随着高能推进剂在火箭发动机中的普遍应用,两相流动对发动机性能的影响越来越受到重视。
高速旋转固体火箭发动机内流场数值计算郝雯;封锋;罗盟;谢爱元【摘要】为计算高速旋转固体火箭发动机内流场特性,采用RNG k-ε湍流模型,对某端面和内孔同时燃烧管状装药固体火箭发动机燃烧室-喷管统一内流场进行了三维流场数值计算.采用用户自定义函数UDF编程给出了质量入口边界.计算结果表明,高速旋转固体火箭发动机内流场切向速度分布较复杂,在燃烧室前封头呈现出了准Rankie组合涡分布,在喷管内呈现明显的准强迫涡特性.燃气切向速度在燃烧室和喷管喉部分别达到峰值.燃烧室前封头处燃气轴向速度变化剧烈,出现了2个峰值,并在近轴线处出现了轴向速度的负值.高速旋转对前封头与喷喉结构强度与热防护带来严重影响.【期刊名称】《固体火箭技术》【年(卷),期】2013(036)004【总页数】5页(P481-485)【关键词】高速旋转;固体火箭发动机;组合涡;切向速度;轴向速度【作者】郝雯;封锋;罗盟;谢爱元【作者单位】南京理工大学机械工程学院,南京210094;南京理工大学机械工程学院,南京210094;南京理工大学机械工程学院,南京210094;南京理工大学机械工程学院,南京210094【正文语种】中文【中图分类】V4350 引言火箭增程及底排-火箭复合增程是当今炮射弹药增程的两种主要方式。
旋转稳定型弹药飞行时的转速可达到20 000 r/min或以上。
这种高速旋转的工作条件显著改变发动机的装药燃烧规律[1]及发动机内弹道参数[2],主要表现在缩短固体火箭发动机的燃烧时间、增加推力与压强、降低比冲[3-5]。
此外,高速旋转所致的强旋流动现象对发动机结构的热防护带来负面效应[6]。
文献[7]中,采用某固体火箭发动机进行高速旋转实验研究,最后拟合出了动态燃速随装药燃面相对半径和转速变化的经验式。
本文参照此文献中拟合的燃速公式,在Fluent软件基础上,使用用户自定义函数UDF,给出质量入口边界,对端面和内孔同时燃烧的管状装药固体火箭发动机进行了三维内流场数值计算。
旋转固体火箭发动机一维内弹道计算
张为华;曹泰岳;郭印诚
【期刊名称】《推进技术》
【年(卷),期】1994()3
【摘要】给出了带"锥-柱形"装药的旋转固体火箭发动机的一维内弹道计算方法,其燃速计算采用稳态分析方法,喷管流量计算采用等熵旋转喷管流的简化方法。
所得结果和旋转试车台上实测的实验结果吻合较好。
计算和实验表明,"锥-柱形"装药呈现较大的加速度效应。
文中提供的计算方法具有较高预示精度,可推广到其它药四旋转发动机的内弹道计算,对旋转发动机设计具有指导作用。
【总页数】5页(P28-32)
【关键词】固体推进剂;火箭发动机;内弹道
【作者】张为华;曹泰岳;郭印诚
【作者单位】国防科技大学航天技术系
【正文语种】中文
【中图分类】V435.11
【相关文献】
1.固体火箭发动机燃烧室结构径向振动与内弹道性能的耦合计算 [J], 李亮;鲍福廷;孟伟;和兴锁
2.串联双燃速固体火箭发动机一维内弹道计算 [J], 鲍福廷;赵飞;李进贤
3.锥柱形装药旋转火箭发动机的零维内弹道计算 [J], 曹泰岳;吴军
4.固体火箭发动机点火过程内弹道计算 [J], 李海涛;鲍福廷
5.固体火箭发动机通用药柱的一维准定常内弹道计算 [J], 田维平;张扬中;顾雪林因版权原因,仅展示原文概要,查看原文内容请购买。
一、引言固体火箭发动机内弹道学方丁酉,是一门深奥而又具有重要意义的科学。
固体火箭发动机是航天领域中不可或缺的一部分,其内弹道学方丁酉的研究对于提高火箭发动机的性能和可靠性具有重要意义。
在本文中,将深入探讨固体火箭发动机内弹道学方丁酉的重要性、原理和应用。
二、固体火箭发动机内弹道学理论的重要性固体火箭发动机内弹道学方丁酉,是对固体火箭发动机内部燃烧过程的运动状态、燃烧特性、瞬态过程等进行研究的学科。
通过对固体火箭发动机内弹道学方丁酉的研究,可以深入了解燃烧室内各种参数的变化规律,确定燃烧室内流场的性质,以及改善发动机的设计和性能等方面提供重要参考。
固体火箭发动机内弹道学方丁酉的研究对提高火箭发动机的推进效率、降低燃料消耗、提高火箭发射成功率等具有重要意义。
通过对固体火箭发动机内弹道学方丁酉的研究,可以进一步提高火箭的响应速度和控制精度,提高火箭的战场生存能力和作战效能。
固体火箭发动机内弹道学方丁酉的研究内容主要包括燃烧室内燃气运动的多相流动、燃烧过程的热力学性质、燃烧产物的组成及其对流场和燃烧室壁面的影响等。
通过对这些内容的研究,可以从微观和宏观两个层面全面了解固体火箭发动机内弹道学方丁酉的运行机理,进而进行针对性的改进和应用。
固体火箭发动机内弹道学方丁酉的原理在航天领域中有着广泛的应用。
在燃烧室内燃气运动的多相流动方面,可以通过优化燃烧室的结构和燃烧室内部的燃烧过程,提高火箭的推进效率;在燃烧产物的组成及其对流场和燃烧室壁面的影响方面,可以减少火箭发动机的燃料消耗,提高火箭的响应速度和控制精度等。
四、固体火箭发动机内弹道学方丁酉的个人观点对于固体火箭发动机内弹道学方丁酉这一学科,我个人认为其研究的深入和广泛应用,对于航天领域的发展和火箭技术的提升有着极其重要的作用。
只有通过对固体火箭发动机内弹道学方丁酉的深入研究,才能够更好地理解火箭发动机内部燃烧过程的运动状态、燃烧特性、瞬态过程等,从而更好地改进火箭的设计和性能,提高火箭的推进效率、降低燃料消耗、提高火箭发射成功率等。
!计?算?星?型í装痢?药?的?几?何?尺?寸?!参?数簓符?号?说μ明÷!-----------------------------------------------------------------------------------!n表括?示?星?角?数簓,num表括?示?将?推?进?剂á沿?肉╝厚?方?向ò分?为a几?等台?分?,m表括?示?选?择?压1力 值μ大洙?小?!d表括?示?外猘径?,len表括?示?长¤度è,thet表括?示?星?边?夹D角?,epsilon表括?示?角?度è系μ数簓!r表括?示?过y度è圆2弧?半?径?,r1表括?示?星?角?圆2弧?半?径?,R0表括?示?通 ?用?气?体?常£数簓!l表括?示?药?柱ù的?特?征÷长¤度è,y0表括?示?初?始?特?征÷参?数簓,y1表括?示?燃?尽?特?征÷参?数簓!I0表括?示?总哩?冲?,F表括?示?推?力 ,Poc表括?示?燃?烧?室酣?的?工¤作痢?压1力!Isp表括?示?比括?冲?,density_p表括?示?密ü度è,k表括?示?比括?热è?比括?rspeed表括?示?燃?速ù,?pn表括?示?压1力 指?数簓!mpeff表括?示?有瓺效§装痢?药?量?,Cf表括?示?推?力 系μ数簓,Ctz表括?示?特?征÷速ù度è,At表括?示?喉í部?面?积y,a表括?示?燃?速ù系μ数簓!S表括?示?平?均ù燃?烧?面?积y,e1表括?示?平?均ù肉╝厚?,epsilon1表括?示?减?面?比括?epsilon2表括?示?增?面?比括?!foresmax表括?示?前°段?最?大洙?相à对?周ü边?长¤,backsmax表括?示?后ó段?最?大洙?相à对?周ü边?长¤!smin表括?示?最?小?相à对?周ü边?长¤,thet1表括?示?周ü边?长¤取?得?最?小?值μ时骸?的?星?边?夹D角?!Ap表括?示?初?始?通 ?气?面?积y,J表括?示?初?始?通 ?气?参?量?,eta表括?示?装痢?填?系μ数簓,Af表括?示?剩骸?药?面?积y,etaf表括?示?剩骸?药?系μ数簓!de表括?示?每?一?份 肉╝厚?的?长¤度è,Sa表括?示?燃?烧?面?积y数簓组哩?Apa表括?示?通 ?气?面?积y数簓组哩?!------------------------------------------------------------------------------------program mainimplicit nonereal(kind=8),parameter::Pi=3.14integer::n,num,i,mreal(kind=8)::d,len,thet,epsilon,r,r1real(kind=8)::l,y0,y1real(kind=8)::I0,F,Poc,Pe,R0real(kind=8)::Isp,density_p,k,rspeed,Pnreal(kind=8)::mpeff,Cf,Ctz,At,areal(kind=8)::S,e1,epsilon1,epsilon2real(kind=8)::foresmax,backsmax,smin,thet1real(kind=8)::Ap,J,eta,Af,etafreal(kind=8)::error,dereal(kind=8),allocatable::Sa(:),Apa(:)!读á入?所ù需è要癮所ù用?参?数簓值μopen(3,file="design_parameter.dat")read(3,*)I0!总哩?冲?read(3,*)F!推?力read(3,*)Poc!燃?烧?室酣?压1力read(3,*)Isp!比括?冲?read(3,*)density_p!推?进?剂á的?密ü度èread(3,*)k!比括?热è?比括?read(3,*)a!燃?速ù系μ数簓read(3,*)Ctz!读á取?特?征÷速ù度èCtzread(3,*)Pn!压1力 指?数簓read(3,*)Pe!读á取?喷?管ü出?口ú的?压1力 Peread(3,*)R0!读á取?通 ?用?气?体?常£数簓R0read(3,*)d,r!读á取?装痢?药?直±径?d和í过y度è圆2弧?半?径?rread(3,*)epsilon2!读á取?增?面?比括╡psilon2read(3,*)r1!读á取?星?角?圆2弧?半?径?r1read(3,*)thet1!试?取?周ü边?长¤取?得?最?小?值μ时骸?的?星?边?夹D角?thet1read(3,*)thet!试?取?初?始?时骸?的?星?边?夹D角?thetread(3,*)n!读á取?星?角?数簓n(辍?,4,5,6,7,8)?read(3,*)epsilon!试?取?角?度è系μ数簓epsilonread(3,*)num!读á取?等台?分?肉╝厚?的?等台?分?数簓numread(3,*)m!m=1表括?示?最?小?压1力 给?定¨,m=2表括?示?最?大洙?压1力 给?定¨,m=0表括?示?平?均ù压1力 给?定¨close(3)allocate(Sa(0:num),Apa(0:num))!根ù据Y规?定¨的?总哩?冲?计?算?有瓺效§装痢?药?量?mpeffmpeff=1.02*I0/Isp!计?算?推?力 系μ数簓和í喉í部?面?积yCf=sqrt(k)*(2/(k+1))**((k+1)/(k-1)/2)*sqrt(2*k*(1-(Pe/Poc)**((k-1)/k))/(k-1))At=F/Cf/(Poc*101325.0)!计?算?平?均ù肉╝厚?e1和í计?算?特?征÷长¤度è和í平?均ù燃?烧?面?积yS!若?给?定¨的?是?最?大洙?压1力 则ò需è根ù据Y增?面?比括?计?算?最?小?燃?面?再ù计?算?平?均ù燃?面?值μ!y1一?般?取?值μ在ú1附?近ü,(0.8-1.2)if(m==0)then!平?均ù压1力 给?定¨S=At*(Poc*101325.0)**(1-Pn)/(Ctz*density_p*a)e1=mpeff/(density_p*S)else if(m==1)then!最?小?压1力 给?定¨S=At*(Poc*101325.0)**(1-Pn)/(Ctz*density_p*a)S=(epsilon2+1.0)*S/2.0e1=mpeff/(density_p*S)else!最?大洙?压1力 给?定¨S=At*(Poc*101325.0)**(1-Pn)/(Ctz*density_p*a)S=(1.0/epsilon2+1.0)*S/2.0e1=mpeff/(density_p*S)end ifl=d/2-e1-ry0=r/ly1=(e1+r)/l!根ù据Yepsilon2选?择?角?度è系μ数簓epsilon!计?算?取?得?最?小?周ü长¤时骸?的?thet1do while(.true.)error=thet1/2+cotan(thet1/2)-Pi/n-Pi/2if(error>=0.0001)thenthet1=thet1+0.00001else if(error<=-0.0001)thenthet1=thet1-0.00001elseexitend ifend do!计?算?角?度è系μ数簓epsilon的?准?确ā?值μdo while(.true.)backsmax=2*n*((1-epsilon)*Pi/n+y1*(Pi/n+asin(sin(epsilon*Pi/n)/y1))) smin=2*n*(sin(epsilon*Pi/n)/sin(thet1/2)+(1-epsilon)*Pi/n)if(backSmax/Smin>epsilon2)thenepsilon=epsilon+0.001elseexitend ifend do!根ù据Y减?面?比括╡psilon1,?计?算?thetepsilon1=1/epsilon2foresmax=smin/epsilon1do while(.true.)error=2*n*(sin(epsilon*Pi/n)/sin(thet/2)+(1-epsilon)*Pi/n+&(r1+r)*(Pi/n+Pi/2-thet/2-cos(thet/2)/sin(thet/2))/l)-foresmaxif(error>=0.0001)thenthet=thet+0.00001else if(error<=-0.0001)thenthet=thet-0.00001elseexitend do!计?算?初?始?通 ?气?面?积y,?通 ?气?参?量?,?装痢?填?系μ数簓,药?柱ù长¤度èAp=(n*((1-epsilon)*Pi/n+sin(epsilon*Pi/n)*(cos(epsilon*Pi/n)-&sin(epsilon*Pi/n)*cotan(thet/2)))+2*n*r*(sin(epsilon*Pi/n)/sin(thet/2)+&(1-epsilon)*Pi/n)/l+n*r**2*(Pi/n+Pi/2-thet/2-cotan(thet/2))/l**2+&n*r1**2*(thet/2+cotan(thet/2)-Pi/2)/l**2)*l**2J=At/Apeta=4*(Pi*d**2/4-Ap)/(Pi*d**2)Af=(epsilon*Pi*(1+y1)**2-n*(sin(epsilon*Pi/n)*(sqrt(y1**2-sin(epsilon*Pi/n)**2)+&cos(epsilon*Pi/n)))-n*y1**2*(epsilon*Pi/n+asin(sin(epsilon*Pi/n)/y1)))*l**2etaf=4*Af/(pi*d**2)len=mpeff/density_p/(Pi*d**2/4-Ap-Af)!计?算?推?进?剂á的?燃?面?变?化ˉ规?律é并¢输?出?结á果?de=e1/numopen(10,file="export_burnS.dat")open(20,file="export_Ap.dat")do i=0,num,1if((i*de)<=r1)thenSa(i)=2*n*(sin(epsilon*Pi/n)/sin(thet/2)+(1-epsilon)*Pi/n+&(r1+r)*(Pi/n+Pi/2-thet/2-cotan(thet/2))/l-(r1-i*de)*Pi/n/l)*l*len Apa(i)=(n*((1-epsilon)*Pi/n+sin(epsilon*Pi/n)*(cos(epsilon*Pi/n)-sin(epsilon*Pi/n)* cotan(thet/2)))+&2*n*(r+i*de)*(sin(epsilon*Pi/n)/sin(thet/2)+(1-epsilon)*Pi/n)/l+&n*(r+i*de)**2*(Pi/n+Pi/2-thet/2-cotan(thet/2))/l**2+n*(r1-i*de)**2*(thet/2+cotan(thet/2) -Pi/2)/l**2)*l**2else if((i*de)>r1.and.(i*de)<=(l*sin(epsilon*Pi/n)/cos(thet/2)-r))thenSa(i)=2*n*(sin(epsilon*Pi/n)/sin(thet/2)+(1-epsilon)*Pi/n+&(i*de+r)*(Pi/n+Pi/2-thet/2-cotan(thet/2))/l)*l*lenApa(i)=(n*((1-epsilon)*Pi/n+sin(epsilon*Pi/n)*(cos(epsilon*Pi/n)-sin(epsilon*Pi/n)* cotan(thet/2)))+&2*n*(r+i*de)*(sin(epsilon*Pi/n)/sin(thet/2)+(1-epsilon)*Pi/n)/l+&n*(r+i*de)**2*(Pi/n+Pi/2-thet/2-cotan(thet/2))/l**2)*l**2else if((i*de)>(l*sin(epsilon*Pi/n)/cos(thet/2)-r).and.(i*de)<=e1)thenSa(i)=2*n*((1-epsilon)*Pi/n+(r+i*de)*(Pi/n+asin(l*sin(epsilon*Pi/n)/(i*de+r)))/l)*l *lenApa(i)=n*((1-epsilon)*Pi*(1+(r+i*de)/l)**2/n+sin(epsilon*Pi/n)*(sqrt((r+i*de)**2/l* *2-&sin(epsilon*Pi/n)**2)+cos(epsilon*Pi/n))+(r+i*de)**2*(epsilon*Pi/n+&asin(l*sin(epsilon*Pi/n)/(r+i*de)))/l**2)*l**2write(10,"(f10.5,2X,f15.5)")i*de,Sa(i)write(20,"(f10.5,2X,f15.5)")i*de,Apa(i)end doclose(10)close(20)!调獭?用?子哩?程ì序ò计?算?装痢?药?内ú弹獭?道台?曲ú线?callinternal_ballistics0(d,len,e1,n,thet,epsilon,r,r1,l,a,pn,Ctz,At,k,density_p,mpeff,Isp) !输?出?星?型í装痢?药?的?几?何?参?数簓open(30,file="export_star_geometry.dat")write(30,"(A6,f11.5,A2)")"mpeff=",mpeff,"Kg"write(30,"(A2,f8.5,A2)")"D=",d,"m"write(30,"(A2,f8.5,A2)")"L=",len,"m"write(30,"(A2,I2)")"n=",nwrite(30,"(A3,f8.5)")"θè=",thetwrite(30,"(A3,f7.5)")"ε?=",epsilonwrite(30,"(A2,f7.5,A2)")"r=",r,"m"write(30,"(A3,f7.5,A2)")"r1=",r1,"m"write(30,"(A3,f8.5,A2)")"e1=",e1,"m"write(30,"(A3,f10.5,A2)")"Ap=",Ap,"㎡O"write(30,"(A2,f10.5)")"J=",Jwrite(30,"(A3,f7.5)")"η?=",etawrite(30,"(A3,f10.5,A2)")"Af=",Af,"㎡O"write(30,"(A4,f10.5)")"η?f=",etafwrite(30,"(A2,f10.5)")"l=",lwrite(30,"(A3,f10.5)")"y0=",y0write(30,"(A3,f10.5)")"y1=",y1close(30)stopend!该?子哩?程ì序ò用?于 ?计?算?零?维?变?截?面?燃?烧?装痢?药?的?内ú弹獭?道台?!利?用?龙ⅷ?格?-库a塔t法ぁ?计?算?内ú弹獭?道台?曲ú线?!可é用?于 ?计?算?侵?蚀骸?燃?烧?效§应畖下?的?内ú弹獭?道台?曲ú线?!--------------------------------------------------------------------------!d外猘径?,len长¤度è,e1平?均ù肉╝厚?,n星?角?数簓,thet星?边?夹D角?,epsilon角?度è系μ数簓!r过y度è圆2弧?半?径?,r1星?角?圆2弧?半?径?和íl特?征÷尺?寸?,key表括?示?是?否?考?虑?侵?蚀骸?燃?烧?!density_p推?进?剂á的?密ü度è,k比括?热è?比括?ga系μ数簓,C特?征÷速ù度è,mpeff药?柱ù质ê量?!a速ù度è系μ数簓,pn压1强?指?数簓,At喉í部?面?积y,Pc燃?烧?室酣?的?设Θ?计?压1力 ,Isp表括?示?理え?论?比括?冲?!dt时骸?间?步?长¤,e燃?层?厚?度è,time时骸?间?,Sa同?一?时骸?刻ì轴á向ò的?各÷节ú点?的?燃?面?,Apa通 ?气?面?积y!Poc各÷节ú点?的?压1力 ,Pe喷?管ü出?口ú压1力 ,P_I压1力 冲?量?,ep侵?蚀骸?比括?!P_av平?均ù压1力 ,F_av平?均ù推?力 ,I0总哩?冲?,Im重?量?比括?冲?Iv体?积y比括?冲?,Cf 推?力 系μ数簓!--------------------------------------------------------------------------subroutineinternal_ballistics0(d,len,e1,n,thet,epsilon,r,r1,l,a,pn,Ctz,At,k,density_p,mpeff,Isp) implicit nonereal,parameter::Pi=3.14integer::n,keyreal(kind=8),intent(in)::d,len,e1,thet,epsilon,r,r1,lreal(kind=8),intent(in)::density_p,k,Ctz,mpeff,Ispreal(kind=8)::a,pn,At,dt,e,time,Sa,Apa,Poc,Pe,P_I,Freal(kind=8)::ep,P_av,F_av,I0,Im,Iv,Cf,gareal(kind=8)::f1,f2,f3,f4Pe=101325.0P_I=0.0!药?柱ù的?能ü量?特?性?参?量?读á入?open(3,file="inter_parameter.dat")read(3,*)keyread(3,*)dtread(3,*)Pocclose(3)!计?算?点?火e压1强?时骸?的?推?力 值μPoc=Poc*101325.0ga=sqrt(k)*(2/(k+1))**((k+1)/(k-1)/2)Cf=ga*sqrt(2*k*(1-(Pe/Poc)**((k-1)/k))/(k-1))F=Cf*Poc*Ate=0.0time=0.0open(10,file="export_Poc.dat")open(20,file="export_F.dat")write(10,"(f8.5,2X,f10.5)")time,Poc/101325.0write(20,"(f8.5,2X,f15.5)")time,Fdo while(.true.)!计?算?燃?面?和í通 ?气?面?积yif(e<=r1)thenSa=2*n*(sin(epsilon*Pi/n)/sin(thet/2)+(1-epsilon)*Pi/n+&(r1+r)*(Pi/n+Pi/2-thet/2-cotan(thet/2))/l-(r1-e)*Pi/n/l)*l*len Apa=(n*((1-epsilon)*Pi/n+sin(epsilon*Pi/n)*(cos(epsilon*Pi/n)-&sin(epsilon*Pi/n)*cotan(thet/2)))+2*n*(r+e)*(sin(epsilon*Pi/n)/&sin(thet/2)+(1-epsilon)*Pi/n)/l+n*(r+e)**2*(Pi/n+Pi/2-thet/2-&cotan(thet/2))/l**2+n*(r1-e)**2*(thet/2+cotan(thet/2)-Pi/2)/l**2)*l**2else if(e>r1.and.e<=(l*sin(epsilon*Pi/n)/cos(thet/2)-r))thenSa=2*n*(sin(epsilon*Pi/n)/sin(thet/2)+(1-epsilon)*Pi/n+&(e+r)*(Pi/n+Pi/2-thet/2-cotan(thet/2))/l)*l*lenApa=(n*((1-epsilon)*Pi/n+sin(epsilon*Pi/n)*(cos(epsilon*Pi/n)-&sin(epsilon*Pi/n)*cotan(thet/2)))+2*n*(r+e)*(sin(epsilon*Pi/n)/&sin(thet/2)+(1-epsilon)*Pi/n)/l+n*(r+e)**2*(Pi/n+Pi/2-thet/2-cotan(thet/2))/l**2)*l**2 else if(e>(l*sin(epsilon*Pi/n)/cos(thet/2)-r).and.e<=e1)thenSa=2*n*((1-epsilon)*Pi/n+(r+e)*(Pi/n+asin(l*sin(epsilon*Pi/n)/(e+r)))/l)*l*lenApa=n*((1-epsilon)*Pi*(1+(r+e)/l)**2/n+sin(epsilon*Pi/n)*(sqrt((r+e)**2/l**2-& sin(epsilon*Pi/n)**2)+cos(epsilon*Pi/n))+(r+e)**2*&(epsilon*Pi/n+asin(l*sin(epsilon*Pi/n)/(r+e)))/l**2)*l**2elseSa=0.0end if!计?算?侵?蚀骸?比括?if(key==1)thenif(Sa/Apa<=72.9)thenep=1.0elseep=1.3128-1.3249e-2*Sa/Apa+1.5527e-4*(Sa/Apa)**2-4.3868e-7*(Sa/Apa)**3 end ifelseep=1.0end if!采é用?4阶×的?龙ⅷ?格?库a塔t法ぁ?计?算?f1=ga**2*Ctz**2*(density_p*Sa*ep*a*Poc**pn-Poc*At/Ctz)/(Apa*len)f2=ga**2*Ctz**2*(density_p*Sa*ep*a*(Poc+dt*f1/2)**pn-(Poc+dt*f1/2)*At/Ctz)/(Apa*len )f3=ga**2*Ctz**2*(density_p*Sa*ep*a*(Poc+dt*f2/2)**pn-(Poc+dt*f2/2)*At/Ctz)/(Apa*len )f4=ga**2*Ctz**2*(density_p*Sa*ep*a*(Poc+dt*f3)**pn-(Poc+dt*f3)*At/Ctz)/(Apa*len) Poc=Poc+dt*(f1+2.0*f2+2.0*f3+f4)/6.0Cf=ga*sqrt(2*k*(1-(Pe/Poc)**((k-1)/k))/(k-1))F=Cf*Poc*Attime=time+dte=e+dt*ep*a*Poc**pnif(e>e1.and.Poc<101325.0)exitwrite(10,"(f8.5,2X,f10.5)")time,Poc/101325.0write(20,"(f8.5,2X,f15.5)")time,FP_I=P_I+Poc*dtend doclose(10)close(20)!循-环·结á束?P_av=P_I/timeCf=ga*sqrt(2*k*(1-(Pe/P_av)**((k-1)/k))/(k-1))F_av=Cf*P_av*AtI0=F_av*timeIm=I0/mpeff/9.81Iv=I0/(Apa*len)!输?出?固ì体?火e箭y发ぁ?动ˉ机ú的?工¤作痢?参?数簓值μopen(30,file="export_propellant.dat")write(30,"('平?均ù压1力 p_av=',f9.5)")p_av/101325.0write(30,"('平?均ù推?力 F_av=',f10.5,1X,'KN')")F_av/1000.0 write(30,"('总哩?冲?I0=',f10.5,1X,'KN.s')")I0/1000.0write(30,"('重?量?比括?冲?Im=',f10.5,1X,'s')")Imwrite(30,"('体?积y比括?冲?Iv=',f10.5,1X,'KN.s/m3')")Iv/1000.0 write(30,"('工¤作痢?时骸?间?time=',f8.5,1X,'s')")timeclose(30)returnend subroutine。
!计?算?星?型³装痢?药?的?几?何?尺?寸?!参?数簓符?号?说μ明¶!-----------------------------------------------------------------------------------!n表括?示?星?角?数簓,num表括?示?将?推?进?剂®沿?肉╝厚?方?向´分?为a几?等台?分?,m表括?示?选?择?压1力 值μ大洙?小?!d表括?示?外猘径?,len表括?示?长¤度¯,thet表括?示?星?边?夹D角?,epsilon表括?示?角?度¯系μ数簓!r表括?示?过y度¯圆2弧?半?径?,r1表括?示?星?角?圆2弧?半?径?,R0表括?示?通 ?用?气?体?常£数簓!l表括?示?药?柱·的?特?征¶长¤度¯,y0表括?示?初?始?特?征¶参?数簓,y1表括?示?燃?尽?特?征¶参?数簓!I0表括?示?总哩?冲?,F表括?示?推?力 ,Poc表括?示?燃?烧?室酣?的?工¤作痢?压1力!Isp表括?示?比括?冲?,density_p表括?示?密¹度¯,k表括?示?比括?热¯?比括?rspeed表括?示?燃?速·,?pn表括?示?压1力 指?数簓!mpeff表括?示?有瓺效§装痢?药?量?,Cf表括?示?推?力 系μ数簓,Ctz表括?示?特?征¶速·度¯,At表括?示?喉³部?面?积y,a表括?示?燃?速·系μ数簓!S表括?示?平?均·燃?烧?面?积y,e1表括?示?平?均·肉╝厚?,epsilon1表括?示?减?面?比括?epsilon2表括?示?增?面?比括?!foresmax表括?示?前©段?最?大洙?相对?周¹边?长¤,backsmax表括?示?后µ段?最?大洙?相对?周¹边?长¤!smin表括?示?最?小?相对?周¹边?长¤,thet1表括?示?周¹边?长¤取?得?最?小?值μ时骸?的?星?边?夹D角?!Ap表括?示?初?始?通 ?气?面?积y,J表括?示?初?始?通 ?气?参?量?,eta表括?示?装痢?填?系μ数簓,Af表括?示?剩骸?药?面?积y,etaf表括?示?剩骸?药?系μ数簓!de表括?示?每?一?份 肉╝厚?的?长¤度¯,Sa表括?示?燃?烧?面?积y数簓组哩?Apa表括?示?通 ?气?面?积y数簓组哩?!------------------------------------------------------------------------------------program mainimplicit nonereal(kind=8),parameter :: Pi=3.14integer :: n,num,i,mreal(kind=8) :: d,len,thet,epsilon,r,r1real(kind=8) :: l,y0,y1real(kind=8) :: I0,F,Poc,Pe,R0real(kind=8) :: Isp,density_p,k,rspeed,Pnreal(kind=8) :: mpeff,Cf,Ctz,At,areal(kind=8) :: S,e1,epsilon1,epsilon2real(kind=8) :: foresmax,backsmax,smin,thet1real(kind=8) :: Ap,J,eta,Af,etafreal(kind=8) :: error,dereal(kind=8),allocatable :: Sa(:),Apa(:)!读®入?所·需¯要癮所·用?参?数簓值μopen(3,file="design_parameter.dat")read(3,*) I0 !总哩?冲?read(3,*) F !推?力read(3,*) Poc !燃?烧?室酣?压1力read(3,*) Isp !比括?冲?read(3,*) density_p !推?进?剂®的?密¹度¯read(3,*) k !比括?热¯?比括?read(3,*) a !燃?速·系μ数簓read(3,*) Ctz !读®取?特?征¶速·度¯Ctzread(3,*) Pn !压1力 指?数簓read(3,*) Pe !读®取?喷?管¹出?口¸的?压1力 Peread(3,*) R0 !读®取?通 ?用?气?体?常£数簓R0read(3,*) d,r !读®取?装痢?药?直ª径?d和³过y度¯圆2弧?半?径?rread(3,*) epsilon2 !读®取?增?面?比括╡psilon2read(3,*) r1 !读®取?星?角?圆2弧?半?径?r1read(3,*) thet1 !试?取?周¹边?长¤取?得?最?小?值μ时骸?的?星?边?夹D角?thet1read(3,*) thet !试?取?初?始?时骸?的?星?边?夹D角?thetread(3,*) n !读®取?星?角?数簓n(辍?,4,5,6,7,8)?read(3,*) epsilon !试?取?角?度¯系μ数簓epsilonread(3,*) num !读®取?等台?分?肉╝厚?的?等台?分?数簓numread(3,*) m !m=1表括?示?最?小?压1力 给?定¨,m=2表括?示?最?大洙?压1力 给?定¨,m=0表括?示?平?均·压1力 给?定¨close(3)allocate(Sa(0:num),Apa(0:num))!根·据Y规?定¨的?总哩?冲?计?算?有瓺效§装痢?药?量?mpeffmpeff=1.02*I0/Isp!计?算?推?力 系μ数簓和³喉³部?面?积yCf=sqrt(k)*(2/(k+1))**((k+1)/(k-1)/2)*sqrt(2*k*(1-(Pe/Poc)**((k-1)/k))/(k-1))At=F/Cf/(Poc*101325.0)!计?算?平?均·肉╝厚?e1和³计?算?特?征¶长¤度¯和³平?均·燃?烧?面?积yS!若?给?定¨的?是?最?大洙?压1力 则´需¯根·据Y增?面?比括?计?算?最?小?燃?面?再·计?算?平?均·燃?面?值μ!y1一?般?取?值μ在¸1附?近¹,(0.8-1.2)if(m==0) then!平?均·压1力 给?定¨S=At*(Poc*101325.0)**(1-Pn)/(Ctz*density_p*a)e1=mpeff/(density_p*S)else if(m==1) then!最?小?压1力 给?定¨S=At*(Poc*101325.0)**(1-Pn)/(Ctz*density_p*a)S=(epsilon2+1.0)*S/2.0e1=mpeff/(density_p*S)else!最?大洙?压1力 给?定¨S=At*(Poc*101325.0)**(1-Pn)/(Ctz*density_p*a)S=(1.0/epsilon2+1.0)*S/2.0e1=mpeff/(density_p*S)end ifl=d/2-e1-ry0=r/ly1=(e1+r)/l!根·据Yepsilon2选?择?角?度¯系μ数簓epsilon!计?算?取?得?最?小?周¹长¤时骸?的?thet1do while(.true.)error=thet1/2+cotan(thet1/2)-Pi/n-Pi/2if(error>=0.0001) thenthet1=thet1+0.00001else if(error<=-0.0001) thenthet1=thet1-0.00001elseexitend ifend do!计?算?角?度¯系μ数簓epsilon的?准?确º?值μdo while(.true.)backsmax=2*n*((1-epsilon)*Pi/n+y1*(Pi/n+asin(sin(epsilon*Pi/n)/y1))) smin=2*n*(sin(epsilon*Pi/n)/sin(thet1/2)+(1-epsilon)*Pi/n) if(backSmax/Smin>epsilon2) thenepsilon=epsilon+0.001elseexitend ifend do!根·据Y减?面?比括╡psilon1,?计?算?thetepsilon1=1/epsilon2foresmax=smin/epsilon1do while(.true.)error=2*n*(sin(epsilon*Pi/n)/sin(thet/2)+(1-epsilon)*Pi/n+&(r1+r)*(Pi/n+Pi/2-thet/2-cos(thet/2)/sin(thet/2))/l)-foresmax if(error>=0.0001) thenthet=thet+0.00001else if(error<=-0.0001) thenthet=thet-0.00001elseexitend do!计?算?初?始?通 ?气?面?积y,?通 ?气?参?量?,?装痢?填?系μ数簓,药?柱·长¤度¯Ap=(n*((1-epsilon)*Pi/n+sin(epsilon*Pi/n)*(cos(epsilon*Pi/n)-&sin(epsilon*Pi/n)*cotan(thet/2)))+2*n*r*(sin(epsilon*Pi/n)/sin(thet/2)+&(1-epsilon)*Pi/n)/l+n*r**2*(Pi/n+Pi/2-thet/2-cotan(thet/2))/l**2+&n*r1**2*(thet/2+cotan(thet/2)-Pi/2)/l**2)*l**2J=At/Apeta=4*(Pi*d**2/4-Ap)/(Pi*d**2)Af=(epsilon*Pi*(1+y1)**2-n*(sin(epsilon*Pi/n)*(sqrt(y1**2-sin(epsilon*Pi/n)**2)+& cos(epsilon*Pi/n)))-n*y1**2*(epsilon*Pi/n+asin(sin(epsilon*Pi/n)/y1)))*l**2etaf=4*Af/(pi*d**2)len=mpeff/density_p/(Pi*d**2/4-Ap-Af)!计?算?推?进?剂®的?燃?面?变?化ˉ规?律°并¢输?出?结®果?de=e1/numopen(10,file="export_burnS.dat")open(20,file="export_Ap.dat")do i=0,num,1if((i*de)<=r1) thenSa(i)=2*n*(sin(epsilon*Pi/n)/sin(thet/2)+(1-epsilon)*Pi/n+&(r1+r)*(Pi/n+Pi/2-thet/2-cotan(thet/2))/l-(r1-i*de)*Pi/n/l)*l*len Apa(i)=(n*((1-epsilon)*Pi/n+sin(epsilon*Pi/n)*(cos(epsilon*Pi/n)-sin(epsilon*Pi/n)* cotan(thet/2)))+&2*n*(r+i*de)*(sin(epsilon*Pi/n)/sin(thet/2)+(1-epsilon)*Pi/n)/l+&n*(r+i*de)**2*(Pi/n+Pi/2-thet/2-cotan(thet/2))/l**2+n*(r1-i*de)**2*(thet/2+cotan(thet/2 )-Pi/2)/l**2)*l**2else if((i*de)>r1.and.(i*de)<=(l*sin(epsilon*Pi/n)/cos(thet/2)-r)) thenSa(i)=2*n*(sin(epsilon*Pi/n)/sin(thet/2)+(1-epsilon)*Pi/n+&(i*de+r)*(Pi/n+Pi/2-thet/2-cotan(thet/2))/l)*l*lenApa(i)=(n*((1-epsilon)*Pi/n+sin(epsilon*Pi/n)*(cos(epsilon*Pi/n)-sin(epsilon*Pi/n)* cotan(thet/2)))+&2*n*(r+i*de)*(sin(epsilon*Pi/n)/sin(thet/2)+(1-epsilon)*Pi/n)/l+&n*(r+i*de)**2*(Pi/n+Pi/2-thet/2-cotan(thet/2))/l**2)*l**2else if((i*de)>(l*sin(epsilon*Pi/n)/cos(thet/2)-r).and.(i*de)<=e1) thenSa(i)=2*n*((1-epsilon)*Pi/n+(r+i*de)*(Pi/n+asin(l*sin(epsilon*Pi/n)/(i*de+r)))/l)*l *lenApa(i)=n*((1-epsilon)*Pi*(1+(r+i*de)/l)**2/n+sin(epsilon*Pi/n)*(sqrt((r+i*de)**2/l* *2-&sin(epsilon*Pi/n)**2)+cos(epsilon*Pi/n))+(r+i*de)**2*(epsilon*Pi/n+&asin(l*sin(epsilon*Pi/n)/(r+i*de)))/l**2)*l**2write(10,"(f10.5,2X,f15.5)") i*de,Sa(i)write(20,"(f10.5,2X,f15.5)") i*de,Apa(i)end doclose(10)close(20)!调獭?用?子哩?程²序´计?算?装痢?药?内¸弹獭?道台?曲¸线?callinternal_ballistics0(d,len,e1,n,thet,epsilon,r,r1,l,a,pn,Ctz,At,k,density_p,mpeff,Isp) !输?出?星?型³装痢?药?的?几?何?参?数簓open(30,file="export_star_geometry.dat")write(30,"(A6,f11.5,A2)") "mpeff=",mpeff,"Kg"write(30,"(A2,f8.5,A2)") "D=",d,"m"write(30,"(A2,f8.5,A2)") "L=",len,"m"write(30,"(A2,I2)") "n=",nwrite(30,"(A3,f8.5)") "θ¯=",thetwrite(30,"(A3,f7.5)") "ε?=",epsilonwrite(30,"(A2,f7.5,A2)") "r=",r,"m"write(30,"(A3,f7.5,A2)") "r1=",r1,"m"write(30,"(A3,f8.5,A2)") "e1=",e1,"m"write(30,"(A3,f10.5,A2)") "Ap=",Ap,"㎡O"write(30,"(A2,f10.5)") "J=",Jwrite(30,"(A3,f7.5)") "η?=",etawrite(30,"(A3,f10.5,A2)") "Af=",Af,"㎡O"write(30,"(A4,f10.5)") "η?f=",etafwrite(30,"(A2,f10.5)") "l=",lwrite(30,"(A3,f10.5)") "y0=",y0write(30,"(A3,f10.5)") "y1=",y1close(30)stopend!该?子哩?程²序´用?于 ?计?算?零?维?变?截?面?燃?烧?装痢?药?的?内¸弹獭?道台?!利?用?龙ⅷ?格?-库a塔t法ぁ?计?算?内¸弹獭?道台?曲¸线?!可°用?于 ?计?算?侵?蚀骸?燃?烧?效§应畖下?的?内¸弹獭?道台?曲¸线?!--------------------------------------------------------------------------!d外猘径?,len长¤度¯,e1平?均·肉╝厚?,n星?角?数簓,thet星?边?夹D角?,epsilon角?度¯系μ数簓!r过y度¯圆2弧?半?径?,r1星?角?圆2弧?半?径?和³l特?征¶尺?寸?,key表括?示?是?否?考?虑?侵?蚀骸?燃?烧?!density_p推?进?剂®的?密¹度¯,k比括?热¯?比括?ga系μ数簓,C特?征¶速·度¯,mpeff药?柱·质±量?!a速·度¯系μ数簓,pn压1强?指?数簓,At喉³部?面?积y,Pc燃?烧?室酣?的?设Θ?计?压1力 ,Isp表括?示?理え?论?比括?冲?!dt时骸?间?步?长¤,e燃?层?厚?度¯,time时骸?间?,Sa同?一?时骸?刻²轴®向´的?各¶节¸点?的?燃?面?,Apa通 ?气?面?积y!Poc各¶节¸点?的?压1力 ,Pe喷?管¹出?口¸压1力 ,P_I压1力 冲?量?,ep侵?蚀骸?比括?!P_av平?均·压1力 ,F_av平?均·推?力 ,I0总哩?冲?,Im重?量?比括?冲?Iv体?积y比括?冲?,Cf推?力 系μ数簓!--------------------------------------------------------------------------subroutineinternal_ballistics0(d,len,e1,n,thet,epsilon,r,r1,l,a,pn,Ctz,At,k,density_p,mpeff,Isp)implicit nonereal,parameter::Pi=3.14integer :: n,keyreal(kind=8),intent(in) :: d,len,e1,thet,epsilon,r,r1,lreal(kind=8),intent(in) :: density_p,k,Ctz,mpeff,Ispreal(kind=8) :: a,pn,At,dt,e,time,Sa,Apa,Poc,Pe,P_I,F real(kind=8) :: ep,P_av,F_av,I0,Im,Iv,Cf,gareal(kind=8) :: f1,f2,f3,f4Pe=101325.0P_I=0.0!药?柱·的?能¹量?特?性?参?量?读®入?open(3,file="inter_parameter.dat")read(3,*) keyread(3,*) dtread(3,*) Pocclose(3)!计?算?点?火e压1强?时骸?的?推?力 值μPoc=Poc*101325.0ga=sqrt(k)*(2/(k+1))**((k+1)/(k-1)/2)Cf=ga*sqrt(2*k*(1-(Pe/Poc)**((k-1)/k))/(k-1))F=Cf*Poc*Ate=0.0time=0.0open(10,file="export_Poc.dat")open(20,file="export_F.dat")write(10,"(f8.5,2X,f10.5)") time,Poc/101325.0write(20,"(f8.5,2X,f15.5)") time,Fdo while(.true.)!计?算?燃?面?和³通 ?气?面?积yif(e<=r1) thenSa=2*n*(sin(epsilon*Pi/n)/sin(thet/2)+(1-epsilon)*Pi/n+&(r1+r)*(Pi/n+Pi/2-thet/2-cotan(thet/2))/l-(r1-e)*Pi/n/l)*l*len Apa=(n*((1-epsilon)*Pi/n+sin(epsilon*Pi/n)*(cos(epsilon*Pi/n)-&sin(epsilon*Pi/n)*cotan(thet/2)))+2*n*(r+e)*(sin(epsilon*Pi/n)/&sin(thet/2)+(1-epsilon)*Pi/n)/l+n*(r+e)**2*(Pi/n+Pi/2-thet/2-&cotan(thet/2))/l**2+n*(r1-e)**2*(thet/2+cotan(thet/2)-Pi/2)/l**2)*l**2else if(e>r1.and.e<=(l*sin(epsilon*Pi/n)/cos(thet/2)-r)) thenSa=2*n*(sin(epsilon*Pi/n)/sin(thet/2)+(1-epsilon)*Pi/n+&(e+r)*(Pi/n+Pi/2-thet/2-cotan(thet/2))/l)*l*lenApa=(n*((1-epsilon)*Pi/n+sin(epsilon*Pi/n)*(cos(epsilon*Pi/n)-&sin(epsilon*Pi/n)*cotan(thet/2)))+2*n*(r+e)*(sin(epsilon*Pi/n)/&sin(thet/2)+(1-epsilon)*Pi/n)/l+n*(r+e)**2*(Pi/n+Pi/2-thet/2-cotan(thet/2))/l**2)*l**2 else if(e>(l*sin(epsilon*Pi/n)/cos(thet/2)-r).and.e<=e1) thenSa=2*n*((1-epsilon)*Pi/n+(r+e)*(Pi/n+asin(l*sin(epsilon*Pi/n)/(e+r)))/l)*l*len Apa=n*((1-epsilon)*Pi*(1+(r+e)/l)**2/n+sin(epsilon*Pi/n)*(sqrt((r+e)**2/l**2-& sin(epsilon*Pi/n)**2)+cos(epsilon*Pi/n))+(r+e)**2*&(epsilon*Pi/n+asin(l*sin(epsilon*Pi/n)/(r+e)))/l**2)*l**2elseSa=0.0end if!计?算?侵?蚀骸?比括?if(key==1) thenif(Sa/Apa<=72.9) thenep=1.0elseep=1.3128-1.3249e-2*Sa/Apa+1.5527e-4*(Sa/Apa)**2-4.3868e-7*(Sa/Apa)**3end ifelseep=1.0end if!采°用?4阶¬的?龙ⅷ?格?库a塔t法ぁ?计?算?f1=ga**2*Ctz**2*(density_p*Sa*ep*a*Poc**pn-Poc*At/Ctz)/(Apa*len)f2=ga**2*Ctz**2*(density_p*Sa*ep*a*(Poc+dt*f1/2)**pn-(Poc+dt*f1/2)*At/Ctz)/(Apa*len )f3=ga**2*Ctz**2*(density_p*Sa*ep*a*(Poc+dt*f2/2)**pn-(Poc+dt*f2/2)*At/Ctz)/(Apa*len )f4=ga**2*Ctz**2*(density_p*Sa*ep*a*(Poc+dt*f3)**pn-(Poc+dt*f3)*At/Ctz)/(Apa*len) Poc=Poc+dt*(f1+2.0*f2+2.0*f3+f4)/6.0Cf=ga*sqrt(2*k*(1-(Pe/Poc)**((k-1)/k))/(k-1))F=Cf*Poc*Attime=time+dte=e+dt*ep*a*Poc**pnif(e>e1.and.Poc<101325.0) exitwrite(10,"(f8.5,2X,f10.5)") time,Poc/101325.0write(20,"(f8.5,2X,f15.5)") time,FP_I=P_I+Poc*dtend doclose(10)close(20)!循-环«结®束?P_av=P_I/timeCf=ga*sqrt(2*k*(1-(Pe/P_av)**((k-1)/k))/(k-1))F_av=Cf*P_av*AtI0=F_av*timeIm=I0/mpeff/9.81Iv=I0/(Apa*len)!输?出?固²体?火e箭y发ぁ?动ˉ机¸的?工¤作痢?参?数簓值μopen(30,file="export_propellant.dat")write(30,"('平?均·压1力 p_av=',f9.5)") p_av/101325.0write(30,"('平?均·推?力 F_av=',f10.5,1X,'KN')") F_av/1000.0 write(30,"('总哩?冲?I0=',f10.5,1X,'KN.s')") I0/1000.0write(30,"('重?量?比括?冲?Im=',f10.5,1X,'s')") Imwrite(30,"('体?积y比括?冲?Iv=',f10.5,1X,'KN.s/m3')") Iv/1000.0 write(30,"('工¤作痢?时骸?间?time=',f8.5,1X,'s')") timeclose(30)returnend subroutine。