【问题标题】:how to avoid singularity in numerical integration?如何避免数值积分中的奇异性?
【发布时间】:2015-07-01 16:22:54
【问题描述】:

我想在 Maple18 中绘制 E2(t),但它在 t=0 左右爆炸。如何对这个积分应用渐近评估来避免奇异性并解决这个问题?

restart: with(plots):
with(Student[NumericalAnalysis]):
g1:=(x,t)->(-sqrt(t)/(2*sqrt(Pi*r^3)))*(sin((r*(x-1)^(2)/(4*t))+(Pi/4))-sin((r*(x+1)^(2)/(4*t))+(Pi/4))):

g2:=(x,t)->int(g1(x,t),r=1..infinity);

g3:=(x,t) -> (diff(g2(x,t),t)):


g4:=(x,t) -> (diff(g2(x,t),x,x)):


g5:=(x,t) -> ((1/2)*(g3(x,t)^2+g4(x,t)^2)):


E2:=t->(int(g5(x,t),x=0..100)):

evalf(E2(0));
Error, (in g1) numeric exception: division by zero
evalf(E2(1));
Error, (in g3) invalid input: diff received 1, which is not valid for its 2nd argument
plot(E2(t),t=0..20);

最好的问候,

【问题讨论】:

  • 我刚刚在 MaplePrimes 上回答了这个问题。当人们同时在不同的论坛上问同样的问题时,这很烦人。
  • 卡尔,在这里复制粘贴你的答案。
  • @Carl 感谢您的回答并对这种情况表示抱歉。这对我来说是一个紧迫的问题,我想增加获得答案的机会。
  • @2che:您对我在 MaplePrimes 上的回答没有任何回应令人不安。
  • @lolmaus - Andrey Mikhaylov:好的,完成。

标签: plot integration maple expansion


【解决方案1】:

我不知道您所说的“避免奇点”是什么意思。你的函数在 t=0 时是奇异的,没有办法改变它。只是不要试图在 0 处评估它。

您的另一个问题是由于尝试对数字进行区分而引起的。您需要做的是使用 unapply 创建过程 E2。然后将象征性地进行区分。本质上,这是您的代码。我删除了多余的括号、(x,t)-> 和不需要的包加载。

restart:
g1:= -sqrt(t/Pi/r^3)/2 * 
     (sin(r*(x-1)^2/4/t+Pi/4)-sin(r*(x+1)^2/4/t+Pi/4)):
g2:= int(g1, r= 1..infinity):
g3:= diff(g2,t):
g4:= diff(g2,x$2):
g5:= (g3^2+g4^2)/2:

这里是E2

E2:= unapply(Int(g5, x= 0..100, epsilon= 1e-4, digits= 7), t):

注意 Int 中的大写 I。这将避免浪费时间尝试符号积分,并在应用 evalf 时直接进入数值积分。 epsilondigits 参数会降低精度并提高速度。对于绘图,精度仍然足够。绘图命令是

plot(E2, 0..20, numpoints= 50, labels= [t, ``]);

这需要两到三分钟。如果您增加 numpoints 以获得更准确的绘图,时间将成比例增加。

【讨论】:

    猜你喜欢
    • 2020-03-22
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 2014-04-17
    • 1970-01-01
    • 2013-05-19
    相关资源
    最近更新 更多