1一、问题的重述考虑航天器在仅受到地球万有引力、航天器自身发动机作用力的作用下作平面运动,将地球和航天器视为质点,建立航天器运动的数学模型。显然这样的数学模型在精度上是远远不能满足实际需要的,在其他要求精确制导等有关高科技的实际问题中,我们都面临着类似的问题:我们必须建立高精度的数学模型,必须高精度地估计模型中的大批参数,因为只有这样的数学模型才能解决实际问题,而不会出现差之毫厘,结果却失之千里的情况。由于航天器的问题太复杂,本题仅考虑较简单的确定高精度参数问题。假设有一个生态系统,其中含有两种生物,即:A生物和B生物,其中A生物是捕食者,B生物是被捕食者。假设t时刻捕食者A的数目为xt,被捕食者B数目为yt,它们之间满足以下变化规律:1234xtxtytytytxt初始条件为:0506xtyt其中16kk为模型的待定参数。通过对此生态系统的观测,可以得到相关的观测数据。要利用有关数据,解决以下问题:1)在观测数据无误差的情况下,若已知2,求其它5个参数1,3,4,5,6kk?2)若2也未知,至少需要多少组观测数据,才能确定参数16kk?3)在观测资料有误差(时间变量不含有误差)的情况下,确定参数16kk在某种意义下的最优解,并与仿真结果比较,进而改进数学模型。4)假设连观测资料的时间变量也含有误差,确定参数k在某种意义下的最优解。二、航天器运动模型的建立考虑航天器在仅受到地球万有引力、航天器自身发动机作用力的作用下作平面运动,将地球和航天器视为质点,由理论力学可知,一个刚体在空间的运动可以看作质心的移动,因此可以应用质心运动定理来研究刚体质心的移动规律。以地球中心为原点,建立直角坐标系,航天器绕地球飞行,可以出现在该直角坐标系中四个象限的任意一个之内。平面直角坐标系如图1。符号说明如下:xv——航天器在x方向的速度2yv——航天器在y方向的速度1F——万有引力,1222MmMmFGGrxy2F——航天器发动机作用力,为控制变量——万有引力与x轴正方向的夹角——航天器发动机作用力与x轴正方向的夹角0t——初始时刻0x,0y——航天器初始位置0xv——航天器x方向初速度0yv——航天器y方向初速度航天器受的万有引力1F方向指向地球中心(原点),航天器受推力2F的方向与x轴正方向成角。将1F和2F投影到该直角坐标系上,见图1图1航天器受力分解图其中,2122222coscoscosxMmxFFFFGxtytxtyt2122222sinsinsinyMmyFFFFGxtytxtyt初始条件为300000000()()()xxyyxtxytyvtvvtvxt,yt都是关于时间的位置函数,航天器在x方向的分速度即xt对时间t求导,航天器在y方向的分速度即yt对时间t求导,航天器在x方向的加速度即xt对时间t求二阶导,航天器在y方向的加速度即yt对时间t求二阶导,根据牛顿第二定律有方程(3)和(4)。由此建立的航天器模型如下:2cosxyxxxtvytvFFvtm222222222(1)sinyyxtMGmxtytxtytFytFMvtGmmxtytxtyt显然这样的数学模型在精度上是远远不能满足实际需要的,在其他要求精确制导等有关高科技的实际问题中,我们都面临着类似的问题:我们必须建立高精度的数学模型,必须高精度地估计模型中的大批参数,因为只有这样的数学模型才能解决实际问题,而不会出现差之毫厘,结果却失之千里的情况。这时所建立数学模型的精度就成了数学模型的生命线。例如上述问题中的航天器还要受到地球质量分布不均匀所引起的摄动力,大气阻力,日、月及其它星球的摄动引力的影响,以及航天器发动机为调整航天器自身姿态运作时作用力的影响。这样不但数学模型十分复杂,而且在这些数学模型中还要涉及到许多重要的参数,如地球的引力场模型就有许多待定参数。不仅如此,在对航天器进行测量时,还涉及到观测站的地理位置以及设备的系统误差等参数。为此人们要设法利用长期积累的丰富的观测资料,高精度确定这些重要的参数。由于航天器的问题太复杂,下面本题仅考虑较简单的确定高精度参数问题。三、捕食者与被捕食者生态系统问题的分析题中假设有一个生态系统,含有两种生物,A生物和B生物,A生物是捕食者,B生物是被捕食者。假设t时刻捕食者A的数目为xt,被捕食者B数目为yt,它们之间满足以下变化规律:1234xtxtytytytxt4初始条件为:0506xtyt该模型中,1捕食者独自存在时死亡率,10;2被捕食者对捕食者的供养能力20;3是被捕食者的独立生存增长率,30;4是捕食者掠取被捕食者的能力,40。[2]这个方程就是生态系统中被捕食者与捕食者的volterra模型,16kk为模型的待定参数。对于该模型理论上不存在解析解,因此我们不能通过参数拟合确定模型的参数。Volterra模型在给定参数和初始值的情形下可以采用数值积分获得任意时间点的数值解。根据volterra模型进行一些公式推导如下:(2)dxt=xt+ytdtdyt=ytxtdt两个方程相除得:yt+xtdyt=dxtxt+yt移项得:34ytytdyt=dxtytxt两边积分:(3)00yxyxytytdyt=dxtytxt得到相轨方程:(lnln)()(lnln)()0(4)0000y-y+y-y-x-x-x-x=移项得:lnlnlnln0000y+y-x-x=y+y-x-x5该式右边为lnln0000y+y-x-x只与系统初始状态有关,令lnln0000y+y-x-x=C(5)易知0C,将(5)式代入(4)式得到lnln(6)y+y-x-x=C式中方程两边同除以C,得:lnln1y+y-x-x=CCCC(7)在(7)式中,令1C=,22C=,33C=,44C=得到1234lnln1y+y-x-x=这一方程体现了Volterra模型中两个变量之间的变化关系,我们称此方程为相轨方程。进一步研究相轨方程,可以发现Volterra模型中两个变量呈现周期性变化。第一问,对k来说,相轨方程是一个4未知数的方程,DATA1中有6组数据,用6组数据确定4个k可以采用极小范数最小二乘解。又因为2已知,C可求,从而k14k可求,由于各观测值真实准确,5,6可取DATA1中任意一组数据,不失一般性,我们取第一组数据为初始值。第二问,我们可以证明参数C与系统周期成反比,由参考资料可以知道volterra模型中k的含义,从而确定C的正负性,在2未知的情况下,求C可以从C与时间的关系入手,我们先在DATA1的6组数据中取4组算出k,然后设计一个1C的新生态系统,以无误差的观测数据DATA1为准,设计搜索算法找到与DATA1中x,y值极为接近的数值点,找到对应的观测时间,得到观测间隔,这个观测间隔与DATA1中已知的观测间隔一起可以求出C,从而得6到k14k,5,6同第一问。第三问,用所有数据求得极小范数最小二乘解,可以确定k。经过与第二问类似的方法获得C。进一步求出k14k,5,6可取DATA2中观测初始时间的值,这一套k16k就是我们所求的最小二乘意义下的最优解。将这一组参数带入volterra模型,获得各观测点上的仿真结果。通过与观测结果比较,我们发现误差普遍较大。于是我们改进了参数估计模型,改为求取均方误差意义下的最优解。获得了较好的效果。第四问,我们采取了第三问中的改良算法,以求取使x,y,t三者均方误差最小的参数组为目标,进行计算,然而误差较大。经过判断,我们认为这是由于时间和数据均存在误差导致搜索结果不够精确,我们改进搜索算法,结果大为改善。四、模型的建立及求解问题一,2已知,求其它5个参数1,3,4,5,6kk结合DATA1.TXT中6组无误差的观测数据(包括了观测时刻jt、A生物数目)(jtx、B生物数目)(jty),(7)式含有4个未知数,而题中提供了6组数据,写为矩阵形式即:1111222213333234444455556666lnlnlnln1lnln1lnln11lnlnlnlnyyxxyyxxyyxxyyxxyyxxyyxx记为矩阵形式AX=b,其中111122223333444455556666lnlnlnlnlnlnlnlnlnlnlnlnyyxxyyxxyyxxyyxxyyxxyyxxA,1234X,1111b,矛盾方程组AX=b的惟一极小范数最小二乘解为XAb,采用极小范数最小二乘解,得到kk(14)的值,如表1所示7表1最小二乘的kk(14)的值12340.14471560254457-0.01447146313960-0.868317924024250.07236037259773因为22C=,由题意215=,而从volterra模型本身出发,2是捕食者掠取被捕食者的能力,所以利用DATA1中数据算出的20,所以C0,这一问中C=13.82030262390574-,把C代入1C=,33C=,44C=得到表22已知,1,3,4kk的值12.00001342156688=-312.00041648377672=41.00004224727917=-对于5,6来说,05xt,06yt,由观测数据)(jtx,)(jty已知,而0t未知,所以,56,,可以是观测数据中任意一组)(jtx,)(jty的值。不失一般性我们取00t,由DATA1知道,0()10xt=0()60yt=,所以510=,660=。问题二,2未知,至少需要几组数据才能确定(16)kk的值1.k的求解由问题一的分析可知,为了求取k,至少需要四个方程,即四组观测数据。列为线性方程11111222223333344444lnln1lnln1lnln11lnlnyyxxyyxxyyxxyyxx可以解得1C=,22C=,33C=,44C=2.观测系统时间变换分析8对于同样系统的观测,当选取的观测时间起点和观测时间单位不同时,得到同一物理系统的两种不同时间观测结果。不失一般性,我们将两种观测之间的时间变换关系表示为21t=at+b,其中,12tt表示两种观测系统的时间单位数量,b反映了两个观测系统时间起点上的差异,而a反映了两个观测系统观测时间单位的比例,线性关系是由时间的均匀性确定的。定理一:对于参数1234