численное решение интеграла от экспоненты в Mathematica

[salang]

Leierkastenmann в сообщении #812770 писал(а):

Это график именно того, с чем мы мучаемся? Слабо верится, значения не похожи. В то, что это какой-то аналог, поверить могу. В том же конкретном случае, что есть у нас тут, не будет такого плавного спуска, будет тоненький пик в нанообласти, который дальше резко падает в ноль. Замените LogListPlot на ListPlot, можете даже домножить на коэффициент и полюбуйтесь, что получится.

Да, это нужная форма сигнала. Должен быть короткий пик с крутым фронтом и резким спадом, а дальше очень плавное уменьшение до 1...2мкс

[salang]

не получилось с дискретным временем

(Оффтоп)

9. nb

In[23]:= c = 3*10^8;
H = 7.2*10^5;
\[Tau] = 47*10^-6;
Fd = 7200;
Tp = 50*10^-6;
f = 3.5*10^8;
v = 7500;
Q = 18000;
L = 100;
s2 = 0.003;
s1 = 2;
\[CapitalTheta]0 = 0.018;
m = 256;
\[Lambda] = 0.022;
n = 3;
lm = 25;
d = 1.2;
\[Lambda] = 0.022;
n = 3;
t = List[-10^-8, -20^-9, -10^-9, 0, 10^-11, 10^-10, 10^-9, 10^-8, 
   20*10^-7, 2*10^-6, 10^-6];
p := Sum[NIntegrate[
    Exp[(\[Pi]*Fd^2*Tp^2*(m - Abs[k])^2*x^2)/
       H^2 - (n^2*\[Pi]^2*d^2*y^2)/(\[CapitalTheta]0^2*\[Lambda]^2*
         H^2*141^2) - (100*lm^2*x^2*s1^2)/(H^2*L^2*\[Pi]^2*
         s2^2) - (100*lm^2*y^2*s1^2)/(H^2*L^2*\[Pi]^2*
         s2^2) - (5.55*(y^2 + 
           x^2))/(H^2*\[CapitalTheta]0^2) - \[Pi]*(f^2*(tt - k*Tp - 
             x^2/(c*H) - y^2/(c*H))^2 + 
         2*Q*(v*k*Tp)/(H*
             d) + \[Tau]^2*((v*k*Tp)/(H*d))^2)], {x, -Infinity, 
     Infinity}, {y, -Infinity, Infinity}, MaxRecursion -> 40, 
    AccuracyGoal -> 60, Method -> "AdaptiveMonteCarlo"], {k, -255, 
    255}];
ListPlot[Transpose[{t, p}]]



During evaluation of In[23]:= NIntegrate::inumr: The integrand E^(-5.42915*10^-7 x^2-6.22059*10^-7 y^2-3.30433*10^-8 (x^2+y^2)-\[Pi] (-3.98438+1.225*10^17 (51/4000+tt+Times[<<2>>]+Times[<<2>>])^2)) has evaluated to non-numerical values for all sampling points in the region with boundaries {{0,1},{0,1}}. >>

During evaluation of In[23]:= NIntegrate::inumr: The integrand E^(-5.42913*10^-7 x^2-6.22059*10^-7 y^2-3.30433*10^-8 (x^2+y^2)-\[Pi] (-3.96875+1.225*10^17 (127/10000+tt+Times[<<2>>]+Times[<<2>>])^2)) has evaluated to non-numerical values for all sampling points in the region with boundaries {{0,1},{0,1}}. >>

During evaluation of In[23]:= NIntegrate::inumr: The integrand E^(-5.42909*10^-7 x^2-6.22059*10^-7 y^2-3.30433*10^-8 (x^2+y^2)-\[Pi] (-3.95313+1.225*10^17 (253/20000+tt+Times[<<2>>]+Times[<<2>>])^2)) has evaluated to non-numerical values for all sampling points in the region with boundaries {{0,1},{0,1}}. >>

During evaluation of In[23]:= General::stop: Further output of NIntegrate::inumr will be suppressed during this calculation. >>

During evaluation of In[23]:= Transpose::nmtx: The first two levels of the one-dimensional list {{-(1/100000000),-(1/512000000000),-(1/1000000000),0,1/100000000000,1/10000000000,1/1000000000,1/100000000,1/500000,1/500000,1/1000000},511 NIntegrate[Exp[(\[Pi] Power[<<2>>] Power[<<2>>] Power[<<2>>] Power[<<2>>])/Power[<<2>>]-<<4>> Power[<<2>>]-<<4>> Power[<<2>>]-<<4>> Power[<<2>>]-<<2>> Power[<<2>>]-\[Pi] Plus[<<3>>]],{x,-\[Infinity],\[Infinity]},{y,-\[Infinity],\[Infinity]},MaxRecursion->40,AccuracyGoal->60,Method->AdaptiveMonteCarlo]} cannot be transposed. >>

During evaluation of In[23]:= Transpose::nmtx: The first two levels of the one-dimensional list {{-(1/100000000),-(1/512000000000),-(1/1000000000),0,1/100000000000,1/10000000000,1/1000000000,1/100000000,1/500000,1/500000,1/1000000},511 NIntegrate[Exp[(\[Pi] Power[<<2>>] Power[<<2>>] Power[<<2>>] Power[<<2>>])/Power[<<2>>]-<<4>> Power[<<2>>]-<<4>> Power[<<2>>]-<<4>> Power[<<2>>]-<<2>> Power[<<2>>]-\[Pi] Plus[<<3>>]],{x,-\[Infinity],\[Infinity]},{y,-\[Infinity],\[Infinity]},MaxRecursion->40,AccuracyGoal->60,Method->AdaptiveMonteCarlo]} cannot be transposed. >>

During evaluation of In[23]:= Transpose::nmtx: The first two levels of the one-dimensional list {{-(1/100000000),-(1/512000000000),-(1/1000000000),0,1/100000000000,1/10000000000,1/1000000000,1/100000000,1/500000,1/500000,1/1000000},511 NIntegrate[Exp[(\[Pi] Power[<<2>>] Power[<<2>>] Power[<<2>>] Power[<<2>>])/Power[<<2>>]-<<4>> Power[<<2>>]-<<4>> Power[<<2>>]-<<4>> Power[<<2>>]-<<2>> Power[<<2>>]-\[Pi] Plus[<<3>>]],{x,-\[Infinity],\[Infinity]},{y,-\[Infinity],\[Infinity]},MaxRecursion->40,AccuracyGoal->60,Method->AdaptiveMonteCarlo]} cannot be transposed. >>

During evaluation of In[23]:= General::stop: Further output of Transpose::nmtx will be suppressed during this calculation. >>

During evaluation of In[23]:= ListPlot::lpn: Transpose[{{-(1/100000000),-(1/512000000000),-(1/1000000000),0,1/100000000000,1/10000000000,1/1000000000,1/100000000,1/500000,1/500000,1/1000000},511 NIntegrate[Exp[<<5>> Power[<<2>>]-Times[<<2>>]-Times[<<2>>]-Times[<<2>>]-Times[<<2>>]-Times[<<2>>]],{x,-\[Infinity],\[Infinity]},{y,-\[Infinity],\[Infinity]},MaxRecursion->40,AccuracyGoal->60,Method->AdaptiveMonteCarlo]}] is not a list of numbers or pairs of numbers. >>

Out[44]= ListPlot[
 Transpose[{{-(1/100000000), -(1/512000000000), -(1/1000000000), 0, 1/
    100000000000, 1/10000000000, 1/1000000000, 1/100000000, 1/500000, 
    1/500000, 1/1000000}, 
   511 NIntegrate[
     Exp[(\[Pi] Fd^2 Tp^2 (m - Abs[k])^2 x^2)/H^2 - (
       n^2 \[Pi]^2 d^2 y^2)/(\[CapitalTheta]0^2 \[Lambda]^2 H^2 \
141^2) - (100 lm^2 x^2 s1^2)/(H^2 L^2 \[Pi]^2 s2^2) - (
       100 lm^2 y^2 s1^2)/(H^2 L^2 \[Pi]^2 s2^2) - (
       5.55 (y^2 + x^2))/(
       H^2 \[CapitalTheta]0^2) - \[Pi] (f^2 (tt - k Tp - x^2/(c H) - 
             y^2/(c H))^2 + (2 Q (v k Tp))/(
          H d) + \[Tau]^2 ((v k Tp)/(
            H d))^2)], {x, -\[Infinity], \[Infinity]}, {y, -\
\[Infinity], \[Infinity]}, MaxRecursion -> 40, AccuracyGoal -> 60, 
     Method -> "AdaptiveMonteCarlo"]}]]

,
с непререрывным тоже

(Оффтоп)

7. nb
In[66]:= c = 3*10^8;
H = 7.2*10^5;
\[Tau] = 47*10^-6;
Fd = 7200;
Tp = 50*10^-6;
f = 3.5*10^8;
v = 7500;
Q = 18000;
L = 100;
s2 = 0.003;
s1 = 2;
\[CapitalTheta]0 = 0.018;
m = 256;
\[Lambda] = 0.022;
n = 3;
lm = 25;
d = 1.2;
\[Lambda] = 0.022;
n = 3;
p[t] := Sum[
   NIntegrate[
    Exp[(\[Pi]*Fd^2*Tp^2*(m - Abs[k])^2*x^2)/
       H^2 - (n^2*\[Pi]^2*d^2*y^2)/(\[CapitalTheta]0^2*\[Lambda]^2*
         H^2*141^2) - (100*lm^2*x^2*s1^2)/(H^2*L^2*\[Pi]^2*
         s2^2) - (100*lm^2*y^2*s1^2)/(H^2*L^2*\[Pi]^2*
         s2^2) - (5.55*(y^2 + 
           x^2))/(H^2*\[CapitalTheta]0^2) - \[Pi]*(f^2*(tt - k*Tp - 
             x^2/(c*H) - y^2/(c*H))^2 + 
         2*Q*(v*k*Tp)/(H*
             d) + \[Tau]^2*((v*k*Tp)/(H*d))^2)], {x, -Infinity, 
     Infinity}, {y, -Infinity, Infinity}, MaxRecursion -> 40, 
    AccuracyGoal -> 60, Method -> "AdaptiveMonteCarlo"], {k, -255, 
    255}];
Plot[p[t], {-10^(-8), 10^(-6)}, PlotPoints -> 30, MaxRecursion -> 0]



During evaluation of In[66]:= NIntegrate::inumr: The integrand E^(-5.42915*10^-7 x^2-6.22059*10^-7 y^2-3.30433*10^-8 (x^2+y^2)-\[Pi] (-3.98438+1.225*10^17 (51/4000+tt+Times[<<2>>]+Times[<<2>>])^2)) has evaluated to non-numerical values for all sampling points in the region with boundaries {{0,1},{0,1}}. >>

During evaluation of In[66]:= NIntegrate::inumr: The integrand E^(-5.42913*10^-7 x^2-6.22059*10^-7 y^2-3.30433*10^-8 (x^2+y^2)-\[Pi] (-3.96875+1.225*10^17 (127/10000+tt+Times[<<2>>]+Times[<<2>>])^2)) has evaluated to non-numerical values for all sampling points in the region with boundaries {{0,1},{0,1}}. >>

During evaluation of In[66]:= NIntegrate::inumr: The integrand E^(-5.42909*10^-7 x^2-6.22059*10^-7 y^2-3.30433*10^-8 (x^2+y^2)-\[Pi] (-3.95313+1.225*10^17 (253/20000+tt+Times[<<2>>]+Times[<<2>>])^2)) has evaluated to non-numerical values for all sampling points in the region with boundaries {{0,1},{0,1}}. >>

During evaluation of In[66]:= General::stop: Further output of NIntegrate::inumr will be suppressed during this calculation. >>
During evaluation of In[66]:= SetDelayed::write: Tag Times in (511 NIntegrate[Exp[(\[Pi] Fd^2 Tp^2 Plus[<<2>>]^2 x^2)/H^2-(Power[<<2>>] Power[<<2>>] Power[<<2>>] Power[<<2>>])/Times[<<4>>]-(100 Power[<<2>>] Power[<<2>>] Power[<<2>>])/Times[<<4>>]-(100 Power[<<2>>] Power[<<2>>] Power[<<2>>])/Times[<<4>>]-(5.55 Plus[<<2>>])/Times[<<2>>]-\[Pi] (Times[<<2>>]+Times[<<3>>]+Times[<<2>>])],{x,-\[Infinity],\[Infinity]},{y,-\[Infinity],\[Infinity]},MaxRecursion->40,AccuracyGoal->60,Method->AdaptiveMonteCarlo])[{-(1/100000000),-(1/512000000000),-(1/1000000000),0,1/100000000000,1/10000000000,1/1000000000,1/100000000,1/500000,1/500000,1/1000000}] is Protected. >>

During evaluation of In[66]:= Plot::pllim: Range specification {-(1/10^8),1/10^6} is not of the form {x, xmin, xmax}. >>

Out[86]= Plot[p[t], {-(1/10^8), 1/10^6}, PlotPoints -> 30, 
 MaxRecursion -> 0]

. Для ускорения последнего расчета явно указал 100 отсчетов в PlotPoints.

[Aritaborian]

Во втором случае у вас p[t]:=..., когда должно быть p[t_]:=.... И в Plot неправильно указаны пределы: должно быть {t, -10^(-8), 10^(-6)}, а у вас {-10^(-8), 10^(-6)}.
Это первое, что я заметил; может, ещё посмотрю на код повнимательнее.

[salang]

Aritaborian в сообщении #813849 писал(а):

Во втором случае у вас p[t]:=..., когда должно быть p[t_]:=.... И в Plot неправильно указаны пределы

c этими исправлениями ошибок нет, но график пустой и что удивительно диапазон значений p от -1 до 1

[Aritaborian]

salang в сообщении #813870 писал(а):

что удивительно диапазон значений p от -1 до 1

На пустом графике это диапазон значений по умолчанию. Насчёт остального не знаю пока; может, я дал вам неправильные советы ;-)

[salang]

с нетерпением жду результата

[Aritaborian]

Ну вот снова у вас эта элементарная невнимательность.
Почему сначала t, а потом вдруг tt? И величины в t почему-то не по порядку возрастания...
Итак, дискретный случай.

c = 3*10^8; H = 7.2*10^5; \[Tau] = 47*10^-6; Fd = 7200; Tp = 50*10^-6; f = 3.5*10^8; v = 7500; Q = 18000; L = 100; s2 = 0.003; s1 = 2; \[CapitalTheta]0 = 0.018; m = 256; \[Lambda] = 0.022; n = 3; lm = 25; d = 1.2; \[Lambda] = 0.022; n = 3;
t = Sort@{-10^-8, -20^-9, -10^-9, 0, 10^-11, 10^-10, 10^-9, 10^-8, 
    20*10^-7, 2*10^-6, 10^-6};

p = Sum[NIntegrate[Exp[(\[Pi] Fd^2 Tp^2 (m - Abs[k])^2 x^2)/
       H^2 - (n^2 \[Pi]^2 d^2 y^2)/(\[CapitalTheta]0^2 \[Lambda]^2 H^2*141^2) - (100 lm^2 x^2 s1^2)/(H^2 L^2 \[Pi]^2 s2^2) - (100 lm^2 y^2*s1^2)/(H^2 L^2 \[Pi]^2 s2^2) - (5.55 (y^2 + x^2))/(H^2 \[CapitalTheta]0^2) - \[Pi] (f^2 (t - k Tp - x^2/(c H) - y^2/(c H))^2 + 2 Q (v k Tp)/(H d) + \[Tau]^2 ((v k Tp)/(H d))^2)], {x, -Infinity, Infinity}, {y, -Infinity, Infinity}, MaxRecursion -> 40, AccuracyGoal -> 60, Method -> "AdaptiveMonteCarlo"], {k, -255, 255}];

ListLogPlot[Transpose[{t, p}], Joined -> True, PlotRange -> All]

[salang]

Aritaborian в сообщении #813904 писал(а):

Почему сначала t, а потом вдруг tt? И величины в t почему-то не по порядку возрастания...Итак, дискретный случай.

С ранжировкой согласен, проглядел. А где именно tt в моем дискретном варианте?
Самое главное- график с ListPlot как надо, т.е. после среза с высокой крутизной идет очень медленное уменьшение, но это из-за масштаба очень плохо видно. Как растянуть график по вертикали в окрестностях значения р 1000...10000 при линейном масштабе?
А с непрерывным временем что-то возможно сделать?

[Aritaborian]

salang в сообщении #813974 писал(а):

А где именно tt в моем дискретном варианте?

Возьмите да просмотрите внимательно.
Над остальным пока думаю. Надеюсь, кто-нибудь присоединится.
Кстати, чем вас не устраивает мой последний график?

[salang]

Aritaborian в сообщении #813990 писал(а):

salang в сообщении #813974 писал(а):

А где именно tt в моем дискретном варианте?

Возьмите да просмотрите внимательно

не смог найти

Aritaborian в сообщении #813974 писал(а):

чем вас не устраивает мой последний график?

График замечательный, но ListLogPlot показывает только вершину импульса, ListPlot еще и срез. Нужен график на интервале времени 10нс-1мкс c соответствующим диапазоном p (1000 ?). Если это невозможно на одном графике отобразить, может на нескольких получится?

[salang]

а почему при максимальном указанном в исходном массиве времени 1мкс график строится аж до 20мкс (а там точно ничего не может быть)? Не исправляется даже прямым указанием

PlotRange ->{t, -10^(-8), 10^(-6)}

. Как исправить это?
еще не получается простой расчет константы

(Оффтоп)

In[53]:= H = 7.2*10^5;
\[Tau] = 47*10^-6;
f = 3.5*10^8;
d = 1.2;
\[Lambda] = 0.022;
n = 3;
P0 = 25;
\[Eta] = 0.8;
KIP = 0.4;
kf = Sqrt[0.5];

In[63]:= (P0*(
 4*\[Pi] (\[Pi]*
   d^2)/^4)/\[Lambda]^2*\[Lambda]^2*f*\[Tau]*KIP*\[Eta]*kf^2)/(64*\
\[Pi]^2*H^2)

Out[63]= 1.14236*10^-8/Null^4

[salang]

еще не получился интеграл, который точно берется

(Оффтоп)

c = 3*10^8;
H = 7.2*10^5;
\[Tau] = 47*10^-6;
Fd = 7200;
Tp = 50*10^-6;
f = 3.5*10^8;
v = 7500;
Q = 18000;
L = 100;
s2 = 0.003;
s1 = 2;
\[CapitalTheta]0 = 0.018;
m = 256;
\[Lambda] = 0.022;
n = 3;
lm = 25;
d = 1.2;
\[Lambda] = 0.022;
n = 3;
P0 = 25;
\[Eta] = 0.8;
KIP = 0.4;
kf = Sqrt[0.5];
\[CapitalTheta]s = 1;
t = Sort@{-10^-8, -20^-9, -10^-9, 0, 10^-11, 10^-10, 10^-9, 10^-8, 
    20*10^-7, 200*10^-7, 10^-6};
z = Sum[NIntegrate[
    Exp[-(5.55* r^2)/(H^2 *\[CapitalTheta]0^2) - (11.1*
         r*\[CapitalTheta]s^2*Cos[\[Phi]])/(H *\[CapitalTheta]0^2) - 
      r^2/(H^2 *s1^2) - \[Pi]*
       f^2*(t - k Tp - r^2/(c H) - \[Eta])^2], {\[Eta], -Infinity, 
     Infinity}, {r, 0, Infinity}, {\[Phi], 0, 2*\[Pi]}, 
    MaxRecursion -> 40, AccuracyGoal -> 60, 
    Method -> "AdaptiveMonteCarlo"], {k, -255, 255}];

ListPlot[Transpose[{t, z}], Joined -> True, PlotRange -> All]

[Leierkastenmann]

1. Какой именно график строится до 20 мкс? У меня все графики строго до 2 мкс строятся в соответствии с тем, как задано t.
2. Расчет константы не получается, потому что формула набрана с ошибками, там какое-то очень странное возведение в четвертую степень.
3. Интеграл точно берется в каком смысле? Вы точно знаете, что он берется или он может быть вычислен аналитически? В общем-то он и в примере тоже берется, другое дело, что значение получается 0, видимо это Вас не устраивает, а чего ему не быть нулем, когда там у эскпоненты просто какой-то безумный показатель степени, причем "константно-безумный", такой что никаких надежд на отличие интеграла от 0 нет.
4. Про непрерывный случай забудьте. Сглаживайте график просто добавлением точек по t до той поры, пока Вам не надоест ждать пока расчет завершится.

[salang]

Leierkastenmann в сообщении #817124 писал(а):

2. Расчет константы не получается, потому что формула набрана с ошибками, там какое-то очень странное возведение в четвертую степень

а что именно странно? Площадь круга возводится в квадрат.

Leierkastenmann в сообщении #817124 писал(а):

3. Интеграл точно берется в каком смысле? Вы точно знаете, что он берется или он может быть вычислен аналитически? В общем-то он и в примере тоже берется, другое дело, что значение получается 0, видимо это Вас не устраивает, а чего ему не быть нулем, когда там у эскпоненты просто какой-то безумный показатель степени, причем "константно-безумный", такой что никаких надежд на отличие интеграла от 0 нет

он аналитически вычисляется (интегралы табличные). Там есть четкий максимум и плавный спад

Leierkastenmann в сообщении #817124 писал(а):

4. Сглаживайте график просто добавлением точек по t до той поры, пока Вам не надоест ждать пока расчет завершится.

имеется в виду указание

PlotPoints -> 30

?

[Aritaborian]

salang в сообщении #817144 писал(а):

имеется в виду указание

PlotPoints -> 30

?

Нет. Сказано же было:

Leierkastenmann в сообщении #817124 писал(а):

Про непрерывный случай забудьте.

Имеется в виду увеличение количества точек в дискретном случае.