How can I plot b[a]/b[0.0001] vs. a using Mathematica?

  • Thread starter Thread starter kptsilva
  • Start date Start date
  • Tags Tags
    Mathematica Plot
AI Thread Summary
The discussion revolves around plotting the function b[a]/b[0.0001] versus a using Mathematica, based on a differential equation from a research paper. The equation involves a boundary condition b'[0.0001]=0 and a specific function w[a] defined in terms of a and a variable c, which is set to 1 for this instance. Participants clarify that the equation should be solved using NDSolve, not NIntegrate, as it is a second-order differential equation requiring two integrations. A critical correction was made regarding a negative sign in the definition of w[a], which significantly impacted the results. The provided Mathematica code successfully generates the desired plot, although some issues with different values of c remain unresolved.
kptsilva
Messages
33
Reaction score
0
hey guys,
I came across a research paper stating to numerically integrate the following equation.

2/3 a^2 b''[a] + (1 - w[a]) a b'[a] - (1 + w[a]) (1 - 3 c w[a]) b[a]

A Boundary condition is given b'[0.0001]=0
Where;

w[a_] := 2*a^(3*(1 + c))/(1 + 2*a^(3*(1 + c)));

(c=1) (c is a variable but let's consider a particular instance c=1)

'a' goes from 10^-4 to 1000 in a log scale.

I want to plot b[a]/b[0.0001] vs. a

how can i plot this using mathematica?
 
Physics news on Phys.org
2/3 a^2 b''[a] + (1 - w[a]) a b'[a] - (1 + w[a]) (1 - 3 c w[a]) b[a] == ?
 
?==0
 
See the attached notebook. The second plot on a log scale is probably more instructive.
 

Attachments

the paper says that they have numerically integrated that equation.I can't understand that statement. The program should be written with NIntegrate rather than NDSolve right?
 
Since it's a second order equation, you would need to integrate it twice. Integrating along step-by-step is really what NDSolve is doing. I don't see a simple way to write a solution using the NIntegarate command. Does it matter?
 
Figures 2 and 3 are the same graph, just plotted on different scales. They've numerically solved the differential equation just like I did. Also, in your original post, you forgot the minus sign in w[a]. This makes a big difference! Now it looks like the paper. See attached.
 

Attachments

phyzguy thanks alot, in my program i have forgotten the negative sign it gives me the graphs perfectly.
thhis is the code i wrote,
c = 1;
w[a_] := (2*a^(3*(1 + c)))/(1 + (2*a^(3*(1 + c))));
fun = (2/3)*(a^2)*b''[a] + (1 - w[a])*a*
b'[a] - (1 + w[a])*(1 - 3*c*w[a])*b[a]; sol =
DSolve [{fun == 0, b'[10^-4] == 0}, b, a];
A = LogLinearPlot[Evaluate[b[a]/b[10^-4] /. sol], {a, 10^-4, 10^3},
PlotRange -> {-500, 12000}]
 
  • #10
There are still some problems, i think there is abit of an error for some values of c(=alpha)
 
  • #11
for c=1 the graph is perfect
 
  • #12
Any news to my problem?
 
Back
Top