Take the 2-minute tour ×
Mathematica Stack Exchange is a question and answer site for users of Mathematica. It's 100% free, no registration required.

Possible Duplicate:
Kramers-Kronig in Mathematica

I am trying to get epsilon2 (imaginary part) from known epsilon1 (real part) data with the respective energy values (w), using the Kramers-Kronig relations. My code is:

{epsilon1, w} = ToExpression@Import["http://pastebin.com/raw.php?i=Z08GwMa8"];
epsilon2[w_]:= 2/Pi NIntegrate[a epsilon1[a]/(a^2-w^2), {a, 0, 5.5793}, 
    Method-> "PrincipalValue", MaxRecursion->20, Exclusions->{(a^2-w^2)==0}]//Quiet

Plot[epsilon2[w],{w,0,5.57932},AxesOrigin->{0,0},PlotPoints->400]

However, by using the code above, I managed only to get a list of irrelevant numbers. The graph of epsilon1 looks like a sine curve, while that of epsilon2 looks like an "M" shape. Any help is greatly appreciated! =)

P.S. If pastebin is down, you can get the original data from revision 1.

share|improve this question
add comment

marked as duplicate by whuber, Oleksandr R., rcollyer, rm -rf Oct 27 '12 at 3:17

This question has been asked before and already has an answer. If those answers do not fully address your question, please ask a new question.

1 Answer

up vote 3 down vote accepted

epsilon1as you defined it, isn't a function but a list.
You may try something like:

{epsilon1, w} = ToExpression@Import["http://pastebin.com/raw.php?i=Z08GwMa8"];
f = Interpolation[Transpose[{Flatten[w], Flatten[epsilon1]}]]
epsilon2[w_] := epsilon2[w] = 
                2/Pi NIntegrate[a f[a]/(a^2 - w^2), {a, 0, 5.5793}, 
                Method -> "PrincipalValue", Exclusions -> {(a^2 - w^2) == 0}] //  Quiet
DiscretePlot[epsilon2[w], {w, 0.5, 5.5, .25}, AxesOrigin -> {0, 0}, Joined -> True]

Mathematica graphics

share|improve this answer
    
Hi @belisarius! Thank you =) 2 things: 1) The shape of epsilon2 is there, but it is supposed to be above the x-axis (i.e. flipping the negative region up to the positive region). 2) From my epsilon2 data, the curve should start off from the origin {0,0} rather than from the high values. Is there a problem or something? =) –  weesiang Oct 16 '12 at 2:27
    
@weesiang You may have to check them. That curve is a reasonable approach to your integrand. –  belisarius Oct 16 '12 at 10:22
1  
Oh yes!! I miss out a parameter in the equation! Now I got it! Thank you very much @belisarius! –  weesiang Oct 16 '12 at 10:47
    
Hi @belisarius, If I were to integrate from 0 to infinity, how would you alter the syntax? And also, I would like to ask you if there are alternative commands to DiscretePlot? =) –  weesiang Dec 6 '12 at 11:20
add comment

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