I have the following long functions f1and f2 dependent on parameters ω and g1:

f1 = (1.92`*^8 g1^2)/(10000 (1 + 40000 ω^2) + (400 + 4 ω^2) (1 + 160000 g1^4 + 40004 ω^2 + 160000 ω^4 - 800 g1^2 (-1 + 400 ω^2)) + 200 (-4 ω^2 (1 - 40000 g1^2 + 40000 ω^2) + 20 (1 + 400 g1^2 + 40000 ω^2))) + (96000.` g1^2 (400 + 4 ω^2))/(10000 (1 + 40000 ω^2) + (400 + 4 ω^2) (1 + 160000 g1^4 + 40004 ω^2 + 160000 ω^4 - 800 g1^2 (-1 + 400 ω^2)) + 200 (-4 ω^2 (1 - 40000 g1^2 + 40000 ω^2) + 20 (1 + 400 g1^2 + 40000 ω^2))) + (120200.` (10000 + 200 (20 - 4 ω^2) + (1 + 4 ω^2) (400 + 4 ω^2)))/(10000 (1 + 40000 ω^2) + (400 + 4 ω^2) (1 + 160000 g1^4 + 40004 ω^2 + 160000 ω^4 - 800 g1^2 (-1 + 400 ω^2)) + 200 (-4 ω^2 (1 - 40000 g1^2 + 40000 ω^2) + 20 (1 + 400 g1^2 + 40000 ω^2)))
f2 = (7.68`*^8 g1^2)/(160000 (1 + 40000 ω^2) + (400 + 4 ω^2) (1 + 160000 g1^4 + 40004 ω^2 + 160000 ω^4 - 800 g1^2 (-1 + 400 ω^2)) + 800 (-4 ω^2 (1 - 40000 g1^2 + 40000 ω^2) + 20 (1 + 400 g1^2 + 40000 ω^2))) + (96000.` g1^2 (400 + 4 ω^2))/(160000 (1 + 40000 ω^2) + (400 + 4 ω^2) (1 + 160000 g1^4 + 40004 ω^2 + 160000 ω^4 - 800 g1^2 (-1 + 400 ω^2)) + 800 (-4 ω^2 (1 - 40000 g1^2 + 40000 ω^2) + 20 (1 + 400 g1^2 + 40000 ω^2))) + (120200.` (160000 + 800 (20 - 4 ω^2) + (1 + 4 ω^2) (400 + 4 ω^2)))/(160000 (1 + 40000 ω^2) + (400 + 4 ω^2) (1 + 160000 g1^4 + 40004 ω^2 + 160000 ω^4 - 800 g1^2 (-1 + 400 ω^2)) + 800 (-4 ω^2 (1 - 40000 g1^2 + 40000 ω^2) + 20 (1 + 400 g1^2 + 40000 ω^2)))

I intend to find the integral of f1 and f2 with respect to ω for varying g1 values. Then doing a Log-Log plot of the resulting integrand against g1:

test1 = Table[{g1, 1/(2 π)*NIntegrate[Evaluate[f1], {ω, -100, 100}]}, {g1, 0, 150}];
test2 = Table[{g1, 1/(2 π)*NIntegrate[Evaluate[f2], {ω, -100, 100}]}, {g1, 0, 150}];

ListLogLogPlot[{test1, test2}, PlotRange -> All, Joined -> True, PlotLegends -> {"test1", "test2"}]

I am returned with: enter image description here

First of all it is important to me that the graph shows that they both start at the same point (which they really do). But it is clear that when the x-axis (g1) starts at 1, the two curves starts at two different y values. I was thinking maybe this is because it is not possible to show Log[0]. I tried to remedy the issue by picking a g1 value close to 0 and replot it:

enter image description here Clearly they do start from the same value when g1=0.01 but there appears to be a kink around g1=1. In general, I just want a Log-Log graph that shows that both curves starts at the same point when the x-axis is 0 How should I approach this?

Second of all, it appears that beyond g1=100 the curve drops suddenly to about 0.01. I think it might be numerical instabilities but I'm not sure what's going on here. If it were numerics, I should have had warnings but I didn't see any. What is going on here?

  • I notice I can change the position of the "kink" at the right end of your graphs by changing your {ω, -100, 100} to {ω, -80, 80} or to {ω, -140, 140}. There are some warnings about convergence, but the "kink" appears to track that ω bound exactly. Maybe this is enough of a hint to figure out why the "kink" is there. – Bill Oct 10 at 22:20
  • @Bill Thanks for the help. I still need to show both curves to be starting from the same value smoothly in the Log plot. Any ideas? – kowalski Oct 11 at 0:53
  • I can't be certain about the instability, but Plot3D[f1, {ω, -100, 100}, {g1, 0, 150}, PlotRange -> {0, 10}] seems to show that your expressions appear to blow up for certain values of g1 and ω. I'm guessing that is because of your large complicated denominators going to zero. Simplify[Solve[ Denominator[Simplify[f1]] == 0, {ω, g1}]] might give you some hints about why that is happening. – Bill Oct 11 at 1:12
up vote 2 down vote accepted

You can get a smoother result by calculating the integrals analytically (though not in closed form) first:

aux1[g1_, \[Omega]_] = Integrate[f1, \[Omega]];
aux2[g1_, \[Omega]_] = Integrate[f2, \[Omega]];

int1[g1_] = aux1[g1, 100] - aux1[g1, -100];
int2[g1_] = aux2[g1, 100] - aux2[g1, -100];

LogLogPlot[{1/(2 \[Pi]) int1[g], 1/(2 \[Pi]) int2[g]}, {g, 0.01, 150}]

enter image description here

Your Answer

By clicking "Post Your Answer", you acknowledge that you have read our updated terms of service, privacy policy and cookie policy, and that your continued use of the website is subject to these policies.

Not the answer you're looking for? Browse other questions tagged or ask your own question.