# How to solve a first order coupled system of ODEs?

How can I solve this coupled differential equation?

ClearAll[x, y, g, w1, w2]
n = 1;b = 1;c = 1;
w1 = -(1/3) - 2 Sqrt[g[x]]/3;
w2 = -(b*b)*(1 + y[x])^n/(g[x])^(n - 1)*(1 - g[x] + y[x]);
g'[x] == -3*g[x]*(1 - g[x] + y[x])*(w1 - w2) + y[x]*g[x]*(1 + 3*w1);
y'[x] == -3*y[x]*(1 - g[x] + y[x])*(w1 - w2) +
y[x]*(1 + y[x])*(1 + 3*w1);

• What are the initial conditions?
– zhk
Commented Jan 8, 2017 at 9:39
• with the initial conditions g[0] == 0.72, y[0] == 0.01 Commented Jan 8, 2017 at 9:41
• you should make that comment about initial conditions a part of your question. You should be able to edit your question by clicking on the edit link right below it... Commented Jan 8, 2017 at 9:48

You can use NDSolve to solve your system of odes numerically,

n = 1; b = 1; c = 1;
w1 = -(1/3) - 2 Sqrt[g[x]]/3;
w2 = -(b*b)*(1 + y[x])^n/(g[x])^(n - 1)*(1 - g[x] + y[x]);
Eq1 = g'[x] == -3*g[x]*(1 - g[x] + y[x])*(w1 - w2) +
y[x]*g[x]*(1 + 3*w1)
Eq2 = y'[x] == -3*y[x]*(1 - g[x] + y[x])*(w1 - w2) +
y[x]*(1 + y[x])*(1 + 3*w1)
sol = NDSolve[{Eq1, Eq2, g[0] == 0.72, y[0] == 0.01}, {g, y}, {x, -3,
3}]
Plot[Evaluate[{g[x], y[x]} /. sol], {x, -3, 10}, PlotStyle -> Thick]


To plot g[x] vs y[x]

ParametricPlot[{g[x], y[x]} /. sol, {x, -3, 3}, PlotRange -> All,
MaxRecursion -> 8, AxesLabel -> {"g", "y"}]


• thank you so much for your help Commented Jan 8, 2017 at 9:50
• @merve Sure. ParametricPlot[{g[x], y[x]} /. sol, {x, 0, 1}, PlotRange -> All, MaxRecursion -> 8, AxesLabel -> {"g", "y"}]
– zhk
Commented Jan 8, 2017 at 16:05
• @merve Plot[Evaluate[{g[x], y[x]} /. sol], {x, -3, 10}, PlotStyle -> Thick]
– zhk
Commented Jan 8, 2017 at 16:45
• @merve PlotLegends -> {"g", "y"}.
– zhk
Commented Jan 8, 2017 at 17:04
• @merve data = Plot[Evaluate[{g[x], y[x]} /. sol], {x, -3, 10}, PlotStyle -> Thick]points = Cases[Normal@data, Line[pts_, ___] :> Flatten[pts, Depth[pts] - 3], Infinity]; ifuncts = Interpolation[#, Method -> "Spline", InterpolationOrder -> 2][x] & /@ points; data = Table[Prepend[ifuncts, x], {x, -3, 10, .1}]; TableForm[data, TableHeadings -> {None, Prepend[Array["X", Length[points]], x]}]Export["C:/tcdata/myfile.txt", data, "Table"]
– zhk
Commented Jan 28, 2017 at 18:27