■ 问题提出
在博文 逆Laplace数值逆变换[1] 给出了数值计算Laplace逆变换的简易程序。其中存在以下几个问题需要讨论:
- 问题1: 程序实现过程原理以及优化;
- 问题2: 运算参数:,,对于积分数值的影响。从上篇博文中明显看到有些计算出的结果大大偏离的实际值。比如其中的函数所计算出来的幅值超过了1.
![4153450ba0499fbb8e16f5b4bd145645.png](https://img-blog.csdnimg.cn/img_convert/4153450ba0499fbb8e16f5b4bd145645.png)
01程序实现原理以及优化
1.Laplace逆运算
![742bcae4ede319f3aa977a14e429caeb.png](https://img-blog.csdnimg.cn/img_convert/742bcae4ede319f3aa977a14e429caeb.png)
所以,算法的核心是对进行傅里叶反变换,然后在乘以。
由于确认变换后的函数是实函数,因此,为了节省计算时间,只对傅里叶反变换的积分,进行正半轴的积分,同时积分的上限有参数决定。对积分的数值取齐实部,再乘以2便可以得到。
![b3aa8922274e2aa9989a069fabdd158e.png](https://img-blog.csdnimg.cn/img_convert/b3aa8922274e2aa9989a069fabdd158e.png)
2.Python积分程序实现优化
使用梯形积分来实现函数的积分,可以获得更精确的积分值。理论分析可知:
![459d966188c6668b93dcd381c97774c4.png](https://img-blog.csdnimg.cn/img_convert/459d966188c6668b93dcd381c97774c4.png)
其中 以及 。那么积分误差上限为:
![279832f01c44dd6a6487862fc0935a63.png](https://img-blog.csdnimg.cn/img_convert/279832f01c44dd6a6487862fc0935a63.png)
其中,, 。
![f36e9b88a1c9950dcf259cbdc8c5a0ad.png](https://img-blog.csdnimg.cn/img_convert/f36e9b88a1c9950dcf259cbdc8c5a0ad.png)
▲ 图2.1 梯形积分方法示意图
计算可以由两种方式:
方式1: 对左、右黎曼积分加权平均:
![6df29bbfb88285acf6413844e46f19c8.png](https://img-blog.csdnimg.cn/img_convert/6df29bbfb88285acf6413844e46f19c8.png)
方式2: 利用如下的公式计算:
![51ae1eb0c25d09fa3870e1b9c2d00f01.png](https://img-blog.csdnimg.cn/img_convert/51ae1eb0c25d09fa3870e1b9c2d00f01.png)
最后实现的代码为: 利用上述公式,对进行积分测试:
- ounter(line
printf(trapz(sin, 0, pi/2, 1000))
所得到结果为:0.9999997943832332。 可以看到使用n=1000对应的结果后面的积分精度达到了小数点后面6位9的小数点的位数。
3.实现Laplace逆运算
通过上面梯形积分方法,实现Laplace数值逆变换,具体的子程序如下面所示。
02一些基本函数的实验
下面通过对一些基本常见函数的laplace变换,来测试一下上述程序的性能。
Ⅰ.sin(t)
![9c129aea5bb35bde3e725e4895d8887a.png](https://img-blog.csdnimg.cn/img_convert/9c129aea5bb35bde3e725e4895d8887a.png)
sigma=0.2, omiga=200, nint=omiga*50
![573d789a8d8888a511935162cfab0aca.png](https://img-blog.csdnimg.cn/img_convert/573d789a8d8888a511935162cfab0aca.png)
Ⅱ.exp(-t)
![b23225988229cae234f2269ef2ea99fb.png](https://img-blog.csdnimg.cn/img_convert/b23225988229cae234f2269ef2ea99fb.png)
sigam=-1+0.1, omiga=200, nint=omiga*50
![22931baf24e14948b9cb969d150b9b91.png](https://img-blog.csdnimg.cn/img_convert/22931baf24e14948b9cb969d150b9b91.png)
Ⅲ.u(t)
![a8a09c5ad5c42f2c0df88fd733498c7e.png](https://img-blog.csdnimg.cn/img_convert/a8a09c5ad5c42f2c0df88fd733498c7e.png)
![85036d32f92142a3e6beab8d8c6fc621.png](https://img-blog.csdnimg.cn/img_convert/85036d32f92142a3e6beab8d8c6fc621.png)
Ⅳ.u(t-1)
![fff0d5e10b7be797ac23477aa6b39a6e.png](https://img-blog.csdnimg.cn/img_convert/fff0d5e10b7be797ac23477aa6b39a6e.png)
![56faab20e9a2744c6604be4cb3265ba7.png](https://img-blog.csdnimg.cn/img_convert/56faab20e9a2744c6604be4cb3265ba7.png)
Ⅴ.周期化脉冲信号
![87232598653f7a23e093d35d39333b19.png](https://img-blog.csdnimg.cn/img_convert/87232598653f7a23e093d35d39333b19.png)
![c076ac6d166015c8bef450fd42432869.png](https://img-blog.csdnimg.cn/img_convert/c076ac6d166015c8bef450fd42432869.png)
※ 结论
通过原理分析,可以获得建议的Laplace数值逆运算的正确的PYTHON算法程序。
这个程序是直接对Laplace反变换公式利用梯形积分方法获得计算结果。通过对几种常见的信号Laplace的反变换,验证了这个算法的正确性。
通过在此过程,可以看到,对于参数sigma, omiga, nint对于计算结果还是有很大的影响。另外,对于时间t,只能在比较小的范围内有效,当t超过一定长度,前面所计算的结果都会出现比较大的误差。
![b306a5fa7e229909f07d172a81d722a8.png](https://img-blog.csdnimg.cn/img_convert/b306a5fa7e229909f07d172a81d722a8.png)
▲ 对于方波进行Laplace数值反变换的结果
参考资料
[1]
逆Laplace数值逆变换: https://zhuoqing.blog.csdn.net/article/details/107241738